# Evidence for a mass-dependent AGN Eddington ratio distribution via the flat relationship between SFR and AGN luminosity

Evidence for a mass-dependent AGN Eddington ratio distribution via the flat relationship between... Abstract The lack of a strong correlation between AGN X-ray luminosity (LX; a proxy for AGN power) and the star formation rate (SFR) of their host galaxies has recently been attributed to stochastic AGN variability. Studies using population synthesis models have incorporated this by assuming a broad, universal (i.e. does not depend on the host galaxy properties) probability distribution for AGN specific X-ray luminosities (i.e. the ratio of LX to host stellar mass; a common proxy for Eddington ratio). However, recent studies have demonstrated that this universal Eddington ratio distribution fails to reproduce the observed X-ray luminosity functions beyond z ∼ 1.2. Furthermore, empirical studies have recently shown that the Eddington ratio distribution may instead depend upon host galaxy properties, such as SFR and/or stellar mass. To investigate this further, we develop a population synthesis model in which the Eddington ratio distribution is different for star-forming and quiescent host galaxies. We show that, although this model is able to reproduce the observed X-ray luminosity functions out to z ∼ 2, it fails to simultaneously reproduce the observed flat relationship between SFR and X-ray luminosity. We can solve this, however, by incorporating a mass dependency in the AGN Eddington ratio distribution for star-forming host galaxies. Overall, our models indicate that a relative suppression of low Eddington ratios (λEdd ≲ 0.1) in lower mass galaxies (M* ≲ 1010 − 11 M⊙) is required to reproduce both the observed X-ray luminosity functions and the observed flat SFR/X-ray relationship. galaxies: active, galaxies: evolution, galaxies: luminosity function, mass function, galaxies: statistics 1 INTRODUCTION Most galaxies host a central super massive black hole (SMBH), the masses of which display a tight correlation with their host bulge masses, implying a co-evolution between SMBHs and their host galaxies (e.g. Kormendy & Richstone 1995; Magorrian et al. 1998; Merritt 2000; Kormendy, Bender & Cornell 2011). Furthermore, the similar redshift evolution of the total SMBH accretion rate density and the total star formation rate (SFR) density strongly suggests the existence of a link between SMBH accretion and star formation at galactic scales (e.g. Heckman et al. 2004; Merloni, Rudnick & Di Matteo 2004; Silverman et al. 2008; Aird et al. 2010; Assef et al. 2011). However, the precise connection between the two processes remains poorly understood. One means of connecting SMBH and galaxy growth that has gained popular support over the past two decades is for accreting SMBHs (observed as active galactic nuclei, AGNs) to directly influence the SFRs of their hosts via a variety of feedback mechanisms (see the review of Fabian 2012, for details on the feedback mechanisms). However, a major difficulty in identifying the precise role that SMBH accretion has on influencing SFR is the stochastic nature of AGNs (e.g. Aird et al. 2013; Hickox et al. 2014, see also section 4.3 of Stanley et al. 2015). It is argued that this randomness is the reason why most studies that have explored the relationship between SFR and X-ray luminosity (a proxy for SMBH accretion power) have reported no evidence of a strong correlation between these parameters, at least for moderate luminosity AGNs, which form the majority of the population (i.e. 1042 < LX < 1045 erg s−1; e.g. Lutz et al. 2010; Harrison et al. 2012; Mullaney et al. 2012; Rosario et al. 2012; Rovilos et al. 2012; Santini et al. 2012; Azadi et al. 2015; Stanley et al. 2015; Lanzuisi et al. 2017). Recently, using optically selected AGNs, Stanley et al. (2017) reported an enhancement of SFR among the highest luminosity AGNs (i.e. LX >1045 erg s−1). However, they demonstrated that this is a direct consequence of the most luminous AGNs residing in more massive host galaxies, meaning the enhanced SFR is a consequence of the relationship between stellar mass and SFR (hereafter referred to as the main sequence, MS; e.g. Daddi et al. 2007; Rodighiero et al. 2011; Schreiber et al. 2015). A major problem in investigating the relationship between AGN and SFR is that observations are often subjected to biases (e.g. flux limited samples, rarity of high accretion rate AGNs). Therefore, many studies investigating the connection between SMBH and galaxy growth have adopted a modelling approach (e.g. Aird et al. 2013; Conroy & White 2013; Veale, White & Conroy 2014; Caplar, Lilly & Trakhtenbrot 2015; Jones et al. 2017; Weigel et al. 2017). In these models, the stochastic nature of AGNs is often implemented via a broad specific AGN X-ray luminosity (i.e. the X-ray luminosity relative to the host stellar mass; LX/M*) distribution that, once combined with a mass function, reproduces the observed X-ray luminosity function of AGNs. Overall, these studies find that below z ∼ 1 the specific X-ray luminosity distribution is well-represented by a broad, universal (i.e. the same distribution, irrespective of SFR or stellar mass) distribution, described by a broken power-law (or similar; e.g. a Schechter 1976 function), the normalization of which increases with increasing redshifts (e.g. Aird et al. 2013; Veale et al. 2014; Jones et al. 2017). These models are also often successful at reproducing the observed flat relationship between SFR and X-ray luminosity (Hickox et al. 2014; Veale et al. 2014; Stanley et al. 2015). However, Jones et al. (2017) recently found that using a universal broad distribution for the specific X-ray luminosity distribution of AGNs cannot reproduce the X-ray luminosity functions of Aird et al. (2010) beyond z ∼ 1.2. They argue for the need of a more complicated specific X-ray luminosity distribution of AGNs beyond that redshift. Furthermore, recent studies have reported that the specific X-ray luminosity distribution for AGNs in star-forming hosts significantly differs from that of AGNs in quiescent hosts (e.g. Georgakakis et al. 2014; Wang et al. 2017; Aird et al. 2018). Overall, they report a higher contribution from star-forming AGN hosts to the total specific X-ray luminosity distribution out to z ∼ 2. Motivated by the above findings, we explore a model for the specific X-ray luminosity distribution split between star-forming and quiescent galaxies. After demonstrating that this simple model cannot simultaneously reproduce both the X-ray luminosity functions and the flat relationship between SFR and X-ray luminosity out to z ∼ 2, we introduce a mass dependency for the specific X-ray luminosity distribution of star-forming AGN hosts, as recently observed (e.g. Bongiorno et al. 2016; Georgakakis et al. 2017; Aird et al. 2018). Hereafter, we refer to our model specific X-ray luminosity distribution split between star-forming and quiescent galaxies as ‘mass-independent’, and to our model that includes the mass dependency for the specific X-ray luminosity distribution of AGNs in star-forming galaxies as ‘mass-dependent’. In Section 2, we describe our method to infer the specific X-ray luminosity distribution for both the mass-independent and the mass-dependent model. Subsequently, in the same section, we present how we incorporate host galaxy properties. We then report the results for the mass-independent model in Section 3.1, and for the mass-dependent model in Section 3.2. We discuss the general implications of our results in Section 4, and conclude in Section 5. Throughout, we assume the specific X-ray luminosity (i.e. LX/M*) to be proportional to Eddington ratio (i.e. the ratio between bolometric AGN luminosity and Eddington luminosity, LAGN/LEdd), converting the X-ray 2–10 keV luminosity into bolometric AGN luminosity using a universal conversion factor of 22.4 (median value found in Vasudevan & Fabian 2007, and based on a local AGN sample with LX = 1041 − 46 erg s−1), and assuming a ratio of 0.002 between SMBH masses and stellar masses (Marconi & Hunt 2003). Therefore, we have,   $${\rm \lambda _{\rm Edd}} = \frac{22.4 \times {\rm L_{\rm X}}}{10^{38.1}\,{\rm erg\,s^{-1}} \times 0.002\ \frac{{\rm M_\ast }}{{\rm \mathrm{M}_{\odot }}}},$$ (1)where λEdd is the Eddington ratio, and LX is the intrinsic (absorption-corrected) 2–10 keV X-ray luminosity. We stress that the important parameter we investigate is the specific X-ray luminosity (LX/M*), and that the Eddington ratio is solely used for its greater familiarity. Throughout, we assume a WMAP–7 year cosmology (Larson et al. 2011) and a Salpeter (1955) initial mass function (IMF) when calculating galaxy stellar masses and SFRs. 2 METHODOLOGY Many recent studies have used Population Synthesis Models (PSMs) to investigate the connection between AGNs and host galaxies. PSMs randomly draw properties from probability distribution functions (e.g. Aird et al. 2013; Hickox et al. 2014; Stanley et al. 2015; Jones et al. 2017). Using PSMs, we are able to generate populations of N massive mock galaxies, with stellar masses M*, drawn from a mass distribution, for which we randomly allocate Eddington ratios, λEdd, following an Eddington ratio distribution, and, using equation (1), derive X-ray luminosities, LX, the histogram of which is proportional to the model X-ray luminosity function. The key point is to optimize the Eddington ratio distribution until the model X-ray luminosity function matches the observed one (e.g. Aird et al. 2013; Jones et al. 2017). This can be repeated at different redshifts to investigate any redshift dependence of the Eddington ratio distribution. Because it implies generating a new population of galaxies while optimizing the Eddington ratio distribution, it can be computationally expensive. Therefore, we will use an analytical approach (similar to those of e.g. Veale et al. 2014; Caplar et al. 2015) by directly considering the various distributions as probability density functions. In this study, we fit to the X-ray luminosity functions of Aird et al. (2015) and use the stellar mass functions of Davidzon et al. (2017). The latter are shown in Fig. 1 for various redshifts, out to z = 2.5. Figure 1. View largeDownload slide Mass functions employed in our models. Top panel: mass functions for star-forming galaxies as reported in Davidzon et al. (2017). Each colour corresponds to a different redshift bin, out to z = 2.5 (see key). Bottom panel: same as top panel but for the quiescent galaxies. Figure 1. View largeDownload slide Mass functions employed in our models. Top panel: mass functions for star-forming galaxies as reported in Davidzon et al. (2017). Each colour corresponds to a different redshift bin, out to z = 2.5 (see key). Bottom panel: same as top panel but for the quiescent galaxies. In our analytical approach, we first assume a model λEdd distribution, which is described by a set of parameters θ (for example, for a power law θ = {N, α}, where N is the normalization and α the slope). Our aim is to iterate over θ to identify the solution that best fits the observed X-ray luminosity function. The model X-ray luminosity function is generated by analytically combining the observed mass function with the λEdd distribution given by the set of parameters θ. To achieve this, we treat the X-ray luminosity function as a probability distribution, p(LX). The probability of observing an AGN of luminosity LX = X in a galaxy of mass M* is given by1  $$p(L_{\rm X}=X) \equiv p\left(\lambda _{\rm Edd} = \frac{X}{M_\ast }\right),$$ (2)where p(λLedd = X/M*) is the probability of λEdd = X/M* calculated with the current set of model parameters, θ. We derive the total probability of observing LX = X in all galaxies as the sum of these individual probabilities weighted by the galaxy mass function (which we also treat as a probability distribution, p(M*)),   \begin{eqnarray} p(L_{\rm X}=X) = \sum _{Y=M_{\rm Min}}^{M_{\rm Max}} p\left(\lambda _{\rm Edd} = \frac{X}{Y} | M_\ast =Y\right)\times p(M_\ast =Y),\nonumber \\ \end{eqnarray} (3)where p(M* = Y) is the probability of M*=Y (or the mass function evaluated at Y). The total model luminosity function derived from a given set of parameters, θ, is then simply the combination of these probabilities (i.e. for LX = 1041–1046 erg s−1 in bins of 0.1 dex). By taking this approach, we avoid having to model N galaxies (where $$N\sim \mathcal {O}(10^9)$$ to generate sufficient numbers of luminous galaxies), and can simply split our mass function into $$\mathcal {O}(10^3)$$ bins between MMin = 108 M⊙ to MMax = 1014 M⊙ (where MMin and MMax are chosen to cover the full range of stellar masses). Since we adopt an iterative process to identify the best-fitting parameters (see later in this subsection), this factor of $$\mathcal {O}(10^{6})$$ reduction in the number of required calculations per iteration dramatically reduces the time taken to converge. To optimize the parameters that define our Eddington ratio distribution, we performed maximum likelihood estimation. Assuming Gaussian distributions for the uncertainties, the log-likelihood of a set of parameters θ is defined by −0.5χ2, where χ2 is the chi-square between the model and the observed X-ray luminosity functions. We used flat proper (i.e. with finite boundaries) prior distributions for each parameter that defines our model, and checked the posterior distributions to ensure that they are not constrained in any way by the limits of our prior distributions. The parameter space (defined by the flat prior distributions) was explored using a Monte Carlo Markov Chain (MCMC) fully implemented in a python application programming interface, emcee2 (Foreman-Mackey et al. 2013), that uses the affine invariant MCMC ensemble sampler of Goodman & Weare (2010). We ran our MCMC code using 100 walkers for 10 000 steps as a burn-in phase to find the global solution (i.e. these 10 000 steps are then discarded), and for 5000 steps as a production phase (i.e. these 5000 steps are kept for analysis). The convergence of each chain was assessed using the Gelman & Rubin (1992) test. We then extracted the best parameters by fitting a Gaussian function to each sampled posterior distribution, taking the mean and the standard deviation as the best estimate and the 1σ uncertainties, respectively, of each parameter. When the posterior distribution of a parameter was flat below or above a given value, we assumed that it is an upper or a lower limit, respectively (again, we ensured that values of these limits are not affected by our prior distributions). This optimization method was applied to infer the Eddington ratio distribution for both our mass-independent and mass-dependent models. 2.1 A mass-independent λEdd distribution As highlighted in the introduction, a number of studies have recently reported a difference in the λEdd distribution of quiescent and star-forming host galaxies (e.g. Georgakakis et al. 2014; Wang et al. 2017). We therefore take this ‘two-component’ distribution as our starting point (i.e. we do not attempt to model a single universal λEdd distribution for all galaxies). Most studies use a broken power-law (or related functions) to represent the Eddington ratio distribution of AGNs (e.g. Aird et al. 2013; Veale et al. 2014; Jones et al. 2017; Weigel et al. 2017). As such, we take a conservative approach by assuming a broken power-law for both the star-forming and the quiescent component of our λEdd distribution. Each broken power law is defined as   $$p(\lambda _{\rm Edd}) = \frac{A}{\left(\frac{\lambda _{\rm Edd}}{\lambda _{\rm break}} \right)^{\gamma _{1}}+\left(\frac{\lambda _{\rm Edd}}{\lambda _{\rm break}} \right)^{\gamma _{2}}},$$ (4)where $$p({\rm \lambda _{\rm Edd}})$$ is the probability of Eddington ratio, A is the normalization, λbreak is the position of the break, and γ1 and γ2 are the slopes at low λEdd (i.e. below the break) and high λEdd (i.e. above the break), respectively. As such, each of our broken power law components are represented by four parameters, giving a total of eight free parameters to optimize for this mass-independent model (see Section 2 for the optimization method). Following empirical results reported in, for example, Georgakakis et al. (2014), Wang et al. (2017), and Aird et al. (2018), prior to optimization, we assume that both the normalization, ASF, and the position of the break, $$\lambda _{\rm break}^{\rm SF}$$, of the star-forming component are each always higher than their quiescent analogues, AQui and $$\lambda _{\rm break}^{\rm Qui}$$, respectively. This also has the benefit of reducing the degeneracy of the model. However, it is important to stress that we do not assume any specific parameter values for the ratio between the two normalizations and the two positions of the breaks, we simply incorporate this information in the prior distributions by excluding some parts of the parameter space. Incorporating the aforementioned assumptions, we optimized the eight free parameters of our mass-independent model. To investigate any redshift evolution of the Eddington ratio distribution, we performed the optimization at z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.2, z = 1.5, z = 1.7, z = 2.0, z = 2.3, z = 2.5, z = 2.5, z = 2.7, z = 3.0, and z = 3.5, interpolating the mass functions of Davidzon et al. (2017) split between star-forming and quiescent galaxies, and the total X-ray luminosity functions of Aird et al. (2015) at these redshifts. We stress that the Bayesian methodology adopted in Aird et al. (2015) to derive the observed X-ray luminosity functions (see Aird et al. (2015) for details on their method) means that they consider all X-ray detected sources when deriving their X-ray luminosity functions (which our analysis aims to fit). Therefore, there is no implied mass cut beyond the fact that their sources are detected in X-rays. 2.2 A mass-dependent λEdd distribution Motivated by recent studies reporting that the λEdd distribution for star-forming galaxies is mass-dependent (e.g. Bongiorno et al. 2016; Georgakakis et al. 2017; Aird et al. 2018), we develop a second model that can accommodate three different λEdd distributions in three different mass bins for star-forming galaxies. As with our mass-independent model, we still split in terms of quiescent and star-forming galaxies, but for the star-forming component, we assume instead three different broken power-laws (see equation 4), one for each of our mass bins, which we define as: low mass (8 < log(M*/M⊙) < 10), medium mass (10 < log(M*/M⊙) < 11), and high mass (11 < log(M*/M⊙) < 12). The choice of the boundaries for our mass bins is arbitrary and independent of any empirical studies, yet aims to cover a large range of stellar masses (i.e. 108 <M*/M⊙ <1012). Since assuming three broken power-laws generates 12 free parameters, we again rely on some assumptions based on the results of Aird et al. (2018) to help prevent degeneracies. First, we require that each broken power-law breaks at different λEdd values, in such a way that a higher mass bin breaks at a lower value of λEdd than its lower mass neighbour. Secondly, we assume that the three Eddington ratio distributions share the same slope at high Eddington ratios. Finally, we assume the normalization of each of the three broken power-laws is such that the λEdd distribution above the break is always coincident (see Fig. 8). In making these assumptions, we reduce the parameter space to eight free parameters: one normalization, ASF, (from which the others are derived via our third assumption); three break positions, $$\lambda _{\rm break}^{\rm SF}$$, (one for each mass bin and ordered according to our first assumption); three power-law slopes below the break, $$\gamma _{1}^{\rm SF}$$, (which are unconstrained); a single shared power-law slope above the break, $$\gamma _2^{\rm SF}$$. We stress that all these assumptions are inspired by empirical findings reported by Aird et al. (2018). We will demonstrate in Section 3.1.1 that our previous model (i.e. the mass-independent model, see Section 2.1) is able to reproduce the observed X-ray luminosity functions out to z ∼ 2, and recreates both their star-forming and quiescent galaxy components in good agreement with observations, at least out to z ∼ 1. As such, for this mass-dependent model, we assume that the Eddington ratio distribution for quiescent galaxies is unchanged from our previous model, and only update the star-forming component of the model Eddington ratio distribution by implementing the aforementioned mass-dependency. To do this, we optimize the eight free parameters that describe our mass-dependent Eddington ratio distribution for star-forming galaxies. Again, to investigate any redshift evolution of the mass-dependent Eddington ratio distribution, we performed the optimization at z = 0.1, z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.3, z = 1.5, and z = 1.7, interpolating the mass functions for star-forming galaxies of Davidzon et al. (2017) and the star-forming component of the X-ray luminosity functions as extracted from our mass-independent model. This second model relies on the ability of our first, mass-independent, model to split the X-ray luminosity functions in terms of star-forming and quiescent galaxies. As we show in Section 3.1.1, our mass-independent model is unable to reproduce the total X-ray luminosity function beyond z ∼ 2, and thus cannot be used to reliably split the luminosity function between star-forming and quiescent galaxies. As such, we cannot extend our mass-dependent model beyond z ∼ 2. 2.3 Generating galaxies with AGNs As we aim to investigate the connection between AGN accretion and host star formation activity, we must model the star-forming properties of the host galaxies as well as the Eddington ratios of the AGNs. Since the distributions of specific SFRs (i.e. SFR relative to the host stellar mass; sSFR) are now well-defined, we can simply adopt published models to generate our galaxy populations. As such, for this part of the model, we use a PSM to attribute SFRs to our host galaxies. We use the PSM outlined in Bernhard et al. (2014), which generates a population of mock star-forming galaxies using the mass functions of Ilbert et al. (2013) and allocates SFRs using the relationship between stellar masses and SFRs, or MS ( e.g. Daddi et al. 2007; Salim et al. 2007; Elbaz et al. 2011; Rodighiero et al. 2011; Sargent et al. 2012; Schreiber et al. 2015, for the MS) reported in Rodighiero et al. (2011). The SFRs are then randomly scattered around the sSFR distribution of Sargent et al. (2012). This model successfully reproduces the observed infrared and ultraviolet luminosity functions up to z ∼ 6 (see Bernhard et al. 2014, and references therein). However, for this work, we slightly update the empirical relationships used in Bernhard et al. (2014). Briefly, we change the mass functions to the more recent ones reported in Davidzon et al. (2017), we update the MS sSFR distribution to that reported in Schreiber et al. (2015), and lastly, add a population of quiescent galaxies (using the quiescent galaxy mass functions of Davidzon et al. 2017) with sSFRs at least a factor of ten below that of the MS (Ilbert et al. 2013). We confirm that this updated model still reproduces the infrared and ultraviolet luminosity functions. X-ray luminosities are then allocated using our Eddington ratio distributions derived in Section 2. Using this PSM, we generate a population of 33 925 192 galaxies (i.e. 32 939 834 – or 97 per cent – star-forming and 985 358 – or 3 per cent – quiescent galaxies), with stellar masses in the range 108 <M*/M⊙ <1014. This corresponds to a 50 square degrees blank-field survey out to z = 3. We use a very high (arguably un-physical) upper limit on our range of stellar masses to rule-out the possibility of being affected by any boundary effects. The large relative number of star-forming compared to quiescent galaxies in our model is related to the differences of their respective mass functions, particularly at masses M* ≲ 1010 M⊙ and z ≳ 1. As shown in Fig. 1, at these lower masses and higher redshifts, the mass function of star-forming galaxies steeply rises, while that of quiescent galaxies sharply decreases. This leads to a large difference in the absolute number of star-forming galaxies compared to that of quiescent galaxies at lower masses (we have a cut at M* = 108 M⊙ in our model). Above a stellar mass of M* >1010 M⊙, we have a total of 3 309 201 galaxies, among which 2 646 102 (80 per cent) are star-forming and 663 099 (20 per cent) are quiescent. 3 RESULTS In this section, we first present the results of our mass-independent model, followed by the results of our mass-dependent model. We demonstrate that while our mass-independent model can successfully reproduce the observed X-ray luminosity functions out to z ≈ 2, it cannot at the same time reproduce the observed flat relationship between SFR and AGN luminosity (e.g. Rosario et al. 2012; Stanley et al. 2015). Instead, we show that this can be resolved by introducing a mass dependency in the λEdd distribution of AGNs hosted in star-forming galaxies. 3.1 The mass-independent model 3.1.1 X-ray luminosity functions We show in Fig. 2, the fit to the X-ray luminosity functions for our mass-independent model. Up to z = 1.75 our model X-ray luminosity functions are in very good agreement with those of Aird et al. (2015). However, beyond that redshift, our mass-independent model is unable to reproduce the observed X-ray luminosity functions. Fig. 2 also illustrates the model X-ray luminosity functions split between star-forming and quiescent galaxies. Overall, we find that the star-forming galaxies always dominate the X-ray luminosity functions at LX ≳ 1042.5 erg s−1, which corresponds to the bulk of the AGN population. The contribution from quiescent galaxies in our mass-independent model, however, is only significant at lower X-ray luminosities (i.e. LX ≲ 1042 erg s−1). Figure 2. View largeDownload slide Fit to the X-ray 2–10 keV luminosity functions of Aird et al. (2015) when optimizing our mass-independent λEdd distribution model. Each panel represents a different redshift bin out to z = 3, and as indicated at their bottom right-hand side. The downward and upward empty black triangles are the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The thick orange line in each panel is for the model total X-ray luminosity function given by our mass-independent model, the blue dashed line is for that of the star-forming component and the red triple dashed line shows the quiescent component. The blue empty stars and the red empty circles at z < 0.90 are for the star-forming and the quiescent component of the X-ray luminosity functions, respectively, as reported in Georgakakis et al. (2014) out to z ∼ 1. We hatched the two panels at z = 2.25 and z = 3.0 to indicate that our mass-independent λEdd distribution model is unable to reproduce the observed X-ray luminosity functions at these redshifts (see Section 3.1.1). Figure 2. View largeDownload slide Fit to the X-ray 2–10 keV luminosity functions of Aird et al. (2015) when optimizing our mass-independent λEdd distribution model. Each panel represents a different redshift bin out to z = 3, and as indicated at their bottom right-hand side. The downward and upward empty black triangles are the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The thick orange line in each panel is for the model total X-ray luminosity function given by our mass-independent model, the blue dashed line is for that of the star-forming component and the red triple dashed line shows the quiescent component. The blue empty stars and the red empty circles at z < 0.90 are for the star-forming and the quiescent component of the X-ray luminosity functions, respectively, as reported in Georgakakis et al. (2014) out to z ∼ 1. We hatched the two panels at z = 2.25 and z = 3.0 to indicate that our mass-independent λEdd distribution model is unable to reproduce the observed X-ray luminosity functions at these redshifts (see Section 3.1.1). Beyond z ∼ 2, our model cannot reproduce the observed X-ray luminosity functions. At the position of the break, our model X-ray luminosity functions at these redshifts are at least a factor of three below the observed ones. We also find that the contribution from quiescent galaxies is very low (by at least a factor of ten) compared to that of the star-forming galaxies. To explain the reasons for this, we examine the stellar mass functions of Davidzon et al. (2017; see Fig. 1). They report that the contribution from quiescent galaxies to the mass function decreases as the redshift increases (the knee of the mass function for quiescent galaxies is at least a factor of ten in normalization below that of the star-forming galaxies at z ≳ 2; see Fig. 1). In our model, this lower contribution from the quiescent galaxies can be compensated by an increase in the normalization of the model Eddington ratio distribution for quiescent galaxies (as our model relates the mass function to the Eddington ratio distribution). This would increase the contribution from quiescent galaxies to the total model X-ray luminosity functions in order to match the observed ones. However, in our model, the normalization of the quiescent galaxy component cannot exceed that of the star-forming component (see assumptions in Section 2.1). While looking at the normalizations of both the star-forming (ASF) and the quiescent (AQui) component in our model, we find that beyond z ∼ 2, they are similar (i.e. ASF ∼ AQui). As a consequence, the contribution from quiescent galaxies to the total X-ray luminosity functions cannot increase enough to match the observed X-ray luminosity functions. If we relax this assumption, our model finds that the Eddington ratio distribution beyond z ∼ 2 is fully dominated by quiescent galaxies, with the normalization of the Eddington ratio distribution for star-forming galaxies to that of quiescent galaxies as low as ∼0.01 at z = 2.25 and ∼0.001 at z = 2.5. This is at odds with empirical measurements of the Eddington ratio distributions of, for example, Georgakakis et al. (2014) and Wang et al. (2017), that suggest a more equal splitting between quiescent and star-forming galaxies, as well as infrared studies that find large numbers of high-redshift AGNs residing in star-forming galaxies (e.g. Mullaney et al. 2015; Stanley et al. 2015; Netzer et al. 2016). Finally, as a further check, we include in Fig. 2 the measured X-ray luminosity functions separated into star-forming and quiescent galaxies for z < 1 as reported in Georgakakis et al. (2014). Despite not including this information during our optimization, we find good agreement between our model and these observed X-ray luminosity functions of star-forming and quiescent galaxies, increasing our confidence in the model at these lower redshifts. 3.1.2 The Eddington ratio distribution We show in Section 3.1.1 that our mass-independent model is able to reproduce the observed X-ray luminosity function out to z ∼ 2. We now present in Fig. 3 the Eddington ratio distributions at various redshifts after optimization of our mass-independent model (normalized to coincide at λEdd= 10−3). The shape of our total Eddington ratio distribution (see left-hand panel in Fig. 3) is consistent with a broken power-law, as requested by our model and in agreement with previous studies (e.g. Aird et al. 2013; Jones et al. 2017). However, we find a significant difference in the slope below the break (γ1 in equation 4) of the Eddington ratio distribution for star-forming hosts compared to that of quiescent hosts (see central and right-hand panels in Fig. 3 for the Eddington ratio distribution of star-forming and quiescent host galaxies, respectively). We report that below z ∼ 2, the slope at low-λEdd for the star-forming component is rising (although we are only able to place upper limits). As a consequence, our model Eddington ratio distribution of AGNs in star-forming galaxies is consistent with a peaky3 distribution, similar to the light-bulb shaped distribution sometimes explored in past studies (e.g. Veale et al. 2014; Stanley et al. 2015). This suggests that AGNs hosted in star-forming galaxies have typical λEdd higher than AGNs in quiescent hosts, yet the λEdd of AGNs in quiescent galaxies span a broader range of values. Finally, we notice an overall shift of the knee of the Eddington ratio distribution to higher λEdd values as redshift increases. This implies a typical AGN activity that increases as redshift increases. At z ≳ 2, we find that the shape of the Eddington ratio distribution of AGNs in star-forming galaxies becomes consistent with a broken power-law, in contrast to the peaky Eddington ratio distribution reported for lower redshifts. However, as demonstrated in Section 3.1.1 our mass-independent model is unable to reproduce the X-ray luminosity functions at these higher redshifts. We report in Table 1, the redshift evolution of the parameters that define these Eddington ratio distributions. Figure 3. View largeDownload slide The total mass-independent Eddington ratio distribution (left-hand panel), its star-forming component (central panel), and its quiescent component (right-hand panel) evaluated at z = 0.2, z = 0.5, z = 1.0, z = 1.5, z = 2.0, and z = 2.5 (colour coded; see key). The downward arrows indicate the overall effects of the upper and lower limits found in the parameters that describes these distributions. For the star-forming and quiescent panels, we also indicate with faint grey lines the corresponding total Eddington ratio distribution. The normalization is such that the total Eddington ratio distribution at z = 0.2 integrates to unity, applying an arbitrary cut at λEdd = 10−7. For ease of comparison, we then re-normalize all the distributions at higher redshifts such that their values at log(λEdd)= −3.0 coincide with that at z = 0.2, hence the arbitrary unit of the ordinate axis. Figure 3. View largeDownload slide The total mass-independent Eddington ratio distribution (left-hand panel), its star-forming component (central panel), and its quiescent component (right-hand panel) evaluated at z = 0.2, z = 0.5, z = 1.0, z = 1.5, z = 2.0, and z = 2.5 (colour coded; see key). The downward arrows indicate the overall effects of the upper and lower limits found in the parameters that describes these distributions. For the star-forming and quiescent panels, we also indicate with faint grey lines the corresponding total Eddington ratio distribution. The normalization is such that the total Eddington ratio distribution at z = 0.2 integrates to unity, applying an arbitrary cut at λEdd = 10−7. For ease of comparison, we then re-normalize all the distributions at higher redshifts such that their values at log(λEdd)= −3.0 coincide with that at z = 0.2, hence the arbitrary unit of the ordinate axis. Table 1. Redshift evolution of the parameters that describe the Eddington ratio distribution for star-forming and quiescent galaxies, given by equation (4). SF and Qui labels stand for star-forming and quiescent, respectively. The slopes and intercepts are given for an evolution in (1+z). Uncertainties on the parameters are the standard deviation of 1000 Monte Carlo realisations (i.e. 1σ). Parameters  Intercepts  Slopes    log(ASF)  −2.07±0.19  0.13±0.07  for z < 1.7      −1.15±0.03  for z > 1.7  log(AQui)  −2.80±0.23  0.37±0.10  for z < 1.7      −1.20±0.11  for z > 1.7  log$$(\lambda _{\rm break}^{\rm SF})$$  −2.36±0.06  0.65±0.02    log$$(\lambda _{\rm break}^{\rm Qui})$$  −2.62±0.16  0.66±0.06    $$\gamma _1^{\rm SF}$$  ≲ − 3.0  0.0  for z < 1.7    −0.96±0.51  0.25±0.12  for z > 1.7  $$\gamma _1^{\rm Qui}$$  0.28±0.05  0.07±0.02    $$\gamma _2^{\rm SF}$$  2.29±0.04  0.01±0.02  for z > 1.7      0.15±0.03  for z > 1.7  $$\gamma _2^{\rm Qui}$$  ≳2.7  0.0  for z < 1.7    3.75±2.33  −0.53±0.63  for z > 1.7  Parameters  Intercepts  Slopes    log(ASF)  −2.07±0.19  0.13±0.07  for z < 1.7      −1.15±0.03  for z > 1.7  log(AQui)  −2.80±0.23  0.37±0.10  for z < 1.7      −1.20±0.11  for z > 1.7  log$$(\lambda _{\rm break}^{\rm SF})$$  −2.36±0.06  0.65±0.02    log$$(\lambda _{\rm break}^{\rm Qui})$$  −2.62±0.16  0.66±0.06    $$\gamma _1^{\rm SF}$$  ≲ − 3.0  0.0  for z < 1.7    −0.96±0.51  0.25±0.12  for z > 1.7  $$\gamma _1^{\rm Qui}$$  0.28±0.05  0.07±0.02    $$\gamma _2^{\rm SF}$$  2.29±0.04  0.01±0.02  for z > 1.7      0.15±0.03  for z > 1.7  $$\gamma _2^{\rm Qui}$$  ≳2.7  0.0  for z < 1.7    3.75±2.33  −0.53±0.63  for z > 1.7  Notes: Slopes and intercepts are given for an evolution as a function of (1+z). The intercept for z > 1.7, when assuming a break in the z evolution of the parameter, is given by the continuity at z = 1.7 (i.e. [intercept forz > 1.7] = (1 + 1.7) × ([slope forz < 1.7] − [slope forz > 1.7]) + [intercept forz < 1.7]). View Large 3.1.3 SFR in bins of X-ray luminosities As previously highlighted, a key test of any model that describes the Eddington ratio distribution is whether it can reproduce other observed features of the AGN population besides just the (X-ray) luminosity functions. Here, we test our PSM (which incorporates pre-defined sSFR distributions and our own optimized mass-independent Eddington ratio distribution; see Section 2.3) by assessing whether it can reproduce the observed flat relationship between SFR and X-ray luminosity as reported in many observational studies (e.g. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). We do this by calculating the mean-average SFRs of the AGNs in our PSM in bins of 0.5 dex in X-ray luminosity, taking into account both star-forming and quiescent galaxies (see Section 2.3). As shown in Fig. 4, our mass-independent model predicts a strong correlation between these SFRs and X-ray luminosities (with an average slope of ∼0.6), at least up to z ∼ 2. This is in good agreement with empirical results for our lowest redshift bin (i.e. z = 0.35), but strongly contradicts the flat observed relationship between SFR and X-ray luminosity at higher redshifts (e.g. Rosario et al. 2012; Stanley et al. 2015). It is important to note that other studies also find a similar strongly increasing relationship between SFR and X-ray luminosity when employing a peaky (i.e. similar in shape than our peaky distribution; see central panel in Fig. 3) Eddington ratio distribution for star-forming galaxies (e.g. Veale et al. 2014; Stanley et al. 2015). We will show in Section 4 that this strong correlation is due to the AGN host mass distribution and the SFR/stellar-mass relationship (i.e. MS) induced by the peaky Eddington ratio distribution found for our mass-independent model. Figure 4. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosity at z = 0.35, z = 0.65, z = 1.15, and z = 2. The black thick lines shows the correlation predicted by our PSM that incorporates AGNs using our mass-independent Eddington ratio distribution split between star-forming and quiescent galaxies. The error bars show the 1σ uncertainties on the mean for average SFR. The thick black dashed line is identical than the thick black line but for galaxies with M* ≳ 109.5 M⊙. Various symbols are for a compilation of published SFR/X-ray relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line at z = 1.15 shows the relationship reported in Silverman et al. (2009). Finally, the black dashed line represents the correlation displayed by AGN-dominated systems in Netzer (2009). Our mass-independent model predicts a significant correlation between SFR and X-ray luminosity, which is at odds with empirical results for z ≳ 0.5. Figure 4. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosity at z = 0.35, z = 0.65, z = 1.15, and z = 2. The black thick lines shows the correlation predicted by our PSM that incorporates AGNs using our mass-independent Eddington ratio distribution split between star-forming and quiescent galaxies. The error bars show the 1σ uncertainties on the mean for average SFR. The thick black dashed line is identical than the thick black line but for galaxies with M* ≳ 109.5 M⊙. Various symbols are for a compilation of published SFR/X-ray relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line at z = 1.15 shows the relationship reported in Silverman et al. (2009). Finally, the black dashed line represents the correlation displayed by AGN-dominated systems in Netzer (2009). Our mass-independent model predicts a significant correlation between SFR and X-ray luminosity, which is at odds with empirical results for z ≳ 0.5. As SFR is correlated to stellar mass via the MS, it is possible that the empirical flat relationship between SFR and X-ray luminosity is a consequence of an observational bias towards higher mass galaxies. To test this, we recalculate from our model the average SFR in bins of X-ray luminosity, but now excluding galaxies with M* <109.5 M⊙. This mass limit roughly corresponds to that reached in the highest redshift bins of Stanley et al. (2015). We show in Fig. 4 that although this mass cut causes the model SFR/X-ray relationship to rise slightly at z = 2.0 and LX ≲ 1043.5 erg s−1, it is not sufficient to reproduce the observed flat relationship between SFR and X-ray luminosity. 3.2 The mass-dependent model Since our mass-independent model fails to reproduce the X-ray luminosity functions at z ≳ 2 and is unable to reproduce the flat relationship between SFR and X-ray luminosity (at z > 0.5), we consider whether relaxing the mass-independence requirement (as suggested by recent studies, e.g. Bongiorno et al. 2016; Georgakakis et al. 2017; Aird et al. 2018) can help to resolve these issues. As explained in Section 2.2, for this model, we assume that the quiescent component remains unchanged from the mass-independent model meaning that we only fit the star-forming X-ray luminosity functions at z < 2. These are extracted from the results of our mass-independent model, and is justified by the ability of our mass-independent model to reproduce the X-ray luminosity functions of star-forming and quiescent galaxies to z ∼ 1 from Georgakakis et al. (2014). 3.2.1 X-ray luminosity functions We show in Fig. 5, the fit to the star-forming component of the X-ray luminosity functions using our mass-dependent Eddington ratio distribution model to z = 1.75. Again, we do not extend our mass-dependent model beyond z ∼ 2 since it relies on the ability of our previous – mass-independent – model to split the X-ray luminosity functions in terms of star-forming and quiescent components (see Section 2.2). Out to z ∼ 2, our mass-dependent model does a good job at fitting the star-forming X-ray luminosity functions in all our redshift bins. Since our mass-dependent model assumes an λEdd distribution split into three different mass bins (see Section 2.2), we are also able to derive the X-ray luminosity function of each mass bin (see Fig. 5). In doing so, our model predicts that the highest mass galaxies (i.e. M* >1011 M⊙) dominate the star-forming X-ray luminosity functions across almost the full range of X-ray luminosities (i.e. 1041 <LX <1046 erg s−1). The exception being in our lowest redshift bin (i.e. z = 0.11), where 1010 <M*/M⊙ <1011 galaxies contribute marginally more than our highest mass galaxies at LX ≳ 1043.5 erg s−1. In general, however, these medium mass galaxies only contribute significantly (i.e. >10 per cent) at luminosities above the knee of the X-ray luminosity functions (i.e. LX >1044 erg s−1). By contrast, the lowest mass galaxies (i.e. M* <1010 M⊙) provide almost no contribution to the X-ray luminosity functions of star-forming galaxies at LX ≳ 1042.5 erg s−1 (i.e. <10 per cent), although we are only able to place upper limits on this contribution from lower mass galaxies. Figure 5. View largeDownload slide Fit to the star-forming component of the X-ray 2–10 keV luminosity functions using our mass-dependent Eddington ratio distribution model up to z = 1.75. Each panel shows a different redshift bin as indicated on their bottom right-hand side. The black line in each panel shows the star-forming component of the X-ray luminosity function (derived from our mass-independent model), and the orange line shows the best fit given by our mass-dependent model. The dashed, single-dot dashed, and triple-dot dashed lines show the X-ray luminosity functions for our highest (11 < log(M*/M⊙) < 12), medium (10 < log(M*/M⊙) < 11), and lowest (8 < log(M*/M⊙) < 10) mass bin, respectively. The black arrows indicate the overall effects of upper and lower limits on the parameters that define our mass-dependent Eddington ratio distribution. The blue stars and the red circles show the X-ray luminosity functions for the star-forming and quiescent galaxies, respectively, as reported in Georgakakis et al. (2014). The light grey downward and upward triangles show the total X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. Figure 5. View largeDownload slide Fit to the star-forming component of the X-ray 2–10 keV luminosity functions using our mass-dependent Eddington ratio distribution model up to z = 1.75. Each panel shows a different redshift bin as indicated on their bottom right-hand side. The black line in each panel shows the star-forming component of the X-ray luminosity function (derived from our mass-independent model), and the orange line shows the best fit given by our mass-dependent model. The dashed, single-dot dashed, and triple-dot dashed lines show the X-ray luminosity functions for our highest (11 < log(M*/M⊙) < 12), medium (10 < log(M*/M⊙) < 11), and lowest (8 < log(M*/M⊙) < 10) mass bin, respectively. The black arrows indicate the overall effects of upper and lower limits on the parameters that define our mass-dependent Eddington ratio distribution. The blue stars and the red circles show the X-ray luminosity functions for the star-forming and quiescent galaxies, respectively, as reported in Georgakakis et al. (2014). The light grey downward and upward triangles show the total X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. Our mass-dependent model suggests that AGNs hosted by star-forming galaxies with M* >1011 M⊙ dominate the X-ray luminosity functions. This is consistent with Georgakakis et al. (2017), which demonstrates that AGNs with LX >1041 − 42 erg s−1 significantly contribute to the AGN host stellar mass function for masses M* >1011 M⊙. This is also consistent with results from Bongiorno et al. (2016), which show that the host stellar mass function of AGNs with LX >1043 erg s−1 peaks at M* ∼1011 M⊙, out to z = 2. However, the vast majority of AGN host galaxies with LX <1043 erg s−1 display typical stellar masses in the range M* ∼1010 − 11 M⊙ (e.g. Bongiorno et al. 2016; Aird, Coil & Georgakakis 2017), which contrasts with our findings of a model X-ray luminosity function dominated by host galaxies with M* >1011 M⊙. However, our model is limited by the choice for the boundaries of our mass bins. Therefore, should we change our highest mass bin to incorporate galaxies with M* ∼1010.5 M⊙, we would expect to find that galaxies with M* >1010.5 M⊙ dominate the X-ray luminosity functions, in agreement with recent observational studies (e.g. Bongiorno et al. 2016; Aird et al. 2017; Wang et al. 2017). 3.2.2 SFR in bins of X-ray luminosities Having confirmed that our mass-dependent model is able to reproduce the X-ray luminosity functions for star-forming galaxies up to z ∼ 2, we now consider whether the corresponding PSM reproduces the observed flat relationship between SFR and X-ray luminosity (e.g. Rosario et al. 2012; Stanley et al. 2015). We show in Fig. 6, the mean-average SFR of AGNs split into 0.5 dex-wide bins of X-ray luminosity at similar redshifts as those of Stanley et al. (2015). Again, we include quiescent galaxies in our model when calculating these averages. Contrary to our mass-independent model, we find that our mass-dependent model predicts a flat relationship between SFR and X-ray luminosity in very good agreement with observations (see Fig. 6). It is important to stress that the model was not optimized to recreate the flat SFR/X-ray luminosity relationship. Therefore, we conclude that the mass-dependent model is able to reproduce the X-ray luminosity functions for star-forming galaxies (with a good agreement with observations at least up to z ∼ 1), while also independently reproducing the observed flat relationship between SFR and X-ray luminosity out to z ∼ 2. Figure 6. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosities predicted by our mass-dependent model (black thick line) at z = 0.35, z = 0.65, z = 1.15, and z = 2.0. The error bars represent the 1σ uncertainties on the mean for average SFR. The various symbols are from a compilation of published relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line shows the relationship between average SFR and X-ray luminosity found in Silverman et al. (2009), and the black dashed line represents the relationship for AGN-dominated systems, as reported in Netzer (2009). Although it is not implemented in the model, our mass-independent model successfully reproduces the observed flat SFR/X-ray relationship out to z = 2. Figure 6. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosities predicted by our mass-dependent model (black thick line) at z = 0.35, z = 0.65, z = 1.15, and z = 2.0. The error bars represent the 1σ uncertainties on the mean for average SFR. The various symbols are from a compilation of published relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line shows the relationship between average SFR and X-ray luminosity found in Silverman et al. (2009), and the black dashed line represents the relationship for AGN-dominated systems, as reported in Netzer (2009). Although it is not implemented in the model, our mass-independent model successfully reproduces the observed flat SFR/X-ray relationship out to z = 2. 3.2.3 The evolution of the λEdd distribution Having demonstrated that our mass-dependent model is able to reproduce both the X-ray luminosity functions for star-forming galaxies and the flat SFR/X-ray luminosity relationship to z ∼ 2, we now explore how the eight parameters (see Section 2.2) that define our mass-dependent model change with redshift. We show in Fig. 7, the redshift evolution of these parameters. We report that the overall normalization of the star-forming component (ASF) increases very slightly with increasing redshift. The position of the break ($$\lambda _{\rm break}^{\rm SF}$$) in each mass bin shifts towards higher λEdd values. This is consistent with the overall increase of the position of the break observed in our mass-independent model (i.e. all mass bins collapsed together; see Section 3.1.2), suggesting a higher typical accretion rate for AGNs at higher redshifts. Regarding the slopes at lower λEdd, ($$\gamma _1^{\rm SF}$$), we find that it does not evolve with redshift for any of our mass bins. However, although $$\gamma _1^{\rm SF}$$ for our highest mass bin is well constrained, we have large uncertainties associated with $$\gamma _1^{\rm SF}$$ for our medium mass bin, and upper limits for our lowest mass bin (see Fig. 7). These weak constraints on $$\gamma _1^{\rm SF}$$ for our medium and lowest mass bins are a direct consequence of their low contributions to the total Eddington ratio distribution for star-forming galaxies (see Fig. 8), as was hinted by their low contribution to the total X-ray luminosity functions (see Section 3.2). We also find that the shared slope at high λEdd ($$\gamma _2^{\rm SF}$$) does not change significantly with redshift. Overall, our findings suggest that there is a suppression of low-λEdd (i.e. λEdd ≲0.1) for AGNs hosted in lower mass galaxies (i.e. M* ≲ 1010 − 11 M⊙). Finally, we performed a linear fit of the redshift evolution of each parameter, using the standard deviation of 1000 Monte Carlo realisations for the uncertainties on the fitting parameters. We report the results of this fit in Table 2. Figure 7. View largeDownload slide Redshift evolution of the parameters that define our mass-dependent Eddington ratio distribution for star-forming galaxies at z = 0.1, z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.3, z = 1.5, and z = 1.7. Error bars show the 3σ uncertainties. The upward and downward arrows indicate the lower and upper limits, respectively. The two top panels present the shared parameters among our three mass bins, i.e. the normalization, $$\log _{A_{\rm SF}}$$, and the slope at high Eddington ratio, $$\gamma _2^{\rm SF}$$. The left-hand panels are the three different break positions $$\lambda _{\rm break}^{\rm SF}$$, one for each mass bin (as indicated in the bottom left-hand side in each panel), and the right-hand panels are the same but for that of the lower Eddington ratio slope, $$\gamma _1^{\rm SF}$$. The black lines show the best fit with redshift of each parameter. Figure 7. View largeDownload slide Redshift evolution of the parameters that define our mass-dependent Eddington ratio distribution for star-forming galaxies at z = 0.1, z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.3, z = 1.5, and z = 1.7. Error bars show the 3σ uncertainties. The upward and downward arrows indicate the lower and upper limits, respectively. The two top panels present the shared parameters among our three mass bins, i.e. the normalization, $$\log _{A_{\rm SF}}$$, and the slope at high Eddington ratio, $$\gamma _2^{\rm SF}$$. The left-hand panels are the three different break positions $$\lambda _{\rm break}^{\rm SF}$$, one for each mass bin (as indicated in the bottom left-hand side in each panel), and the right-hand panels are the same but for that of the lower Eddington ratio slope, $$\gamma _1^{\rm SF}$$. The black lines show the best fit with redshift of each parameter. Figure 8. View largeDownload slide The Eddington ratio probability distributions of star-forming and quiescent galaxies in our mass-dependent model at z = 0.2, z = 0.5, z = 1.0, z = 1.5, and z = 2.0. The left-hand panel shows both the star-forming (plain lines) and the quiescent components (single-dashed lines) of the Eddington ratio probability distribution. The normalization is such that the total Eddington ratio distribution at each redshift integrates to unity after applying an arbitrary cut at λEdd= 10−7. The right-hand panels are, for each redshift (indicated on their top right-hand corner), the contribution of the different mass bins to the Eddington ratio distributions of AGNs hosted in star-forming galaxies. In each of these panels, the dashed line is for the highest mass bin (11 < log(M*/M⊙) < 12), the single-dotted dashed line is for the medium mass bin (10 < log(M*/M⊙) < 11), and the triple-dotted dashed line is for the lowest mass bin (8 < log(M*/M⊙) < 10). The downward arrows for the slope at low-λEdd in our lower mass bin indicate that we have upper limits. The grey area illustrates the large uncertainties for the slope at low-λEdd in our medium mass bin (see Fig. 7). We normalized each of the total Eddington ratio distribution such that it integrates to unity, applying an arbitrary cut at λEdd = 10−7. Figure 8. View largeDownload slide The Eddington ratio probability distributions of star-forming and quiescent galaxies in our mass-dependent model at z = 0.2, z = 0.5, z = 1.0, z = 1.5, and z = 2.0. The left-hand panel shows both the star-forming (plain lines) and the quiescent components (single-dashed lines) of the Eddington ratio probability distribution. The normalization is such that the total Eddington ratio distribution at each redshift integrates to unity after applying an arbitrary cut at λEdd= 10−7. The right-hand panels are, for each redshift (indicated on their top right-hand corner), the contribution of the different mass bins to the Eddington ratio distributions of AGNs hosted in star-forming galaxies. In each of these panels, the dashed line is for the highest mass bin (11 < log(M*/M⊙) < 12), the single-dotted dashed line is for the medium mass bin (10 < log(M*/M⊙) < 11), and the triple-dotted dashed line is for the lowest mass bin (8 < log(M*/M⊙) < 10). The downward arrows for the slope at low-λEdd in our lower mass bin indicate that we have upper limits. The grey area illustrates the large uncertainties for the slope at low-λEdd in our medium mass bin (see Fig. 7). We normalized each of the total Eddington ratio distribution such that it integrates to unity, applying an arbitrary cut at λEdd = 10−7. Table 2. Redshift evolution of the parameters that describe the Eddington ratio distribution of star-forming galaxies for our mass-dependent model. The slopes and intercepts are given for an evolution as a function of (1+z). Uncertainties on the fitting parameters are extracted via 1000 Monte Carlo realisations (i.e. 1σ). Parameters  Intercepts  Slopes  log(ASF)  −1.57±0.08  0.03±0.04  log($$\lambda _{\rm break}^{\rm low\ mass}$$)  >−2.53  1.27  log($$\lambda _{\rm break}^{\rm medium\ mass}$$)  −2.05±0.15  0.68±0.08  log($$\lambda _{\rm break}^{\rm high\ mass}$$)  −2.60±0.08  0.70±0.04  $$\gamma _1^{\rm low\ mass}$$  <0.5  0.0  $$\gamma _1^{\rm medium\ mass}$$  −1.0±2.67  −0.12±1.30  $$\gamma _1^{\rm high\ mass}$$  0.45±0.05  −0.04±0.02  γ2  2.21±0.06  −0.03±0.03  Parameters  Intercepts  Slopes  log(ASF)  −1.57±0.08  0.03±0.04  log($$\lambda _{\rm break}^{\rm low\ mass}$$)  >−2.53  1.27  log($$\lambda _{\rm break}^{\rm medium\ mass}$$)  −2.05±0.15  0.68±0.08  log($$\lambda _{\rm break}^{\rm high\ mass}$$)  −2.60±0.08  0.70±0.04  $$\gamma _1^{\rm low\ mass}$$  <0.5  0.0  $$\gamma _1^{\rm medium\ mass}$$  −1.0±2.67  −0.12±1.30  $$\gamma _1^{\rm high\ mass}$$  0.45±0.05  −0.04±0.02  γ2  2.21±0.06  −0.03±0.03  Notes: Slopes and intercepts are given for an evolution as a function of (1+z). The intercept for z > 1.7 when assuming a break in the z evolution of the parameters is given by the continuity at z = 1.7 (i.e. [intercept forz > 1.7] = (1 + 1.7) × ([slope forz < 1.7] − [slope forz > 1.7]) + [intercept forz < 1.7]). View Large Each of the aforementioned trends can be seen in the evolution of the overall Eddington ratio distribution, which we plot out to z = 2 in Fig. 8. In this figure, we also plot the Eddington ratio distributions of quiescent host galaxies found in our mass-independent model at similar redshifts. In Fig. 8, for the left-hand panel, each distribution is normalize such that the total Eddington ratio distribution (i.e. the combination of the star-forming and quiescent component) integrates to unity after applying an arbitrary cut at λEdd= 10−7. For the right-hand panels, each total distribution is also normalized such that it integrates to unity after applying an arbitrary cut at λEdd= 10−7. 4 DISCUSSION 4.1 Further checks to our mass-dependent model As a further check of our mass-dependent model, we compared the model Eddington ratio distribution at z = 1 to the empirical one reported in Wang et al. (2017). As our various assumptions for this model are based on observations from Aird et al. (2018), comparing our mass-dependent Eddington ratio distribution to that empirical of Wang et al. (2017) – instead of Aird et al. (2018) – constitutes a more independent test for our model. We show this comparison in Fig. 9. We find very good agreement between our model Eddington ratio distributions to that of Wang et al. (2017) at z = 1 for both star-forming and quiescent component (at least for log(λEdd) < −0.5). This strongly supports our mass-dependent model for the Eddington ratio distribution of star-forming galaxies (since the mass-independent model was predicting a peaky Eddington ratio distribution for star-forming galaxies, which would contrasts with empirical results from Wang et al. 2017). We predict that the low-λEdd end of the Eddington ratio distribution for star-forming galaxies is dominated by galaxies with M* ≳ 1011 M⊙. However, Wang et al. (2017) reported that their sample of AGNs have typical stellar masses between 1010.5 <M*/M⊙ <1011 (once corrected for the differences in the IMFs between Salpeter 1955 for this work and Chabrier 2003 for Wang et al. 2017), which is half a dex below our highest mass bin. As mentioned in Section 3.2.1, this discrepancy in the host stellar masses between our model and observations is a consequence of the chosen boundaries for each of our mass bins. Figure 9. View largeDownload slide Comparison of our mass-dependent model Eddington ratio distribution to empirical results of Wang et al. (2017) at z = 1. The blue thick line is the total model Eddington ratio distribution for star-forming galaxies, while the dashed line, the dotted-dashed line, and the triple-dotted dashed line indicate our highest, medium, and lowest mass bin contributions, respectively. The blue circles are empirical data of the star-forming component from Wang et al. (2017) at z ∼ 1, while the red circles are that of the quiescent component. Figure 9. View largeDownload slide Comparison of our mass-dependent model Eddington ratio distribution to empirical results of Wang et al. (2017) at z = 1. The blue thick line is the total model Eddington ratio distribution for star-forming galaxies, while the dashed line, the dotted-dashed line, and the triple-dotted dashed line indicate our highest, medium, and lowest mass bin contributions, respectively. The blue circles are empirical data of the star-forming component from Wang et al. (2017) at z ∼ 1, while the red circles are that of the quiescent component. As a final check of our mass-dependent model, we show in Fig. 10 how the normalized average SFR (i.e. SFR relative to that of the MS) changes with λEdd in our mass-dependent model. Although it is not incorporated in the optimization, we predict a slight enhancement of average normalize SFR at higher λEdd (i.e. λEdd ≳ 1) compared to that of their lower λEdd counterpart (at least at z ≳ 1.2). This slight enhancement of normalized average SFR at higher λEdd is also observed in Bernhard et al. (2016). In Fig. 10, we only consider star-forming galaxies for our model since we do not have a good prescription for the normalized SFRs of quiescent galaxies. This is a possible reason for the discrepancy at lower λEdd between our model and the empirical results (i.e. that does contain quiescent galaxies). We also stress that in Fig. 10 our highest redshift bin, i.e. 1.8 < z < 2.9 constitutes an extrapolation of our model. However, the results are still consistent with observations (i.e. a slight enhancement of the normalized average SFR at higher λEdd). Figure 10. View largeDownload slide The normalized average SFR versus the Eddington ratio up to z ∼ 3. The blue lines show the prediction for our mass dependent model. the orange circles show the results report in Bernhard et al. (2016). The grey area represent the scatter around the MS as reported in Schreiber et al. (2015), the dashed black line being the MS. Figure 10. View largeDownload slide The normalized average SFR versus the Eddington ratio up to z ∼ 3. The blue lines show the prediction for our mass dependent model. the orange circles show the results report in Bernhard et al. (2016). The grey area represent the scatter around the MS as reported in Schreiber et al. (2015), the dashed black line being the MS. 4.2 Why our mass-independent model predicts a strong SFR/X-ray relationship One important result of our mass-independent model is that the Eddington ratio distribution for star-forming galaxies is better represented by a peaky distribution (i.e. a narrow distribution), rather than the more common broken power-law (i.e. a broader distribution). To explain this, we refer to Caplar et al. (2015; see also Weigel et al. 2017), where they demonstrate that the steepness of the faint-end slope of the X-ray luminosity function is determined by either the low-mass end of the mass function or the low-λEdd slope of the Eddington ratio distribution (since combined together to form the X-ray luminosity function) whichever is the steeper. Since that the faint-end of the observed X-ray luminosity function of Aird et al. (2015) flattens at higher redshifts whilst, in contrast, the low-mass end of the mass function for star-forming galaxies steepens with redshift (see Fig. 1), our mass-independent model attempts to reduce the contribution from star-forming galaxies at lower λEdd, leading the peaky Eddington ratio distribution for star-forming galaxies found by our mass-independent model. However, beyond z ∼ 2, the low-mass end of the mass function for star-forming galaxies of Davidzon et al. (2017) is already too steep compared to the flat faint-end of the X-ray luminosity function, such that our model is unable to reproduce the observed X-ray luminosity function. The steepening of the mass functions combined with the flattening of the X-ray luminosity functions suggests a mass-dependent Eddington ratio distribution for star-forming galaxies. We further saw that our mass-independent model failed to reproduce the observed flat SFR/X-ray relationship (see Fig. 4). To explain why this happens, we show in Fig. 11, which galaxy masses populate different parts of the SFR/X-ray plane. This plot shows that the mass-independent model produces a strong correlation between X-ray luminosity and stellar mass. This is a direct result of the narrower peaky (in contrast to a broader broken power law) Eddington ratio distribution for star-forming galaxies that this model demands in order to reproduce the observed X-ray luminosity functions, at least out to z ∼ 2. The narrowness of this Eddington ratio distribution means that a galaxy of a given mass can produce an AGN only within a limited range of luminosities. When we then include the correlation between SFR and stellar mass for MS galaxies (e.g. Schreiber et al. 2015), the consequence is a correlation between SFR and X-ray luminosity. This correlation between SFR and X-ray luminosity was also noted by Veale et al. (2014) when using their light-bulb Eddington ratio model, which is similar in shape to our peaky distribution for star-forming host galaxies. Figure 11. View largeDownload slide Stellar mass distributions (colour coded) in the SFR/X-ray luminosity parameter space, using our PSM and assuming our mass-independent model for the Eddington ratio distribution. Each panel shows a different redshift, as indicated on their bottom-right hand side. Each discrete colour indicates a stellar mass density in bins of log(M*)=0.5 (see colour-bar). Figure 11. View largeDownload slide Stellar mass distributions (colour coded) in the SFR/X-ray luminosity parameter space, using our PSM and assuming our mass-independent model for the Eddington ratio distribution. Each panel shows a different redshift, as indicated on their bottom-right hand side. Each discrete colour indicates a stellar mass density in bins of log(M*)=0.5 (see colour-bar). Overall, there is a large discrepancy with our mass-independent model that demands a narrow distribution for the Eddington ratio distribution of star-forming galaxies in order to reproduce the X-ray luminosity functions, but also requires a broader Eddington ratio distribution to reproduce the flat SFR/X-ray relationship. This strongly suggests that the Eddington ratio distribution must be somehow mass-dependent, such that low Eddington ratios are suppressed in low-mass galaxies to be able to reproduce the X-ray luminosity functions while still being broad to reproduce the flat relationship between SFR and X-ray luminosity. Our mass-dependent model provides further constrains on how the Eddington ratio distribution for star-forming galaxies should be mass-dependent. 4.3 Caveats to our models While deriving the λEdd distribution, our model is limited by the functional form chosen for each component of the λEdd distribution and the associated assumptions used to avoid degeneracies during fitting (see Section 2). An obvious caveat to our model is that we cannot explore the full range of possible Eddington ratio distributions. However, we find that, should we relax the assumptions made prior to the fit of the X-ray luminosity functions, our solutions are too degenerate to form a set of meaningful parameters for the Eddington ratio distributions, meaning that the X-ray luminosity functions lack the information needed to constrain the Eddington ratio distribution and its mass-dependency. Using these assumptions is, therefore, a requirement in our modelling approach. As stressed in Section 2, all our assumptions are motivated by the findings of recent studies (e.g. Jones et al. 2017; Aird et al. 2018), and the results predicted by the mass-dependent model are in good agreement with many recent observational studies, i.e. the observed X-ray luminosity functions out to z ∼ 2 – including when split in terms of star-forming and quiescent galaxies (e.g. Georgakakis et al. 2014; Aird et al. 2015), the flat relationship between SFR and X-ray luminosity (e.g. Rosario et al. 2012; Stanley et al. 2015), the empirical Eddington ratio distribution split between star-forming and quiescent component at z = 1 (Wang et al. 2017), and the enhancement of the star-forming properties at higher λEdd (Bernhard et al. 2016). As such while we do not claim that our Eddington ratio distributions are universal, they provide a simple means by which to explore the AGN-galaxy connection. One of the main caveats for our mass-dependent model is that we are limited by the choice of the stellar mass bins. As a consequence, we find that our dominant component for the Eddington ratio distribution and the X-ray luminosity function is our highest mass bin with M* >1011 M⊙. This is qualitatively consistent with observations (e.g. Georgakakis et al. 2014; Aird et al. 2015; Wang et al. 2017), yet roughly half a dex above the observed typical host mass for AGNs. Investigating this in more detail would require additional mass bins. However, this would also generates a larger number of free parameters, and therefore degeneracies. We have chosen these mass bins to be broad enough to cover a wide range of masses, yet probing any mass dependency of the Eddington ratio distribution for AGNs in star-forming hosts. 4.4 Extending our mass-dependent model to higher redshifts We recall that we fit the mass-dependent model to the X-ray luminosity functions of star-forming galaxies only, which was isolated using our mass-independent model. The results of our mass-dependent model therefore depends on the reliability of our mass-independent model to separate the X-ray luminosity functions into star-forming and quiescent components (see Section 2). Thankfully, we can verify this up to z ∼ 1 using the results reported in Georgakakis et al. (2014). However, beyond z ∼ 1, we do not have the empirical results to perform this check, while at z ≳ 2 our mass-independent model is actually inconsistent with observations (i.e. too many AGNs in quiescent galaxies). Faced with this situation where we cannot reliably exploit the results of our mass-independent model, we now investigate whether we can fit the observed X-ray luminosity functions at z = 2.25 using both a mass-dependent star-forming component and a mass-independent quiescent component λEdd distributions. This results in a total of 12 free parameters to optimize (i.e. eight parameters for the star-forming component and four parameters for the quiescent one). We use the same assumptions for the λEdd distribution of star-forming galaxies as described in Section 2, and demand that the normalization of the Eddington ratio distribution for quiescent galaxies lies below that of star-forming galaxies (i.e. consistent with the lower redshift bins). We show in Fig. 12, the results of this new model. Contrary to our mass-independent model, we are now able to reproduce the total observed X-ray luminosity functions at z = 2.25 by placing the majority of AGNs in star-forming galaxies. In particular, we find that the X-ray luminosity functions at z = 2.25 is dominated by the highest mass star-forming galaxies, with a smaller contribution from medium mass star-forming galaxies. The contribution from our lowest mass bin is consistent with zero. Similarly, we can only place upper limits on the contribution from quiescent galaxies to the total X-ray luminosity functions at z = 2.25. We conclude, therefore, that while we are able to reproduce the total X-ray luminosity functions using a mass-dependent λEdd distribution for star-forming galaxies (and ensuring that it is dominated by star-forming galaxies), this solution is degenerate with the level of contribution from quiescent galaxies. To break this degeneracy will require the separation of the high-redshift X-ray luminosity functions into quiescent and star-forming components (i.e. as performed by Georgakakis et al. (2014) at z < 1). Figure 12. View largeDownload slide Fit of the X-ray luminosity functions at z = 2.25 of Aird et al. (2015) for our model assuming a mass-dependent Eddington ratio distribution for star-forming galaxies and a mass-independent one for that of quiescent galaxies. The black downward and upward triangles show the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The orange line shows the total X-ray luminosity functions derived using our model, the blue lines are for that of the star-forming component in each different mass bin (see keys), and the quiescent component is shown with a red line. Downward black arrows indicate upper limits. Figure 12. View largeDownload slide Fit of the X-ray luminosity functions at z = 2.25 of Aird et al. (2015) for our model assuming a mass-dependent Eddington ratio distribution for star-forming galaxies and a mass-independent one for that of quiescent galaxies. The black downward and upward triangles show the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The orange line shows the total X-ray luminosity functions derived using our model, the blue lines are for that of the star-forming component in each different mass bin (see keys), and the quiescent component is shown with a red line. Downward black arrows indicate upper limits. 5 CONCLUSION Motivated by recent results reporting a different Eddington ratio distribution for star-forming and quiescent galaxies (Wang et al. 2017; Aird et al. 2018), we attempt to constrain these distributions by using an analytical model to fit the observed X-ray luminosity functions of Aird et al. (2015), assuming the galaxy mass functions of Davidzon et al. (2017). Our first model assumes two mass-independent λEdd distributions: one for star-forming galaxies and another for quiescent galaxies. After optimization, we find that this model is able to reproduce the X-ray luminosity functions out to z ∼ 2 (see Section 3.1.1), but demands a peaky distribution for the Eddington ratio distribution of star-forming galaxies. Despite this, our mass-independent model fails to reproduce the observed flat SFR/X-ray relationship when incorporated into a PSM for galaxies (see Section 3.1.3). We argue that this failure arises because a mass-independent model places too many low-luminosity AGNs in low-mass galaxies (see Section 4). To resolve this problem, we develop a second model in which the Eddington ratio distribution for star-forming galaxies is allowed to be different in three mass bins (motivated by observation from e.g. Georgakakis et al. 2017; Aird et al. 2018). By suppressing low Eddington ratio AGNs in low-mass galaxies, this mass-dependent model is able to reproduce the X-ray luminosity functions for star-forming galaxies while simultaneously reproducing the observed flat relationship between SFR and X-ray luminosity (see Section 3.2.2). Finally, we also find that the mass-dependent model is consistent with empirical Eddington ratio distribution for star-forming and quiescent galaxies at z = 1, and is able to reproduce the slight enhancement of normalized average SFR at higher λEdd, as reported in Bernhard et al. (2016). Overall, we conclude that, in our model, a suppression of lower AGN accretion activity in low mass galaxies is required to reproduce both the observed X-ray luminosity functions and the observed flat SFR/X-ray relationship. Acknowledgements We thank the referee for the useful comments that help to significantly improve the quality of the paper. EB thanks C. Tadhunter and D. Alexander for the helpful discussions. EB thanks the University of Sheffield for the support grant. Footnotes 1 To a constant factor given in equation (1). 2 emcee is available on-line at http://dan.iel.fm/emcee/current/. 3 Hereafter, we refer to the shape of our model Eddington ratio distribution of AGNs hosted by star-forming galaxies as peaky. REFERENCES Aird J. et al.  , 2010, MNRAS , 401, 2531 https://doi.org/10.1111/j.1365-2966.2009.15829.x CrossRef Search ADS   Aird J. et al.  , 2013, ApJ , 775, 41 https://doi.org/10.1088/0004-637X/775/1/41 CrossRef Search ADS   Aird J., Coil A. L., Georgakakis A., Nandra K., Barro G., Pérez-González P. G., 2015, MNRAS , 451, 1892 https://doi.org/10.1093/mnras/stv1062 CrossRef Search ADS   Aird J., Coil A. L., Georgakakis A., 2017, MNRAS , 465, 3390 https://doi.org/10.1093/mnras/stw2932 CrossRef Search ADS   Aird J. Coil A. L. Georgakakis A., 2018, MNRAS , 474, 1225 CrossRef Search ADS   Assef R. J. et al.  , 2011, ApJ , 728, 56 https://doi.org/10.1088/0004-637X/728/1/56 CrossRef Search ADS   Azadi M. et al.  , 2015, ApJ , 806, 187 https://doi.org/10.1088/0004-637X/806/2/187 CrossRef Search ADS   Bernhard E., Béthermin M., Sargent M., Buat V., Mullaney J. R., Pannella M., Heinis S., Daddi E., 2014, MNRAS , 442, 509 https://doi.org/10.1093/mnras/stu896 CrossRef Search ADS   Bernhard E., Mullaney J. R., Daddi E., Ciesla L., Schreiber C., 2016, MNRAS , 460, 902 https://doi.org/10.1093/mnras/stw973 CrossRef Search ADS   Bongiorno A. et al.  , 2016, A&A , 588, A78 CrossRef Search ADS   Caplar N., Lilly S. J., Trakhtenbrot B., 2015, ApJ , 811, 148 https://doi.org/10.1088/0004-637X/811/2/148 CrossRef Search ADS   Chabrier G., 2003, PASP , 115, 763 https://doi.org/10.1086/376392 CrossRef Search ADS   Conroy C., White M., 2013, ApJ , 762, 70 https://doi.org/10.1088/0004-637X/762/2/70 CrossRef Search ADS   Daddi E. et al.  , 2007, ApJ , 670, 156 https://doi.org/10.1086/521818 CrossRef Search ADS   Davidzon I. et al.  , 2017, A&A , 605, A70 CrossRef Search ADS   Elbaz D. et al.  , 2011, A&A , 533, A119 CrossRef Search ADS   Fabian A. C., 2012, ARA&A , 50, 455 CrossRef Search ADS   Foreman-Mackey D., Hogg D. W., Lang D., Goodman J., 2013, PASP , 125, 306 https://doi.org/10.1086/670067 CrossRef Search ADS   Gelman A., Rubin D. B., 1992, Stat. sci. , pp 457– 472 Georgakakis A. et al.  , 2014, MNRAS , 443, 3327 https://doi.org/10.1093/mnras/stu1326 CrossRef Search ADS   Georgakakis A., Aird J., Schulze A., Dwelly T., Salvato M., Nandra K., Merloni A., Schneider D. P., 2017, MNRAS , 471, 1976 https://doi.org/10.1093/mnras/stx1602 CrossRef Search ADS   Goodman J., Weare J., 2010, Commun. Appl. Math. Comput. sci. , 5, 65 https://doi.org/10.2140/camcos.2010.5.65 CrossRef Search ADS   Harrison C. M. et al.  , 2012, ApJ , 760, L15 https://doi.org/10.1088/2041-8205/760/1/L15 CrossRef Search ADS   Heckman T. M., Kauffmann G., Brinchmann J., Charlot S., Tremonti C., White S. D. M., 2004, ApJ , 613, 109 https://doi.org/10.1086/422872 CrossRef Search ADS   Hickox R. C., Mullaney J. R., Alexander D. M., Chen C.-T. J., Civano F. M., Goulding A. D., Hainline K. N., 2014, ApJ , 782, 9 https://doi.org/10.1088/0004-637X/782/1/9 CrossRef Search ADS   Ilbert O. et al.  , 2013, A&A , 556, A55 CrossRef Search ADS   Jones M. L., Hickox R. C., Mutch S. J., Croton D. J., Ptak A. F., DiPompeo M. A., 2017, ApJ , 843, 125 https://doi.org/10.3847/1538-4357/aa7632 CrossRef Search ADS   Kormendy J., Richstone D., 1995, ARA&A , 33, 581 CrossRef Search ADS   Kormendy J., Bender R., Cornell M. E., 2011, Nature , 469, 374 https://doi.org/10.1038/nature09694 CrossRef Search ADS PubMed  Lanzuisi G. et al.  , 2017, A&A , 602, A123 CrossRef Search ADS   Larson D. et al.  , 2011, ApJS , 192, 16 https://doi.org/10.1088/0067-0049/192/2/16 CrossRef Search ADS   Lutz D. et al.  , 2010, ApJ , 712, 1287 https://doi.org/10.1088/0004-637X/712/2/1287 CrossRef Search ADS   Magorrian J. et al.  , 1998, AJ , 115, 2285 https://doi.org/10.1086/300353 CrossRef Search ADS   Marconi A., Hunt L. K., 2003, ApJ , 589, L21 https://doi.org/10.1086/375804 CrossRef Search ADS   Merloni A., Rudnick G., Di Matteo T., 2004, MNRAS , 354, L37 https://doi.org/10.1111/j.1365-2966.2004.08382.x CrossRef Search ADS   Merritt D., 2000, in Combes F. Mamon G. A. Charmandaris V., eds, ASPacific Conf. Ser. Vol. 197, Dynamics of Galaxies: from the Early Universe to the Present . Astron. Soc. Pac., San Francisco, p. 221 Mullaney J. R. et al.  , 2012, ApJ , 753, L30 https://doi.org/10.1088/2041-8205/753/2/L30 CrossRef Search ADS   Mullaney J. R. et al.  , 2015, MNRAS , 453, L83 https://doi.org/10.1093/mnrasl/slv110 CrossRef Search ADS   Netzer H., 2009, MNRAS , 399, 1907 https://doi.org/10.1111/j.1365-2966.2009.15434.x CrossRef Search ADS   Netzer H., Lani C., Nordon R., Trakhtenbrot B., Lira P., Shemmer O., 2016, ApJ , 819, 123 https://doi.org/10.3847/0004-637X/819/2/123 CrossRef Search ADS   Rodighiero G. et al.  , 2011, ApJ , 739, L40 https://doi.org/10.1088/2041-8205/739/2/L40 CrossRef Search ADS   Rosario D. J. et al.  , 2012, A&A , 545, A45 CrossRef Search ADS   Rosario D. J. et al.  , 2013, A&A , 560, A72 CrossRef Search ADS   Rovilos E. et al.  , 2012, A&A , 546, A58 CrossRef Search ADS   Salim S. et al.  , 2007, ApJS , 173, 267 https://doi.org/10.1086/519218 CrossRef Search ADS   Salpeter E. E., 1955, ApJ , 121, 161 https://doi.org/10.1086/145971 CrossRef Search ADS   Santini P. et al.  , 2012, A&A , 540, A109 CrossRef Search ADS   Sargent M. T., Béthermin M., Daddi E., Elbaz D., 2012, ApJ , 747, L31 https://doi.org/10.1088/2041-8205/747/2/L31 CrossRef Search ADS   Schechter P., 1976, ApJ , 203, 297 https://doi.org/10.1086/154079 CrossRef Search ADS   Schreiber C. et al.  , 2015, A&A , 575, A74 CrossRef Search ADS   Silverman J. D. et al.  , 2008, ApJ , 679, 118 https://doi.org/10.1086/529572 CrossRef Search ADS   Silverman J. D. et al.  , 2009, ApJ , 696, 396 https://doi.org/10.1088/0004-637X/696/1/396 CrossRef Search ADS   Stanley F., Harrison C. M., Alexander D. M., Swinbank A. M., Aird J. A., Del Moro A., Hickox R. C., Mullaney J. R., 2015, MNRAS , 453, 591 https://doi.org/10.1093/mnras/stv1678 CrossRef Search ADS   Stanley F. et al.  , 2017, MNRAS , 742, 2221 CrossRef Search ADS   Vasudevan R. V., Fabian A. C., 2007, MNRAS , 381, 1235 https://doi.org/10.1111/j.1365-2966.2007.12328.x CrossRef Search ADS   Veale M., White M., Conroy C., 2014, MNRAS , 445, 1144 https://doi.org/10.1093/mnras/stu1821 CrossRef Search ADS   Wang T. et al.  , 2017, A&A , 601, A63 CrossRef Search ADS   Weigel A. K. et al.  , 2017, ApJ , 845, 145 https://doi.org/10.3847/1538-4357/aa8097 CrossRef Search ADS   © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Monthly Notices of the Royal Astronomical Society Oxford University Press

# Evidence for a mass-dependent AGN Eddington ratio distribution via the flat relationship between SFR and AGN luminosity

, Volume 476 (1) – May 1, 2018
15 pages

/lp/ou_press/evidence-for-a-mass-dependent-agn-eddington-ratio-distribution-via-the-JCmHkDyHhn
Publisher
The Royal Astronomical Society
ISSN
0035-8711
eISSN
1365-2966
D.O.I.
10.1093/mnras/sty219
Publisher site
See Article on Publisher Site

### Abstract

Abstract The lack of a strong correlation between AGN X-ray luminosity (LX; a proxy for AGN power) and the star formation rate (SFR) of their host galaxies has recently been attributed to stochastic AGN variability. Studies using population synthesis models have incorporated this by assuming a broad, universal (i.e. does not depend on the host galaxy properties) probability distribution for AGN specific X-ray luminosities (i.e. the ratio of LX to host stellar mass; a common proxy for Eddington ratio). However, recent studies have demonstrated that this universal Eddington ratio distribution fails to reproduce the observed X-ray luminosity functions beyond z ∼ 1.2. Furthermore, empirical studies have recently shown that the Eddington ratio distribution may instead depend upon host galaxy properties, such as SFR and/or stellar mass. To investigate this further, we develop a population synthesis model in which the Eddington ratio distribution is different for star-forming and quiescent host galaxies. We show that, although this model is able to reproduce the observed X-ray luminosity functions out to z ∼ 2, it fails to simultaneously reproduce the observed flat relationship between SFR and X-ray luminosity. We can solve this, however, by incorporating a mass dependency in the AGN Eddington ratio distribution for star-forming host galaxies. Overall, our models indicate that a relative suppression of low Eddington ratios (λEdd ≲ 0.1) in lower mass galaxies (M* ≲ 1010 − 11 M⊙) is required to reproduce both the observed X-ray luminosity functions and the observed flat SFR/X-ray relationship. galaxies: active, galaxies: evolution, galaxies: luminosity function, mass function, galaxies: statistics 1 INTRODUCTION Most galaxies host a central super massive black hole (SMBH), the masses of which display a tight correlation with their host bulge masses, implying a co-evolution between SMBHs and their host galaxies (e.g. Kormendy & Richstone 1995; Magorrian et al. 1998; Merritt 2000; Kormendy, Bender & Cornell 2011). Furthermore, the similar redshift evolution of the total SMBH accretion rate density and the total star formation rate (SFR) density strongly suggests the existence of a link between SMBH accretion and star formation at galactic scales (e.g. Heckman et al. 2004; Merloni, Rudnick & Di Matteo 2004; Silverman et al. 2008; Aird et al. 2010; Assef et al. 2011). However, the precise connection between the two processes remains poorly understood. One means of connecting SMBH and galaxy growth that has gained popular support over the past two decades is for accreting SMBHs (observed as active galactic nuclei, AGNs) to directly influence the SFRs of their hosts via a variety of feedback mechanisms (see the review of Fabian 2012, for details on the feedback mechanisms). However, a major difficulty in identifying the precise role that SMBH accretion has on influencing SFR is the stochastic nature of AGNs (e.g. Aird et al. 2013; Hickox et al. 2014, see also section 4.3 of Stanley et al. 2015). It is argued that this randomness is the reason why most studies that have explored the relationship between SFR and X-ray luminosity (a proxy for SMBH accretion power) have reported no evidence of a strong correlation between these parameters, at least for moderate luminosity AGNs, which form the majority of the population (i.e. 1042 < LX < 1045 erg s−1; e.g. Lutz et al. 2010; Harrison et al. 2012; Mullaney et al. 2012; Rosario et al. 2012; Rovilos et al. 2012; Santini et al. 2012; Azadi et al. 2015; Stanley et al. 2015; Lanzuisi et al. 2017). Recently, using optically selected AGNs, Stanley et al. (2017) reported an enhancement of SFR among the highest luminosity AGNs (i.e. LX >1045 erg s−1). However, they demonstrated that this is a direct consequence of the most luminous AGNs residing in more massive host galaxies, meaning the enhanced SFR is a consequence of the relationship between stellar mass and SFR (hereafter referred to as the main sequence, MS; e.g. Daddi et al. 2007; Rodighiero et al. 2011; Schreiber et al. 2015). A major problem in investigating the relationship between AGN and SFR is that observations are often subjected to biases (e.g. flux limited samples, rarity of high accretion rate AGNs). Therefore, many studies investigating the connection between SMBH and galaxy growth have adopted a modelling approach (e.g. Aird et al. 2013; Conroy & White 2013; Veale, White & Conroy 2014; Caplar, Lilly & Trakhtenbrot 2015; Jones et al. 2017; Weigel et al. 2017). In these models, the stochastic nature of AGNs is often implemented via a broad specific AGN X-ray luminosity (i.e. the X-ray luminosity relative to the host stellar mass; LX/M*) distribution that, once combined with a mass function, reproduces the observed X-ray luminosity function of AGNs. Overall, these studies find that below z ∼ 1 the specific X-ray luminosity distribution is well-represented by a broad, universal (i.e. the same distribution, irrespective of SFR or stellar mass) distribution, described by a broken power-law (or similar; e.g. a Schechter 1976 function), the normalization of which increases with increasing redshifts (e.g. Aird et al. 2013; Veale et al. 2014; Jones et al. 2017). These models are also often successful at reproducing the observed flat relationship between SFR and X-ray luminosity (Hickox et al. 2014; Veale et al. 2014; Stanley et al. 2015). However, Jones et al. (2017) recently found that using a universal broad distribution for the specific X-ray luminosity distribution of AGNs cannot reproduce the X-ray luminosity functions of Aird et al. (2010) beyond z ∼ 1.2. They argue for the need of a more complicated specific X-ray luminosity distribution of AGNs beyond that redshift. Furthermore, recent studies have reported that the specific X-ray luminosity distribution for AGNs in star-forming hosts significantly differs from that of AGNs in quiescent hosts (e.g. Georgakakis et al. 2014; Wang et al. 2017; Aird et al. 2018). Overall, they report a higher contribution from star-forming AGN hosts to the total specific X-ray luminosity distribution out to z ∼ 2. Motivated by the above findings, we explore a model for the specific X-ray luminosity distribution split between star-forming and quiescent galaxies. After demonstrating that this simple model cannot simultaneously reproduce both the X-ray luminosity functions and the flat relationship between SFR and X-ray luminosity out to z ∼ 2, we introduce a mass dependency for the specific X-ray luminosity distribution of star-forming AGN hosts, as recently observed (e.g. Bongiorno et al. 2016; Georgakakis et al. 2017; Aird et al. 2018). Hereafter, we refer to our model specific X-ray luminosity distribution split between star-forming and quiescent galaxies as ‘mass-independent’, and to our model that includes the mass dependency for the specific X-ray luminosity distribution of AGNs in star-forming galaxies as ‘mass-dependent’. In Section 2, we describe our method to infer the specific X-ray luminosity distribution for both the mass-independent and the mass-dependent model. Subsequently, in the same section, we present how we incorporate host galaxy properties. We then report the results for the mass-independent model in Section 3.1, and for the mass-dependent model in Section 3.2. We discuss the general implications of our results in Section 4, and conclude in Section 5. Throughout, we assume the specific X-ray luminosity (i.e. LX/M*) to be proportional to Eddington ratio (i.e. the ratio between bolometric AGN luminosity and Eddington luminosity, LAGN/LEdd), converting the X-ray 2–10 keV luminosity into bolometric AGN luminosity using a universal conversion factor of 22.4 (median value found in Vasudevan & Fabian 2007, and based on a local AGN sample with LX = 1041 − 46 erg s−1), and assuming a ratio of 0.002 between SMBH masses and stellar masses (Marconi & Hunt 2003). Therefore, we have,   $${\rm \lambda _{\rm Edd}} = \frac{22.4 \times {\rm L_{\rm X}}}{10^{38.1}\,{\rm erg\,s^{-1}} \times 0.002\ \frac{{\rm M_\ast }}{{\rm \mathrm{M}_{\odot }}}},$$ (1)where λEdd is the Eddington ratio, and LX is the intrinsic (absorption-corrected) 2–10 keV X-ray luminosity. We stress that the important parameter we investigate is the specific X-ray luminosity (LX/M*), and that the Eddington ratio is solely used for its greater familiarity. Throughout, we assume a WMAP–7 year cosmology (Larson et al. 2011) and a Salpeter (1955) initial mass function (IMF) when calculating galaxy stellar masses and SFRs. 2 METHODOLOGY Many recent studies have used Population Synthesis Models (PSMs) to investigate the connection between AGNs and host galaxies. PSMs randomly draw properties from probability distribution functions (e.g. Aird et al. 2013; Hickox et al. 2014; Stanley et al. 2015; Jones et al. 2017). Using PSMs, we are able to generate populations of N massive mock galaxies, with stellar masses M*, drawn from a mass distribution, for which we randomly allocate Eddington ratios, λEdd, following an Eddington ratio distribution, and, using equation (1), derive X-ray luminosities, LX, the histogram of which is proportional to the model X-ray luminosity function. The key point is to optimize the Eddington ratio distribution until the model X-ray luminosity function matches the observed one (e.g. Aird et al. 2013; Jones et al. 2017). This can be repeated at different redshifts to investigate any redshift dependence of the Eddington ratio distribution. Because it implies generating a new population of galaxies while optimizing the Eddington ratio distribution, it can be computationally expensive. Therefore, we will use an analytical approach (similar to those of e.g. Veale et al. 2014; Caplar et al. 2015) by directly considering the various distributions as probability density functions. In this study, we fit to the X-ray luminosity functions of Aird et al. (2015) and use the stellar mass functions of Davidzon et al. (2017). The latter are shown in Fig. 1 for various redshifts, out to z = 2.5. Figure 1. View largeDownload slide Mass functions employed in our models. Top panel: mass functions for star-forming galaxies as reported in Davidzon et al. (2017). Each colour corresponds to a different redshift bin, out to z = 2.5 (see key). Bottom panel: same as top panel but for the quiescent galaxies. Figure 1. View largeDownload slide Mass functions employed in our models. Top panel: mass functions for star-forming galaxies as reported in Davidzon et al. (2017). Each colour corresponds to a different redshift bin, out to z = 2.5 (see key). Bottom panel: same as top panel but for the quiescent galaxies. In our analytical approach, we first assume a model λEdd distribution, which is described by a set of parameters θ (for example, for a power law θ = {N, α}, where N is the normalization and α the slope). Our aim is to iterate over θ to identify the solution that best fits the observed X-ray luminosity function. The model X-ray luminosity function is generated by analytically combining the observed mass function with the λEdd distribution given by the set of parameters θ. To achieve this, we treat the X-ray luminosity function as a probability distribution, p(LX). The probability of observing an AGN of luminosity LX = X in a galaxy of mass M* is given by1  $$p(L_{\rm X}=X) \equiv p\left(\lambda _{\rm Edd} = \frac{X}{M_\ast }\right),$$ (2)where p(λLedd = X/M*) is the probability of λEdd = X/M* calculated with the current set of model parameters, θ. We derive the total probability of observing LX = X in all galaxies as the sum of these individual probabilities weighted by the galaxy mass function (which we also treat as a probability distribution, p(M*)),   \begin{eqnarray} p(L_{\rm X}=X) = \sum _{Y=M_{\rm Min}}^{M_{\rm Max}} p\left(\lambda _{\rm Edd} = \frac{X}{Y} | M_\ast =Y\right)\times p(M_\ast =Y),\nonumber \\ \end{eqnarray} (3)where p(M* = Y) is the probability of M*=Y (or the mass function evaluated at Y). The total model luminosity function derived from a given set of parameters, θ, is then simply the combination of these probabilities (i.e. for LX = 1041–1046 erg s−1 in bins of 0.1 dex). By taking this approach, we avoid having to model N galaxies (where $$N\sim \mathcal {O}(10^9)$$ to generate sufficient numbers of luminous galaxies), and can simply split our mass function into $$\mathcal {O}(10^3)$$ bins between MMin = 108 M⊙ to MMax = 1014 M⊙ (where MMin and MMax are chosen to cover the full range of stellar masses). Since we adopt an iterative process to identify the best-fitting parameters (see later in this subsection), this factor of $$\mathcal {O}(10^{6})$$ reduction in the number of required calculations per iteration dramatically reduces the time taken to converge. To optimize the parameters that define our Eddington ratio distribution, we performed maximum likelihood estimation. Assuming Gaussian distributions for the uncertainties, the log-likelihood of a set of parameters θ is defined by −0.5χ2, where χ2 is the chi-square between the model and the observed X-ray luminosity functions. We used flat proper (i.e. with finite boundaries) prior distributions for each parameter that defines our model, and checked the posterior distributions to ensure that they are not constrained in any way by the limits of our prior distributions. The parameter space (defined by the flat prior distributions) was explored using a Monte Carlo Markov Chain (MCMC) fully implemented in a python application programming interface, emcee2 (Foreman-Mackey et al. 2013), that uses the affine invariant MCMC ensemble sampler of Goodman & Weare (2010). We ran our MCMC code using 100 walkers for 10 000 steps as a burn-in phase to find the global solution (i.e. these 10 000 steps are then discarded), and for 5000 steps as a production phase (i.e. these 5000 steps are kept for analysis). The convergence of each chain was assessed using the Gelman & Rubin (1992) test. We then extracted the best parameters by fitting a Gaussian function to each sampled posterior distribution, taking the mean and the standard deviation as the best estimate and the 1σ uncertainties, respectively, of each parameter. When the posterior distribution of a parameter was flat below or above a given value, we assumed that it is an upper or a lower limit, respectively (again, we ensured that values of these limits are not affected by our prior distributions). This optimization method was applied to infer the Eddington ratio distribution for both our mass-independent and mass-dependent models. 2.1 A mass-independent λEdd distribution As highlighted in the introduction, a number of studies have recently reported a difference in the λEdd distribution of quiescent and star-forming host galaxies (e.g. Georgakakis et al. 2014; Wang et al. 2017). We therefore take this ‘two-component’ distribution as our starting point (i.e. we do not attempt to model a single universal λEdd distribution for all galaxies). Most studies use a broken power-law (or related functions) to represent the Eddington ratio distribution of AGNs (e.g. Aird et al. 2013; Veale et al. 2014; Jones et al. 2017; Weigel et al. 2017). As such, we take a conservative approach by assuming a broken power-law for both the star-forming and the quiescent component of our λEdd distribution. Each broken power law is defined as   $$p(\lambda _{\rm Edd}) = \frac{A}{\left(\frac{\lambda _{\rm Edd}}{\lambda _{\rm break}} \right)^{\gamma _{1}}+\left(\frac{\lambda _{\rm Edd}}{\lambda _{\rm break}} \right)^{\gamma _{2}}},$$ (4)where $$p({\rm \lambda _{\rm Edd}})$$ is the probability of Eddington ratio, A is the normalization, λbreak is the position of the break, and γ1 and γ2 are the slopes at low λEdd (i.e. below the break) and high λEdd (i.e. above the break), respectively. As such, each of our broken power law components are represented by four parameters, giving a total of eight free parameters to optimize for this mass-independent model (see Section 2 for the optimization method). Following empirical results reported in, for example, Georgakakis et al. (2014), Wang et al. (2017), and Aird et al. (2018), prior to optimization, we assume that both the normalization, ASF, and the position of the break, $$\lambda _{\rm break}^{\rm SF}$$, of the star-forming component are each always higher than their quiescent analogues, AQui and $$\lambda _{\rm break}^{\rm Qui}$$, respectively. This also has the benefit of reducing the degeneracy of the model. However, it is important to stress that we do not assume any specific parameter values for the ratio between the two normalizations and the two positions of the breaks, we simply incorporate this information in the prior distributions by excluding some parts of the parameter space. Incorporating the aforementioned assumptions, we optimized the eight free parameters of our mass-independent model. To investigate any redshift evolution of the Eddington ratio distribution, we performed the optimization at z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.2, z = 1.5, z = 1.7, z = 2.0, z = 2.3, z = 2.5, z = 2.5, z = 2.7, z = 3.0, and z = 3.5, interpolating the mass functions of Davidzon et al. (2017) split between star-forming and quiescent galaxies, and the total X-ray luminosity functions of Aird et al. (2015) at these redshifts. We stress that the Bayesian methodology adopted in Aird et al. (2015) to derive the observed X-ray luminosity functions (see Aird et al. (2015) for details on their method) means that they consider all X-ray detected sources when deriving their X-ray luminosity functions (which our analysis aims to fit). Therefore, there is no implied mass cut beyond the fact that their sources are detected in X-rays. 2.2 A mass-dependent λEdd distribution Motivated by recent studies reporting that the λEdd distribution for star-forming galaxies is mass-dependent (e.g. Bongiorno et al. 2016; Georgakakis et al. 2017; Aird et al. 2018), we develop a second model that can accommodate three different λEdd distributions in three different mass bins for star-forming galaxies. As with our mass-independent model, we still split in terms of quiescent and star-forming galaxies, but for the star-forming component, we assume instead three different broken power-laws (see equation 4), one for each of our mass bins, which we define as: low mass (8 < log(M*/M⊙) < 10), medium mass (10 < log(M*/M⊙) < 11), and high mass (11 < log(M*/M⊙) < 12). The choice of the boundaries for our mass bins is arbitrary and independent of any empirical studies, yet aims to cover a large range of stellar masses (i.e. 108 <M*/M⊙ <1012). Since assuming three broken power-laws generates 12 free parameters, we again rely on some assumptions based on the results of Aird et al. (2018) to help prevent degeneracies. First, we require that each broken power-law breaks at different λEdd values, in such a way that a higher mass bin breaks at a lower value of λEdd than its lower mass neighbour. Secondly, we assume that the three Eddington ratio distributions share the same slope at high Eddington ratios. Finally, we assume the normalization of each of the three broken power-laws is such that the λEdd distribution above the break is always coincident (see Fig. 8). In making these assumptions, we reduce the parameter space to eight free parameters: one normalization, ASF, (from which the others are derived via our third assumption); three break positions, $$\lambda _{\rm break}^{\rm SF}$$, (one for each mass bin and ordered according to our first assumption); three power-law slopes below the break, $$\gamma _{1}^{\rm SF}$$, (which are unconstrained); a single shared power-law slope above the break, $$\gamma _2^{\rm SF}$$. We stress that all these assumptions are inspired by empirical findings reported by Aird et al. (2018). We will demonstrate in Section 3.1.1 that our previous model (i.e. the mass-independent model, see Section 2.1) is able to reproduce the observed X-ray luminosity functions out to z ∼ 2, and recreates both their star-forming and quiescent galaxy components in good agreement with observations, at least out to z ∼ 1. As such, for this mass-dependent model, we assume that the Eddington ratio distribution for quiescent galaxies is unchanged from our previous model, and only update the star-forming component of the model Eddington ratio distribution by implementing the aforementioned mass-dependency. To do this, we optimize the eight free parameters that describe our mass-dependent Eddington ratio distribution for star-forming galaxies. Again, to investigate any redshift evolution of the mass-dependent Eddington ratio distribution, we performed the optimization at z = 0.1, z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.3, z = 1.5, and z = 1.7, interpolating the mass functions for star-forming galaxies of Davidzon et al. (2017) and the star-forming component of the X-ray luminosity functions as extracted from our mass-independent model. This second model relies on the ability of our first, mass-independent, model to split the X-ray luminosity functions in terms of star-forming and quiescent galaxies. As we show in Section 3.1.1, our mass-independent model is unable to reproduce the total X-ray luminosity function beyond z ∼ 2, and thus cannot be used to reliably split the luminosity function between star-forming and quiescent galaxies. As such, we cannot extend our mass-dependent model beyond z ∼ 2. 2.3 Generating galaxies with AGNs As we aim to investigate the connection between AGN accretion and host star formation activity, we must model the star-forming properties of the host galaxies as well as the Eddington ratios of the AGNs. Since the distributions of specific SFRs (i.e. SFR relative to the host stellar mass; sSFR) are now well-defined, we can simply adopt published models to generate our galaxy populations. As such, for this part of the model, we use a PSM to attribute SFRs to our host galaxies. We use the PSM outlined in Bernhard et al. (2014), which generates a population of mock star-forming galaxies using the mass functions of Ilbert et al. (2013) and allocates SFRs using the relationship between stellar masses and SFRs, or MS ( e.g. Daddi et al. 2007; Salim et al. 2007; Elbaz et al. 2011; Rodighiero et al. 2011; Sargent et al. 2012; Schreiber et al. 2015, for the MS) reported in Rodighiero et al. (2011). The SFRs are then randomly scattered around the sSFR distribution of Sargent et al. (2012). This model successfully reproduces the observed infrared and ultraviolet luminosity functions up to z ∼ 6 (see Bernhard et al. 2014, and references therein). However, for this work, we slightly update the empirical relationships used in Bernhard et al. (2014). Briefly, we change the mass functions to the more recent ones reported in Davidzon et al. (2017), we update the MS sSFR distribution to that reported in Schreiber et al. (2015), and lastly, add a population of quiescent galaxies (using the quiescent galaxy mass functions of Davidzon et al. 2017) with sSFRs at least a factor of ten below that of the MS (Ilbert et al. 2013). We confirm that this updated model still reproduces the infrared and ultraviolet luminosity functions. X-ray luminosities are then allocated using our Eddington ratio distributions derived in Section 2. Using this PSM, we generate a population of 33 925 192 galaxies (i.e. 32 939 834 – or 97 per cent – star-forming and 985 358 – or 3 per cent – quiescent galaxies), with stellar masses in the range 108 <M*/M⊙ <1014. This corresponds to a 50 square degrees blank-field survey out to z = 3. We use a very high (arguably un-physical) upper limit on our range of stellar masses to rule-out the possibility of being affected by any boundary effects. The large relative number of star-forming compared to quiescent galaxies in our model is related to the differences of their respective mass functions, particularly at masses M* ≲ 1010 M⊙ and z ≳ 1. As shown in Fig. 1, at these lower masses and higher redshifts, the mass function of star-forming galaxies steeply rises, while that of quiescent galaxies sharply decreases. This leads to a large difference in the absolute number of star-forming galaxies compared to that of quiescent galaxies at lower masses (we have a cut at M* = 108 M⊙ in our model). Above a stellar mass of M* >1010 M⊙, we have a total of 3 309 201 galaxies, among which 2 646 102 (80 per cent) are star-forming and 663 099 (20 per cent) are quiescent. 3 RESULTS In this section, we first present the results of our mass-independent model, followed by the results of our mass-dependent model. We demonstrate that while our mass-independent model can successfully reproduce the observed X-ray luminosity functions out to z ≈ 2, it cannot at the same time reproduce the observed flat relationship between SFR and AGN luminosity (e.g. Rosario et al. 2012; Stanley et al. 2015). Instead, we show that this can be resolved by introducing a mass dependency in the λEdd distribution of AGNs hosted in star-forming galaxies. 3.1 The mass-independent model 3.1.1 X-ray luminosity functions We show in Fig. 2, the fit to the X-ray luminosity functions for our mass-independent model. Up to z = 1.75 our model X-ray luminosity functions are in very good agreement with those of Aird et al. (2015). However, beyond that redshift, our mass-independent model is unable to reproduce the observed X-ray luminosity functions. Fig. 2 also illustrates the model X-ray luminosity functions split between star-forming and quiescent galaxies. Overall, we find that the star-forming galaxies always dominate the X-ray luminosity functions at LX ≳ 1042.5 erg s−1, which corresponds to the bulk of the AGN population. The contribution from quiescent galaxies in our mass-independent model, however, is only significant at lower X-ray luminosities (i.e. LX ≲ 1042 erg s−1). Figure 2. View largeDownload slide Fit to the X-ray 2–10 keV luminosity functions of Aird et al. (2015) when optimizing our mass-independent λEdd distribution model. Each panel represents a different redshift bin out to z = 3, and as indicated at their bottom right-hand side. The downward and upward empty black triangles are the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The thick orange line in each panel is for the model total X-ray luminosity function given by our mass-independent model, the blue dashed line is for that of the star-forming component and the red triple dashed line shows the quiescent component. The blue empty stars and the red empty circles at z < 0.90 are for the star-forming and the quiescent component of the X-ray luminosity functions, respectively, as reported in Georgakakis et al. (2014) out to z ∼ 1. We hatched the two panels at z = 2.25 and z = 3.0 to indicate that our mass-independent λEdd distribution model is unable to reproduce the observed X-ray luminosity functions at these redshifts (see Section 3.1.1). Figure 2. View largeDownload slide Fit to the X-ray 2–10 keV luminosity functions of Aird et al. (2015) when optimizing our mass-independent λEdd distribution model. Each panel represents a different redshift bin out to z = 3, and as indicated at their bottom right-hand side. The downward and upward empty black triangles are the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The thick orange line in each panel is for the model total X-ray luminosity function given by our mass-independent model, the blue dashed line is for that of the star-forming component and the red triple dashed line shows the quiescent component. The blue empty stars and the red empty circles at z < 0.90 are for the star-forming and the quiescent component of the X-ray luminosity functions, respectively, as reported in Georgakakis et al. (2014) out to z ∼ 1. We hatched the two panels at z = 2.25 and z = 3.0 to indicate that our mass-independent λEdd distribution model is unable to reproduce the observed X-ray luminosity functions at these redshifts (see Section 3.1.1). Beyond z ∼ 2, our model cannot reproduce the observed X-ray luminosity functions. At the position of the break, our model X-ray luminosity functions at these redshifts are at least a factor of three below the observed ones. We also find that the contribution from quiescent galaxies is very low (by at least a factor of ten) compared to that of the star-forming galaxies. To explain the reasons for this, we examine the stellar mass functions of Davidzon et al. (2017; see Fig. 1). They report that the contribution from quiescent galaxies to the mass function decreases as the redshift increases (the knee of the mass function for quiescent galaxies is at least a factor of ten in normalization below that of the star-forming galaxies at z ≳ 2; see Fig. 1). In our model, this lower contribution from the quiescent galaxies can be compensated by an increase in the normalization of the model Eddington ratio distribution for quiescent galaxies (as our model relates the mass function to the Eddington ratio distribution). This would increase the contribution from quiescent galaxies to the total model X-ray luminosity functions in order to match the observed ones. However, in our model, the normalization of the quiescent galaxy component cannot exceed that of the star-forming component (see assumptions in Section 2.1). While looking at the normalizations of both the star-forming (ASF) and the quiescent (AQui) component in our model, we find that beyond z ∼ 2, they are similar (i.e. ASF ∼ AQui). As a consequence, the contribution from quiescent galaxies to the total X-ray luminosity functions cannot increase enough to match the observed X-ray luminosity functions. If we relax this assumption, our model finds that the Eddington ratio distribution beyond z ∼ 2 is fully dominated by quiescent galaxies, with the normalization of the Eddington ratio distribution for star-forming galaxies to that of quiescent galaxies as low as ∼0.01 at z = 2.25 and ∼0.001 at z = 2.5. This is at odds with empirical measurements of the Eddington ratio distributions of, for example, Georgakakis et al. (2014) and Wang et al. (2017), that suggest a more equal splitting between quiescent and star-forming galaxies, as well as infrared studies that find large numbers of high-redshift AGNs residing in star-forming galaxies (e.g. Mullaney et al. 2015; Stanley et al. 2015; Netzer et al. 2016). Finally, as a further check, we include in Fig. 2 the measured X-ray luminosity functions separated into star-forming and quiescent galaxies for z < 1 as reported in Georgakakis et al. (2014). Despite not including this information during our optimization, we find good agreement between our model and these observed X-ray luminosity functions of star-forming and quiescent galaxies, increasing our confidence in the model at these lower redshifts. 3.1.2 The Eddington ratio distribution We show in Section 3.1.1 that our mass-independent model is able to reproduce the observed X-ray luminosity function out to z ∼ 2. We now present in Fig. 3 the Eddington ratio distributions at various redshifts after optimization of our mass-independent model (normalized to coincide at λEdd= 10−3). The shape of our total Eddington ratio distribution (see left-hand panel in Fig. 3) is consistent with a broken power-law, as requested by our model and in agreement with previous studies (e.g. Aird et al. 2013; Jones et al. 2017). However, we find a significant difference in the slope below the break (γ1 in equation 4) of the Eddington ratio distribution for star-forming hosts compared to that of quiescent hosts (see central and right-hand panels in Fig. 3 for the Eddington ratio distribution of star-forming and quiescent host galaxies, respectively). We report that below z ∼ 2, the slope at low-λEdd for the star-forming component is rising (although we are only able to place upper limits). As a consequence, our model Eddington ratio distribution of AGNs in star-forming galaxies is consistent with a peaky3 distribution, similar to the light-bulb shaped distribution sometimes explored in past studies (e.g. Veale et al. 2014; Stanley et al. 2015). This suggests that AGNs hosted in star-forming galaxies have typical λEdd higher than AGNs in quiescent hosts, yet the λEdd of AGNs in quiescent galaxies span a broader range of values. Finally, we notice an overall shift of the knee of the Eddington ratio distribution to higher λEdd values as redshift increases. This implies a typical AGN activity that increases as redshift increases. At z ≳ 2, we find that the shape of the Eddington ratio distribution of AGNs in star-forming galaxies becomes consistent with a broken power-law, in contrast to the peaky Eddington ratio distribution reported for lower redshifts. However, as demonstrated in Section 3.1.1 our mass-independent model is unable to reproduce the X-ray luminosity functions at these higher redshifts. We report in Table 1, the redshift evolution of the parameters that define these Eddington ratio distributions. Figure 3. View largeDownload slide The total mass-independent Eddington ratio distribution (left-hand panel), its star-forming component (central panel), and its quiescent component (right-hand panel) evaluated at z = 0.2, z = 0.5, z = 1.0, z = 1.5, z = 2.0, and z = 2.5 (colour coded; see key). The downward arrows indicate the overall effects of the upper and lower limits found in the parameters that describes these distributions. For the star-forming and quiescent panels, we also indicate with faint grey lines the corresponding total Eddington ratio distribution. The normalization is such that the total Eddington ratio distribution at z = 0.2 integrates to unity, applying an arbitrary cut at λEdd = 10−7. For ease of comparison, we then re-normalize all the distributions at higher redshifts such that their values at log(λEdd)= −3.0 coincide with that at z = 0.2, hence the arbitrary unit of the ordinate axis. Figure 3. View largeDownload slide The total mass-independent Eddington ratio distribution (left-hand panel), its star-forming component (central panel), and its quiescent component (right-hand panel) evaluated at z = 0.2, z = 0.5, z = 1.0, z = 1.5, z = 2.0, and z = 2.5 (colour coded; see key). The downward arrows indicate the overall effects of the upper and lower limits found in the parameters that describes these distributions. For the star-forming and quiescent panels, we also indicate with faint grey lines the corresponding total Eddington ratio distribution. The normalization is such that the total Eddington ratio distribution at z = 0.2 integrates to unity, applying an arbitrary cut at λEdd = 10−7. For ease of comparison, we then re-normalize all the distributions at higher redshifts such that their values at log(λEdd)= −3.0 coincide with that at z = 0.2, hence the arbitrary unit of the ordinate axis. Table 1. Redshift evolution of the parameters that describe the Eddington ratio distribution for star-forming and quiescent galaxies, given by equation (4). SF and Qui labels stand for star-forming and quiescent, respectively. The slopes and intercepts are given for an evolution in (1+z). Uncertainties on the parameters are the standard deviation of 1000 Monte Carlo realisations (i.e. 1σ). Parameters  Intercepts  Slopes    log(ASF)  −2.07±0.19  0.13±0.07  for z < 1.7      −1.15±0.03  for z > 1.7  log(AQui)  −2.80±0.23  0.37±0.10  for z < 1.7      −1.20±0.11  for z > 1.7  log$$(\lambda _{\rm break}^{\rm SF})$$  −2.36±0.06  0.65±0.02    log$$(\lambda _{\rm break}^{\rm Qui})$$  −2.62±0.16  0.66±0.06    $$\gamma _1^{\rm SF}$$  ≲ − 3.0  0.0  for z < 1.7    −0.96±0.51  0.25±0.12  for z > 1.7  $$\gamma _1^{\rm Qui}$$  0.28±0.05  0.07±0.02    $$\gamma _2^{\rm SF}$$  2.29±0.04  0.01±0.02  for z > 1.7      0.15±0.03  for z > 1.7  $$\gamma _2^{\rm Qui}$$  ≳2.7  0.0  for z < 1.7    3.75±2.33  −0.53±0.63  for z > 1.7  Parameters  Intercepts  Slopes    log(ASF)  −2.07±0.19  0.13±0.07  for z < 1.7      −1.15±0.03  for z > 1.7  log(AQui)  −2.80±0.23  0.37±0.10  for z < 1.7      −1.20±0.11  for z > 1.7  log$$(\lambda _{\rm break}^{\rm SF})$$  −2.36±0.06  0.65±0.02    log$$(\lambda _{\rm break}^{\rm Qui})$$  −2.62±0.16  0.66±0.06    $$\gamma _1^{\rm SF}$$  ≲ − 3.0  0.0  for z < 1.7    −0.96±0.51  0.25±0.12  for z > 1.7  $$\gamma _1^{\rm Qui}$$  0.28±0.05  0.07±0.02    $$\gamma _2^{\rm SF}$$  2.29±0.04  0.01±0.02  for z > 1.7      0.15±0.03  for z > 1.7  $$\gamma _2^{\rm Qui}$$  ≳2.7  0.0  for z < 1.7    3.75±2.33  −0.53±0.63  for z > 1.7  Notes: Slopes and intercepts are given for an evolution as a function of (1+z). The intercept for z > 1.7, when assuming a break in the z evolution of the parameter, is given by the continuity at z = 1.7 (i.e. [intercept forz > 1.7] = (1 + 1.7) × ([slope forz < 1.7] − [slope forz > 1.7]) + [intercept forz < 1.7]). View Large 3.1.3 SFR in bins of X-ray luminosities As previously highlighted, a key test of any model that describes the Eddington ratio distribution is whether it can reproduce other observed features of the AGN population besides just the (X-ray) luminosity functions. Here, we test our PSM (which incorporates pre-defined sSFR distributions and our own optimized mass-independent Eddington ratio distribution; see Section 2.3) by assessing whether it can reproduce the observed flat relationship between SFR and X-ray luminosity as reported in many observational studies (e.g. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). We do this by calculating the mean-average SFRs of the AGNs in our PSM in bins of 0.5 dex in X-ray luminosity, taking into account both star-forming and quiescent galaxies (see Section 2.3). As shown in Fig. 4, our mass-independent model predicts a strong correlation between these SFRs and X-ray luminosities (with an average slope of ∼0.6), at least up to z ∼ 2. This is in good agreement with empirical results for our lowest redshift bin (i.e. z = 0.35), but strongly contradicts the flat observed relationship between SFR and X-ray luminosity at higher redshifts (e.g. Rosario et al. 2012; Stanley et al. 2015). It is important to note that other studies also find a similar strongly increasing relationship between SFR and X-ray luminosity when employing a peaky (i.e. similar in shape than our peaky distribution; see central panel in Fig. 3) Eddington ratio distribution for star-forming galaxies (e.g. Veale et al. 2014; Stanley et al. 2015). We will show in Section 4 that this strong correlation is due to the AGN host mass distribution and the SFR/stellar-mass relationship (i.e. MS) induced by the peaky Eddington ratio distribution found for our mass-independent model. Figure 4. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosity at z = 0.35, z = 0.65, z = 1.15, and z = 2. The black thick lines shows the correlation predicted by our PSM that incorporates AGNs using our mass-independent Eddington ratio distribution split between star-forming and quiescent galaxies. The error bars show the 1σ uncertainties on the mean for average SFR. The thick black dashed line is identical than the thick black line but for galaxies with M* ≳ 109.5 M⊙. Various symbols are for a compilation of published SFR/X-ray relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line at z = 1.15 shows the relationship reported in Silverman et al. (2009). Finally, the black dashed line represents the correlation displayed by AGN-dominated systems in Netzer (2009). Our mass-independent model predicts a significant correlation between SFR and X-ray luminosity, which is at odds with empirical results for z ≳ 0.5. Figure 4. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosity at z = 0.35, z = 0.65, z = 1.15, and z = 2. The black thick lines shows the correlation predicted by our PSM that incorporates AGNs using our mass-independent Eddington ratio distribution split between star-forming and quiescent galaxies. The error bars show the 1σ uncertainties on the mean for average SFR. The thick black dashed line is identical than the thick black line but for galaxies with M* ≳ 109.5 M⊙. Various symbols are for a compilation of published SFR/X-ray relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line at z = 1.15 shows the relationship reported in Silverman et al. (2009). Finally, the black dashed line represents the correlation displayed by AGN-dominated systems in Netzer (2009). Our mass-independent model predicts a significant correlation between SFR and X-ray luminosity, which is at odds with empirical results for z ≳ 0.5. As SFR is correlated to stellar mass via the MS, it is possible that the empirical flat relationship between SFR and X-ray luminosity is a consequence of an observational bias towards higher mass galaxies. To test this, we recalculate from our model the average SFR in bins of X-ray luminosity, but now excluding galaxies with M* <109.5 M⊙. This mass limit roughly corresponds to that reached in the highest redshift bins of Stanley et al. (2015). We show in Fig. 4 that although this mass cut causes the model SFR/X-ray relationship to rise slightly at z = 2.0 and LX ≲ 1043.5 erg s−1, it is not sufficient to reproduce the observed flat relationship between SFR and X-ray luminosity. 3.2 The mass-dependent model Since our mass-independent model fails to reproduce the X-ray luminosity functions at z ≳ 2 and is unable to reproduce the flat relationship between SFR and X-ray luminosity (at z > 0.5), we consider whether relaxing the mass-independence requirement (as suggested by recent studies, e.g. Bongiorno et al. 2016; Georgakakis et al. 2017; Aird et al. 2018) can help to resolve these issues. As explained in Section 2.2, for this model, we assume that the quiescent component remains unchanged from the mass-independent model meaning that we only fit the star-forming X-ray luminosity functions at z < 2. These are extracted from the results of our mass-independent model, and is justified by the ability of our mass-independent model to reproduce the X-ray luminosity functions of star-forming and quiescent galaxies to z ∼ 1 from Georgakakis et al. (2014). 3.2.1 X-ray luminosity functions We show in Fig. 5, the fit to the star-forming component of the X-ray luminosity functions using our mass-dependent Eddington ratio distribution model to z = 1.75. Again, we do not extend our mass-dependent model beyond z ∼ 2 since it relies on the ability of our previous – mass-independent – model to split the X-ray luminosity functions in terms of star-forming and quiescent components (see Section 2.2). Out to z ∼ 2, our mass-dependent model does a good job at fitting the star-forming X-ray luminosity functions in all our redshift bins. Since our mass-dependent model assumes an λEdd distribution split into three different mass bins (see Section 2.2), we are also able to derive the X-ray luminosity function of each mass bin (see Fig. 5). In doing so, our model predicts that the highest mass galaxies (i.e. M* >1011 M⊙) dominate the star-forming X-ray luminosity functions across almost the full range of X-ray luminosities (i.e. 1041 <LX <1046 erg s−1). The exception being in our lowest redshift bin (i.e. z = 0.11), where 1010 <M*/M⊙ <1011 galaxies contribute marginally more than our highest mass galaxies at LX ≳ 1043.5 erg s−1. In general, however, these medium mass galaxies only contribute significantly (i.e. >10 per cent) at luminosities above the knee of the X-ray luminosity functions (i.e. LX >1044 erg s−1). By contrast, the lowest mass galaxies (i.e. M* <1010 M⊙) provide almost no contribution to the X-ray luminosity functions of star-forming galaxies at LX ≳ 1042.5 erg s−1 (i.e. <10 per cent), although we are only able to place upper limits on this contribution from lower mass galaxies. Figure 5. View largeDownload slide Fit to the star-forming component of the X-ray 2–10 keV luminosity functions using our mass-dependent Eddington ratio distribution model up to z = 1.75. Each panel shows a different redshift bin as indicated on their bottom right-hand side. The black line in each panel shows the star-forming component of the X-ray luminosity function (derived from our mass-independent model), and the orange line shows the best fit given by our mass-dependent model. The dashed, single-dot dashed, and triple-dot dashed lines show the X-ray luminosity functions for our highest (11 < log(M*/M⊙) < 12), medium (10 < log(M*/M⊙) < 11), and lowest (8 < log(M*/M⊙) < 10) mass bin, respectively. The black arrows indicate the overall effects of upper and lower limits on the parameters that define our mass-dependent Eddington ratio distribution. The blue stars and the red circles show the X-ray luminosity functions for the star-forming and quiescent galaxies, respectively, as reported in Georgakakis et al. (2014). The light grey downward and upward triangles show the total X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. Figure 5. View largeDownload slide Fit to the star-forming component of the X-ray 2–10 keV luminosity functions using our mass-dependent Eddington ratio distribution model up to z = 1.75. Each panel shows a different redshift bin as indicated on their bottom right-hand side. The black line in each panel shows the star-forming component of the X-ray luminosity function (derived from our mass-independent model), and the orange line shows the best fit given by our mass-dependent model. The dashed, single-dot dashed, and triple-dot dashed lines show the X-ray luminosity functions for our highest (11 < log(M*/M⊙) < 12), medium (10 < log(M*/M⊙) < 11), and lowest (8 < log(M*/M⊙) < 10) mass bin, respectively. The black arrows indicate the overall effects of upper and lower limits on the parameters that define our mass-dependent Eddington ratio distribution. The blue stars and the red circles show the X-ray luminosity functions for the star-forming and quiescent galaxies, respectively, as reported in Georgakakis et al. (2014). The light grey downward and upward triangles show the total X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. Our mass-dependent model suggests that AGNs hosted by star-forming galaxies with M* >1011 M⊙ dominate the X-ray luminosity functions. This is consistent with Georgakakis et al. (2017), which demonstrates that AGNs with LX >1041 − 42 erg s−1 significantly contribute to the AGN host stellar mass function for masses M* >1011 M⊙. This is also consistent with results from Bongiorno et al. (2016), which show that the host stellar mass function of AGNs with LX >1043 erg s−1 peaks at M* ∼1011 M⊙, out to z = 2. However, the vast majority of AGN host galaxies with LX <1043 erg s−1 display typical stellar masses in the range M* ∼1010 − 11 M⊙ (e.g. Bongiorno et al. 2016; Aird, Coil & Georgakakis 2017), which contrasts with our findings of a model X-ray luminosity function dominated by host galaxies with M* >1011 M⊙. However, our model is limited by the choice for the boundaries of our mass bins. Therefore, should we change our highest mass bin to incorporate galaxies with M* ∼1010.5 M⊙, we would expect to find that galaxies with M* >1010.5 M⊙ dominate the X-ray luminosity functions, in agreement with recent observational studies (e.g. Bongiorno et al. 2016; Aird et al. 2017; Wang et al. 2017). 3.2.2 SFR in bins of X-ray luminosities Having confirmed that our mass-dependent model is able to reproduce the X-ray luminosity functions for star-forming galaxies up to z ∼ 2, we now consider whether the corresponding PSM reproduces the observed flat relationship between SFR and X-ray luminosity (e.g. Rosario et al. 2012; Stanley et al. 2015). We show in Fig. 6, the mean-average SFR of AGNs split into 0.5 dex-wide bins of X-ray luminosity at similar redshifts as those of Stanley et al. (2015). Again, we include quiescent galaxies in our model when calculating these averages. Contrary to our mass-independent model, we find that our mass-dependent model predicts a flat relationship between SFR and X-ray luminosity in very good agreement with observations (see Fig. 6). It is important to stress that the model was not optimized to recreate the flat SFR/X-ray luminosity relationship. Therefore, we conclude that the mass-dependent model is able to reproduce the X-ray luminosity functions for star-forming galaxies (with a good agreement with observations at least up to z ∼ 1), while also independently reproducing the observed flat relationship between SFR and X-ray luminosity out to z ∼ 2. Figure 6. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosities predicted by our mass-dependent model (black thick line) at z = 0.35, z = 0.65, z = 1.15, and z = 2.0. The error bars represent the 1σ uncertainties on the mean for average SFR. The various symbols are from a compilation of published relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line shows the relationship between average SFR and X-ray luminosity found in Silverman et al. (2009), and the black dashed line represents the relationship for AGN-dominated systems, as reported in Netzer (2009). Although it is not implemented in the model, our mass-independent model successfully reproduces the observed flat SFR/X-ray relationship out to z = 2. Figure 6. View largeDownload slide The mean-average SFR in bins of X-ray 2–10 keV luminosities predicted by our mass-dependent model (black thick line) at z = 0.35, z = 0.65, z = 1.15, and z = 2.0. The error bars represent the 1σ uncertainties on the mean for average SFR. The various symbols are from a compilation of published relationships (see keys; i.e. Rosario et al. 2012, 2013; Stanley et al. 2015; Bernhard et al. 2016). The orange dashed line shows the relationship between average SFR and X-ray luminosity found in Silverman et al. (2009), and the black dashed line represents the relationship for AGN-dominated systems, as reported in Netzer (2009). Although it is not implemented in the model, our mass-independent model successfully reproduces the observed flat SFR/X-ray relationship out to z = 2. 3.2.3 The evolution of the λEdd distribution Having demonstrated that our mass-dependent model is able to reproduce both the X-ray luminosity functions for star-forming galaxies and the flat SFR/X-ray luminosity relationship to z ∼ 2, we now explore how the eight parameters (see Section 2.2) that define our mass-dependent model change with redshift. We show in Fig. 7, the redshift evolution of these parameters. We report that the overall normalization of the star-forming component (ASF) increases very slightly with increasing redshift. The position of the break ($$\lambda _{\rm break}^{\rm SF}$$) in each mass bin shifts towards higher λEdd values. This is consistent with the overall increase of the position of the break observed in our mass-independent model (i.e. all mass bins collapsed together; see Section 3.1.2), suggesting a higher typical accretion rate for AGNs at higher redshifts. Regarding the slopes at lower λEdd, ($$\gamma _1^{\rm SF}$$), we find that it does not evolve with redshift for any of our mass bins. However, although $$\gamma _1^{\rm SF}$$ for our highest mass bin is well constrained, we have large uncertainties associated with $$\gamma _1^{\rm SF}$$ for our medium mass bin, and upper limits for our lowest mass bin (see Fig. 7). These weak constraints on $$\gamma _1^{\rm SF}$$ for our medium and lowest mass bins are a direct consequence of their low contributions to the total Eddington ratio distribution for star-forming galaxies (see Fig. 8), as was hinted by their low contribution to the total X-ray luminosity functions (see Section 3.2). We also find that the shared slope at high λEdd ($$\gamma _2^{\rm SF}$$) does not change significantly with redshift. Overall, our findings suggest that there is a suppression of low-λEdd (i.e. λEdd ≲0.1) for AGNs hosted in lower mass galaxies (i.e. M* ≲ 1010 − 11 M⊙). Finally, we performed a linear fit of the redshift evolution of each parameter, using the standard deviation of 1000 Monte Carlo realisations for the uncertainties on the fitting parameters. We report the results of this fit in Table 2. Figure 7. View largeDownload slide Redshift evolution of the parameters that define our mass-dependent Eddington ratio distribution for star-forming galaxies at z = 0.1, z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.3, z = 1.5, and z = 1.7. Error bars show the 3σ uncertainties. The upward and downward arrows indicate the lower and upper limits, respectively. The two top panels present the shared parameters among our three mass bins, i.e. the normalization, $$\log _{A_{\rm SF}}$$, and the slope at high Eddington ratio, $$\gamma _2^{\rm SF}$$. The left-hand panels are the three different break positions $$\lambda _{\rm break}^{\rm SF}$$, one for each mass bin (as indicated in the bottom left-hand side in each panel), and the right-hand panels are the same but for that of the lower Eddington ratio slope, $$\gamma _1^{\rm SF}$$. The black lines show the best fit with redshift of each parameter. Figure 7. View largeDownload slide Redshift evolution of the parameters that define our mass-dependent Eddington ratio distribution for star-forming galaxies at z = 0.1, z = 0.3, z = 0.5, z = 0.7, z = 1.0, z = 1.3, z = 1.5, and z = 1.7. Error bars show the 3σ uncertainties. The upward and downward arrows indicate the lower and upper limits, respectively. The two top panels present the shared parameters among our three mass bins, i.e. the normalization, $$\log _{A_{\rm SF}}$$, and the slope at high Eddington ratio, $$\gamma _2^{\rm SF}$$. The left-hand panels are the three different break positions $$\lambda _{\rm break}^{\rm SF}$$, one for each mass bin (as indicated in the bottom left-hand side in each panel), and the right-hand panels are the same but for that of the lower Eddington ratio slope, $$\gamma _1^{\rm SF}$$. The black lines show the best fit with redshift of each parameter. Figure 8. View largeDownload slide The Eddington ratio probability distributions of star-forming and quiescent galaxies in our mass-dependent model at z = 0.2, z = 0.5, z = 1.0, z = 1.5, and z = 2.0. The left-hand panel shows both the star-forming (plain lines) and the quiescent components (single-dashed lines) of the Eddington ratio probability distribution. The normalization is such that the total Eddington ratio distribution at each redshift integrates to unity after applying an arbitrary cut at λEdd= 10−7. The right-hand panels are, for each redshift (indicated on their top right-hand corner), the contribution of the different mass bins to the Eddington ratio distributions of AGNs hosted in star-forming galaxies. In each of these panels, the dashed line is for the highest mass bin (11 < log(M*/M⊙) < 12), the single-dotted dashed line is for the medium mass bin (10 < log(M*/M⊙) < 11), and the triple-dotted dashed line is for the lowest mass bin (8 < log(M*/M⊙) < 10). The downward arrows for the slope at low-λEdd in our lower mass bin indicate that we have upper limits. The grey area illustrates the large uncertainties for the slope at low-λEdd in our medium mass bin (see Fig. 7). We normalized each of the total Eddington ratio distribution such that it integrates to unity, applying an arbitrary cut at λEdd = 10−7. Figure 8. View largeDownload slide The Eddington ratio probability distributions of star-forming and quiescent galaxies in our mass-dependent model at z = 0.2, z = 0.5, z = 1.0, z = 1.5, and z = 2.0. The left-hand panel shows both the star-forming (plain lines) and the quiescent components (single-dashed lines) of the Eddington ratio probability distribution. The normalization is such that the total Eddington ratio distribution at each redshift integrates to unity after applying an arbitrary cut at λEdd= 10−7. The right-hand panels are, for each redshift (indicated on their top right-hand corner), the contribution of the different mass bins to the Eddington ratio distributions of AGNs hosted in star-forming galaxies. In each of these panels, the dashed line is for the highest mass bin (11 < log(M*/M⊙) < 12), the single-dotted dashed line is for the medium mass bin (10 < log(M*/M⊙) < 11), and the triple-dotted dashed line is for the lowest mass bin (8 < log(M*/M⊙) < 10). The downward arrows for the slope at low-λEdd in our lower mass bin indicate that we have upper limits. The grey area illustrates the large uncertainties for the slope at low-λEdd in our medium mass bin (see Fig. 7). We normalized each of the total Eddington ratio distribution such that it integrates to unity, applying an arbitrary cut at λEdd = 10−7. Table 2. Redshift evolution of the parameters that describe the Eddington ratio distribution of star-forming galaxies for our mass-dependent model. The slopes and intercepts are given for an evolution as a function of (1+z). Uncertainties on the fitting parameters are extracted via 1000 Monte Carlo realisations (i.e. 1σ). Parameters  Intercepts  Slopes  log(ASF)  −1.57±0.08  0.03±0.04  log($$\lambda _{\rm break}^{\rm low\ mass}$$)  >−2.53  1.27  log($$\lambda _{\rm break}^{\rm medium\ mass}$$)  −2.05±0.15  0.68±0.08  log($$\lambda _{\rm break}^{\rm high\ mass}$$)  −2.60±0.08  0.70±0.04  $$\gamma _1^{\rm low\ mass}$$  <0.5  0.0  $$\gamma _1^{\rm medium\ mass}$$  −1.0±2.67  −0.12±1.30  $$\gamma _1^{\rm high\ mass}$$  0.45±0.05  −0.04±0.02  γ2  2.21±0.06  −0.03±0.03  Parameters  Intercepts  Slopes  log(ASF)  −1.57±0.08  0.03±0.04  log($$\lambda _{\rm break}^{\rm low\ mass}$$)  >−2.53  1.27  log($$\lambda _{\rm break}^{\rm medium\ mass}$$)  −2.05±0.15  0.68±0.08  log($$\lambda _{\rm break}^{\rm high\ mass}$$)  −2.60±0.08  0.70±0.04  $$\gamma _1^{\rm low\ mass}$$  <0.5  0.0  $$\gamma _1^{\rm medium\ mass}$$  −1.0±2.67  −0.12±1.30  $$\gamma _1^{\rm high\ mass}$$  0.45±0.05  −0.04±0.02  γ2  2.21±0.06  −0.03±0.03  Notes: Slopes and intercepts are given for an evolution as a function of (1+z). The intercept for z > 1.7 when assuming a break in the z evolution of the parameters is given by the continuity at z = 1.7 (i.e. [intercept forz > 1.7] = (1 + 1.7) × ([slope forz < 1.7] − [slope forz > 1.7]) + [intercept forz < 1.7]). View Large Each of the aforementioned trends can be seen in the evolution of the overall Eddington ratio distribution, which we plot out to z = 2 in Fig. 8. In this figure, we also plot the Eddington ratio distributions of quiescent host galaxies found in our mass-independent model at similar redshifts. In Fig. 8, for the left-hand panel, each distribution is normalize such that the total Eddington ratio distribution (i.e. the combination of the star-forming and quiescent component) integrates to unity after applying an arbitrary cut at λEdd= 10−7. For the right-hand panels, each total distribution is also normalized such that it integrates to unity after applying an arbitrary cut at λEdd= 10−7. 4 DISCUSSION 4.1 Further checks to our mass-dependent model As a further check of our mass-dependent model, we compared the model Eddington ratio distribution at z = 1 to the empirical one reported in Wang et al. (2017). As our various assumptions for this model are based on observations from Aird et al. (2018), comparing our mass-dependent Eddington ratio distribution to that empirical of Wang et al. (2017) – instead of Aird et al. (2018) – constitutes a more independent test for our model. We show this comparison in Fig. 9. We find very good agreement between our model Eddington ratio distributions to that of Wang et al. (2017) at z = 1 for both star-forming and quiescent component (at least for log(λEdd) < −0.5). This strongly supports our mass-dependent model for the Eddington ratio distribution of star-forming galaxies (since the mass-independent model was predicting a peaky Eddington ratio distribution for star-forming galaxies, which would contrasts with empirical results from Wang et al. 2017). We predict that the low-λEdd end of the Eddington ratio distribution for star-forming galaxies is dominated by galaxies with M* ≳ 1011 M⊙. However, Wang et al. (2017) reported that their sample of AGNs have typical stellar masses between 1010.5 <M*/M⊙ <1011 (once corrected for the differences in the IMFs between Salpeter 1955 for this work and Chabrier 2003 for Wang et al. 2017), which is half a dex below our highest mass bin. As mentioned in Section 3.2.1, this discrepancy in the host stellar masses between our model and observations is a consequence of the chosen boundaries for each of our mass bins. Figure 9. View largeDownload slide Comparison of our mass-dependent model Eddington ratio distribution to empirical results of Wang et al. (2017) at z = 1. The blue thick line is the total model Eddington ratio distribution for star-forming galaxies, while the dashed line, the dotted-dashed line, and the triple-dotted dashed line indicate our highest, medium, and lowest mass bin contributions, respectively. The blue circles are empirical data of the star-forming component from Wang et al. (2017) at z ∼ 1, while the red circles are that of the quiescent component. Figure 9. View largeDownload slide Comparison of our mass-dependent model Eddington ratio distribution to empirical results of Wang et al. (2017) at z = 1. The blue thick line is the total model Eddington ratio distribution for star-forming galaxies, while the dashed line, the dotted-dashed line, and the triple-dotted dashed line indicate our highest, medium, and lowest mass bin contributions, respectively. The blue circles are empirical data of the star-forming component from Wang et al. (2017) at z ∼ 1, while the red circles are that of the quiescent component. As a final check of our mass-dependent model, we show in Fig. 10 how the normalized average SFR (i.e. SFR relative to that of the MS) changes with λEdd in our mass-dependent model. Although it is not incorporated in the optimization, we predict a slight enhancement of average normalize SFR at higher λEdd (i.e. λEdd ≳ 1) compared to that of their lower λEdd counterpart (at least at z ≳ 1.2). This slight enhancement of normalized average SFR at higher λEdd is also observed in Bernhard et al. (2016). In Fig. 10, we only consider star-forming galaxies for our model since we do not have a good prescription for the normalized SFRs of quiescent galaxies. This is a possible reason for the discrepancy at lower λEdd between our model and the empirical results (i.e. that does contain quiescent galaxies). We also stress that in Fig. 10 our highest redshift bin, i.e. 1.8 < z < 2.9 constitutes an extrapolation of our model. However, the results are still consistent with observations (i.e. a slight enhancement of the normalized average SFR at higher λEdd). Figure 10. View largeDownload slide The normalized average SFR versus the Eddington ratio up to z ∼ 3. The blue lines show the prediction for our mass dependent model. the orange circles show the results report in Bernhard et al. (2016). The grey area represent the scatter around the MS as reported in Schreiber et al. (2015), the dashed black line being the MS. Figure 10. View largeDownload slide The normalized average SFR versus the Eddington ratio up to z ∼ 3. The blue lines show the prediction for our mass dependent model. the orange circles show the results report in Bernhard et al. (2016). The grey area represent the scatter around the MS as reported in Schreiber et al. (2015), the dashed black line being the MS. 4.2 Why our mass-independent model predicts a strong SFR/X-ray relationship One important result of our mass-independent model is that the Eddington ratio distribution for star-forming galaxies is better represented by a peaky distribution (i.e. a narrow distribution), rather than the more common broken power-law (i.e. a broader distribution). To explain this, we refer to Caplar et al. (2015; see also Weigel et al. 2017), where they demonstrate that the steepness of the faint-end slope of the X-ray luminosity function is determined by either the low-mass end of the mass function or the low-λEdd slope of the Eddington ratio distribution (since combined together to form the X-ray luminosity function) whichever is the steeper. Since that the faint-end of the observed X-ray luminosity function of Aird et al. (2015) flattens at higher redshifts whilst, in contrast, the low-mass end of the mass function for star-forming galaxies steepens with redshift (see Fig. 1), our mass-independent model attempts to reduce the contribution from star-forming galaxies at lower λEdd, leading the peaky Eddington ratio distribution for star-forming galaxies found by our mass-independent model. However, beyond z ∼ 2, the low-mass end of the mass function for star-forming galaxies of Davidzon et al. (2017) is already too steep compared to the flat faint-end of the X-ray luminosity function, such that our model is unable to reproduce the observed X-ray luminosity function. The steepening of the mass functions combined with the flattening of the X-ray luminosity functions suggests a mass-dependent Eddington ratio distribution for star-forming galaxies. We further saw that our mass-independent model failed to reproduce the observed flat SFR/X-ray relationship (see Fig. 4). To explain why this happens, we show in Fig. 11, which galaxy masses populate different parts of the SFR/X-ray plane. This plot shows that the mass-independent model produces a strong correlation between X-ray luminosity and stellar mass. This is a direct result of the narrower peaky (in contrast to a broader broken power law) Eddington ratio distribution for star-forming galaxies that this model demands in order to reproduce the observed X-ray luminosity functions, at least out to z ∼ 2. The narrowness of this Eddington ratio distribution means that a galaxy of a given mass can produce an AGN only within a limited range of luminosities. When we then include the correlation between SFR and stellar mass for MS galaxies (e.g. Schreiber et al. 2015), the consequence is a correlation between SFR and X-ray luminosity. This correlation between SFR and X-ray luminosity was also noted by Veale et al. (2014) when using their light-bulb Eddington ratio model, which is similar in shape to our peaky distribution for star-forming host galaxies. Figure 11. View largeDownload slide Stellar mass distributions (colour coded) in the SFR/X-ray luminosity parameter space, using our PSM and assuming our mass-independent model for the Eddington ratio distribution. Each panel shows a different redshift, as indicated on their bottom-right hand side. Each discrete colour indicates a stellar mass density in bins of log(M*)=0.5 (see colour-bar). Figure 11. View largeDownload slide Stellar mass distributions (colour coded) in the SFR/X-ray luminosity parameter space, using our PSM and assuming our mass-independent model for the Eddington ratio distribution. Each panel shows a different redshift, as indicated on their bottom-right hand side. Each discrete colour indicates a stellar mass density in bins of log(M*)=0.5 (see colour-bar). Overall, there is a large discrepancy with our mass-independent model that demands a narrow distribution for the Eddington ratio distribution of star-forming galaxies in order to reproduce the X-ray luminosity functions, but also requires a broader Eddington ratio distribution to reproduce the flat SFR/X-ray relationship. This strongly suggests that the Eddington ratio distribution must be somehow mass-dependent, such that low Eddington ratios are suppressed in low-mass galaxies to be able to reproduce the X-ray luminosity functions while still being broad to reproduce the flat relationship between SFR and X-ray luminosity. Our mass-dependent model provides further constrains on how the Eddington ratio distribution for star-forming galaxies should be mass-dependent. 4.3 Caveats to our models While deriving the λEdd distribution, our model is limited by the functional form chosen for each component of the λEdd distribution and the associated assumptions used to avoid degeneracies during fitting (see Section 2). An obvious caveat to our model is that we cannot explore the full range of possible Eddington ratio distributions. However, we find that, should we relax the assumptions made prior to the fit of the X-ray luminosity functions, our solutions are too degenerate to form a set of meaningful parameters for the Eddington ratio distributions, meaning that the X-ray luminosity functions lack the information needed to constrain the Eddington ratio distribution and its mass-dependency. Using these assumptions is, therefore, a requirement in our modelling approach. As stressed in Section 2, all our assumptions are motivated by the findings of recent studies (e.g. Jones et al. 2017; Aird et al. 2018), and the results predicted by the mass-dependent model are in good agreement with many recent observational studies, i.e. the observed X-ray luminosity functions out to z ∼ 2 – including when split in terms of star-forming and quiescent galaxies (e.g. Georgakakis et al. 2014; Aird et al. 2015), the flat relationship between SFR and X-ray luminosity (e.g. Rosario et al. 2012; Stanley et al. 2015), the empirical Eddington ratio distribution split between star-forming and quiescent component at z = 1 (Wang et al. 2017), and the enhancement of the star-forming properties at higher λEdd (Bernhard et al. 2016). As such while we do not claim that our Eddington ratio distributions are universal, they provide a simple means by which to explore the AGN-galaxy connection. One of the main caveats for our mass-dependent model is that we are limited by the choice of the stellar mass bins. As a consequence, we find that our dominant component for the Eddington ratio distribution and the X-ray luminosity function is our highest mass bin with M* >1011 M⊙. This is qualitatively consistent with observations (e.g. Georgakakis et al. 2014; Aird et al. 2015; Wang et al. 2017), yet roughly half a dex above the observed typical host mass for AGNs. Investigating this in more detail would require additional mass bins. However, this would also generates a larger number of free parameters, and therefore degeneracies. We have chosen these mass bins to be broad enough to cover a wide range of masses, yet probing any mass dependency of the Eddington ratio distribution for AGNs in star-forming hosts. 4.4 Extending our mass-dependent model to higher redshifts We recall that we fit the mass-dependent model to the X-ray luminosity functions of star-forming galaxies only, which was isolated using our mass-independent model. The results of our mass-dependent model therefore depends on the reliability of our mass-independent model to separate the X-ray luminosity functions into star-forming and quiescent components (see Section 2). Thankfully, we can verify this up to z ∼ 1 using the results reported in Georgakakis et al. (2014). However, beyond z ∼ 1, we do not have the empirical results to perform this check, while at z ≳ 2 our mass-independent model is actually inconsistent with observations (i.e. too many AGNs in quiescent galaxies). Faced with this situation where we cannot reliably exploit the results of our mass-independent model, we now investigate whether we can fit the observed X-ray luminosity functions at z = 2.25 using both a mass-dependent star-forming component and a mass-independent quiescent component λEdd distributions. This results in a total of 12 free parameters to optimize (i.e. eight parameters for the star-forming component and four parameters for the quiescent one). We use the same assumptions for the λEdd distribution of star-forming galaxies as described in Section 2, and demand that the normalization of the Eddington ratio distribution for quiescent galaxies lies below that of star-forming galaxies (i.e. consistent with the lower redshift bins). We show in Fig. 12, the results of this new model. Contrary to our mass-independent model, we are now able to reproduce the total observed X-ray luminosity functions at z = 2.25 by placing the majority of AGNs in star-forming galaxies. In particular, we find that the X-ray luminosity functions at z = 2.25 is dominated by the highest mass star-forming galaxies, with a smaller contribution from medium mass star-forming galaxies. The contribution from our lowest mass bin is consistent with zero. Similarly, we can only place upper limits on the contribution from quiescent galaxies to the total X-ray luminosity functions at z = 2.25. We conclude, therefore, that while we are able to reproduce the total X-ray luminosity functions using a mass-dependent λEdd distribution for star-forming galaxies (and ensuring that it is dominated by star-forming galaxies), this solution is degenerate with the level of contribution from quiescent galaxies. To break this degeneracy will require the separation of the high-redshift X-ray luminosity functions into quiescent and star-forming components (i.e. as performed by Georgakakis et al. (2014) at z < 1). Figure 12. View largeDownload slide Fit of the X-ray luminosity functions at z = 2.25 of Aird et al. (2015) for our model assuming a mass-dependent Eddington ratio distribution for star-forming galaxies and a mass-independent one for that of quiescent galaxies. The black downward and upward triangles show the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The orange line shows the total X-ray luminosity functions derived using our model, the blue lines are for that of the star-forming component in each different mass bin (see keys), and the quiescent component is shown with a red line. Downward black arrows indicate upper limits. Figure 12. View largeDownload slide Fit of the X-ray luminosity functions at z = 2.25 of Aird et al. (2015) for our model assuming a mass-dependent Eddington ratio distribution for star-forming galaxies and a mass-independent one for that of quiescent galaxies. The black downward and upward triangles show the X-ray luminosity functions of Aird et al. (2015) for the soft and the hard band, respectively. The orange line shows the total X-ray luminosity functions derived using our model, the blue lines are for that of the star-forming component in each different mass bin (see keys), and the quiescent component is shown with a red line. Downward black arrows indicate upper limits. 5 CONCLUSION Motivated by recent results reporting a different Eddington ratio distribution for star-forming and quiescent galaxies (Wang et al. 2017; Aird et al. 2018), we attempt to constrain these distributions by using an analytical model to fit the observed X-ray luminosity functions of Aird et al. (2015), assuming the galaxy mass functions of Davidzon et al. (2017). Our first model assumes two mass-independent λEdd distributions: one for star-forming galaxies and another for quiescent galaxies. After optimization, we find that this model is able to reproduce the X-ray luminosity functions out to z ∼ 2 (see Section 3.1.1), but demands a peaky distribution for the Eddington ratio distribution of star-forming galaxies. Despite this, our mass-independent model fails to reproduce the observed flat SFR/X-ray relationship when incorporated into a PSM for galaxies (see Section 3.1.3). We argue that this failure arises because a mass-independent model places too many low-luminosity AGNs in low-mass galaxies (see Section 4). To resolve this problem, we develop a second model in which the Eddington ratio distribution for star-forming galaxies is allowed to be different in three mass bins (motivated by observation from e.g. Georgakakis et al. 2017; Aird et al. 2018). By suppressing low Eddington ratio AGNs in low-mass galaxies, this mass-dependent model is able to reproduce the X-ray luminosity functions for star-forming galaxies while simultaneously reproducing the observed flat relationship between SFR and X-ray luminosity (see Section 3.2.2). Finally, we also find that the mass-dependent model is consistent with empirical Eddington ratio distribution for star-forming and quiescent galaxies at z = 1, and is able to reproduce the slight enhancement of normalized average SFR at higher λEdd, as reported in Bernhard et al. (2016). Overall, we conclude that, in our model, a suppression of lower AGN accretion activity in low mass galaxies is required to reproduce both the observed X-ray luminosity functions and the observed flat SFR/X-ray relationship. Acknowledgements We thank the referee for the useful comments that help to significantly improve the quality of the paper. EB thanks C. Tadhunter and D. Alexander for the helpful discussions. EB thanks the University of Sheffield for the support grant. Footnotes 1 To a constant factor given in equation (1). 2 emcee is available on-line at http://dan.iel.fm/emcee/current/. 3 Hereafter, we refer to the shape of our model Eddington ratio distribution of AGNs hosted by star-forming galaxies as peaky. REFERENCES Aird J. et al.  , 2010, MNRAS , 401, 2531 https://doi.org/10.1111/j.1365-2966.2009.15829.x CrossRef Search ADS   Aird J. et al.  , 2013, ApJ , 775, 41 https://doi.org/10.1088/0004-637X/775/1/41 CrossRef Search ADS   Aird J., Coil A. L., Georgakakis A., Nandra K., Barro G., Pérez-González P. G., 2015, MNRAS , 451, 1892 https://doi.org/10.1093/mnras/stv1062 CrossRef Search ADS   Aird J., Coil A. L., Georgakakis A., 2017, MNRAS , 465, 3390 https://doi.org/10.1093/mnras/stw2932 CrossRef Search ADS   Aird J. Coil A. L. Georgakakis A., 2018, MNRAS , 474, 1225 CrossRef Search ADS   Assef R. J. et al.  , 2011, ApJ , 728, 56 https://doi.org/10.1088/0004-637X/728/1/56 CrossRef Search ADS   Azadi M. et al.  , 2015, ApJ , 806, 187 https://doi.org/10.1088/0004-637X/806/2/187 CrossRef Search ADS   Bernhard E., Béthermin M., Sargent M., Buat V., Mullaney J. R., Pannella M., Heinis S., Daddi E., 2014, MNRAS , 442, 509 https://doi.org/10.1093/mnras/stu896 CrossRef Search ADS   Bernhard E., Mullaney J. R., Daddi E., Ciesla L., Schreiber C., 2016, MNRAS , 460, 902 https://doi.org/10.1093/mnras/stw973 CrossRef Search ADS   Bongiorno A. et al.  , 2016, A&A , 588, A78 CrossRef Search ADS   Caplar N., Lilly S. J., Trakhtenbrot B., 2015, ApJ , 811, 148 https://doi.org/10.1088/0004-637X/811/2/148 CrossRef Search ADS   Chabrier G., 2003, PASP , 115, 763 https://doi.org/10.1086/376392 CrossRef Search ADS   Conroy C., White M., 2013, ApJ , 762, 70 https://doi.org/10.1088/0004-637X/762/2/70 CrossRef Search ADS   Daddi E. et al.  , 2007, ApJ , 670, 156 https://doi.org/10.1086/521818 CrossRef Search ADS   Davidzon I. et al.  , 2017, A&A , 605, A70 CrossRef Search ADS   Elbaz D. et al.  , 2011, A&A , 533, A119 CrossRef Search ADS   Fabian A. C., 2012, ARA&A , 50, 455 CrossRef Search ADS   Foreman-Mackey D., Hogg D. W., Lang D., Goodman J., 2013, PASP , 125, 306 https://doi.org/10.1086/670067 CrossRef Search ADS   Gelman A., Rubin D. B., 1992, Stat. sci. , pp 457– 472 Georgakakis A. et al.  , 2014, MNRAS , 443, 3327 https://doi.org/10.1093/mnras/stu1326 CrossRef Search ADS   Georgakakis A., Aird J., Schulze A., Dwelly T., Salvato M., Nandra K., Merloni A., Schneider D. P., 2017, MNRAS , 471, 1976 https://doi.org/10.1093/mnras/stx1602 CrossRef Search ADS   Goodman J., Weare J., 2010, Commun. Appl. Math. Comput. sci. , 5, 65 https://doi.org/10.2140/camcos.2010.5.65 CrossRef Search ADS   Harrison C. M. et al.  , 2012, ApJ , 760, L15 https://doi.org/10.1088/2041-8205/760/1/L15 CrossRef Search ADS   Heckman T. M., Kauffmann G., Brinchmann J., Charlot S., Tremonti C., White S. D. M., 2004, ApJ , 613, 109 https://doi.org/10.1086/422872 CrossRef Search ADS   Hickox R. C., Mullaney J. R., Alexander D. M., Chen C.-T. J., Civano F. M., Goulding A. D., Hainline K. N., 2014, ApJ , 782, 9 https://doi.org/10.1088/0004-637X/782/1/9 CrossRef Search ADS   Ilbert O. et al.  , 2013, A&A , 556, A55 CrossRef Search ADS   Jones M. L., Hickox R. C., Mutch S. J., Croton D. J., Ptak A. F., DiPompeo M. A., 2017, ApJ , 843, 125 https://doi.org/10.3847/1538-4357/aa7632 CrossRef Search ADS   Kormendy J., Richstone D., 1995, ARA&A , 33, 581 CrossRef Search ADS   Kormendy J., Bender R., Cornell M. E., 2011, Nature , 469, 374 https://doi.org/10.1038/nature09694 CrossRef Search ADS PubMed  Lanzuisi G. et al.  , 2017, A&A , 602, A123 CrossRef Search ADS   Larson D. et al.  , 2011, ApJS , 192, 16 https://doi.org/10.1088/0067-0049/192/2/16 CrossRef Search ADS   Lutz D. et al.  , 2010, ApJ , 712, 1287 https://doi.org/10.1088/0004-637X/712/2/1287 CrossRef Search ADS   Magorrian J. et al.  , 1998, AJ , 115, 2285 https://doi.org/10.1086/300353 CrossRef Search ADS   Marconi A., Hunt L. K., 2003, ApJ , 589, L21 https://doi.org/10.1086/375804 CrossRef Search ADS   Merloni A., Rudnick G., Di Matteo T., 2004, MNRAS , 354, L37 https://doi.org/10.1111/j.1365-2966.2004.08382.x CrossRef Search ADS   Merritt D., 2000, in Combes F. Mamon G. A. Charmandaris V., eds, ASPacific Conf. Ser. Vol. 197, Dynamics of Galaxies: from the Early Universe to the Present . Astron. Soc. Pac., San Francisco, p. 221 Mullaney J. R. et al.  , 2012, ApJ , 753, L30 https://doi.org/10.1088/2041-8205/753/2/L30 CrossRef Search ADS   Mullaney J. R. et al.  , 2015, MNRAS , 453, L83 https://doi.org/10.1093/mnrasl/slv110 CrossRef Search ADS   Netzer H., 2009, MNRAS , 399, 1907 https://doi.org/10.1111/j.1365-2966.2009.15434.x CrossRef Search ADS   Netzer H., Lani C., Nordon R., Trakhtenbrot B., Lira P., Shemmer O., 2016, ApJ , 819, 123 https://doi.org/10.3847/0004-637X/819/2/123 CrossRef Search ADS   Rodighiero G. et al.  , 2011, ApJ , 739, L40 https://doi.org/10.1088/2041-8205/739/2/L40 CrossRef Search ADS   Rosario D. J. et al.  , 2012, A&A , 545, A45 CrossRef Search ADS   Rosario D. J. et al.  , 2013, A&A , 560, A72 CrossRef Search ADS   Rovilos E. et al.  , 2012, A&A , 546, A58 CrossRef Search ADS   Salim S. et al.  , 2007, ApJS , 173, 267 https://doi.org/10.1086/519218 CrossRef Search ADS   Salpeter E. E., 1955, ApJ , 121, 161 https://doi.org/10.1086/145971 CrossRef Search ADS   Santini P. et al.  , 2012, A&A , 540, A109 CrossRef Search ADS   Sargent M. T., Béthermin M., Daddi E., Elbaz D., 2012, ApJ , 747, L31 https://doi.org/10.1088/2041-8205/747/2/L31 CrossRef Search ADS   Schechter P., 1976, ApJ , 203, 297 https://doi.org/10.1086/154079 CrossRef Search ADS   Schreiber C. et al.  , 2015, A&A , 575, A74 CrossRef Search ADS   Silverman J. D. et al.  , 2008, ApJ , 679, 118 https://doi.org/10.1086/529572 CrossRef Search ADS   Silverman J. D. et al.  , 2009, ApJ , 696, 396 https://doi.org/10.1088/0004-637X/696/1/396 CrossRef Search ADS   Stanley F., Harrison C. M., Alexander D. M., Swinbank A. M., Aird J. A., Del Moro A., Hickox R. C., Mullaney J. R., 2015, MNRAS , 453, 591 https://doi.org/10.1093/mnras/stv1678 CrossRef Search ADS   Stanley F. et al.  , 2017, MNRAS , 742, 2221 CrossRef Search ADS   Vasudevan R. V., Fabian A. C., 2007, MNRAS , 381, 1235 https://doi.org/10.1111/j.1365-2966.2007.12328.x CrossRef Search ADS   Veale M., White M., Conroy C., 2014, MNRAS , 445, 1144 https://doi.org/10.1093/mnras/stu1821 CrossRef Search ADS   Wang T. et al.  , 2017, A&A , 601, A63 CrossRef Search ADS   Weigel A. K. et al.  , 2017, ApJ , 845, 145 https://doi.org/10.3847/1538-4357/aa8097 CrossRef Search ADS   © 2018 The Author(s) Published by Oxford University Press on behalf of the Royal Astronomical Society

### Journal

Monthly Notices of the Royal Astronomical SocietyOxford University Press

Published: May 1, 2018

## You’re reading a free preview. Subscribe to read the entire article.

### DeepDyve is your personal research library

It’s your single place to instantly
that matters to you.

over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month ### Explore the DeepDyve Library ### Search Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly ### Organize Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place. ### Access Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals. ### Your journals are on DeepDyve Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more. All the latest content is available, no embargo periods. DeepDyve ### Freelancer DeepDyve ### Pro Price FREE$49/month
\$360/year

Save searches from
PubMed

Create lists to

Export lists, citations