ABSTRACT

JWST is unveiling for the first time accreting black holes (BHs) with masses of 106 − 107 M at z > 4, with the most distant residing in GNz11 at z = 10.6. Are we really surprised to find them in the nuclei of z ≃ 5 − 11 galaxies? Here we predict the properties of 4 < z < 11 BHs and their host galaxies considering an Eddington-limited (EL) and a super-Eddington (SE) BH accretion scenario, using the Cosmic Archaeology Tool (CAT) semi-analytical model. We calculate the transmitted spectral energy distribution of CAT synthetic candidates, representative of the BH/galaxy properties of GNz11. We also examine the possibility that the z = 8.7 galaxy CEERS-1019 could host an active BH. We find that the luminosity of high-z JWST detected BHs are better reproduced by the SE model, where BHs descend from efficiently growing light and heavy seeds. Conversely, the host galaxy stellar masses are better matched in the EL model, in which all the systems detectable with JWST surveys JADES and CEERS descend from heavy BH seeds. We support the interpretation that the central point source of GNz11 could be powered by a SE (λEdd ≃ 2 − 3) accreting BH with mass 1.5 × 106 M, while the emission from CEERS-1019 is dominated by the host galaxy; if it harbours an active BH, we find it to have a mass of MBH ≃ 107 M, and to be accreting at sub-Eddington rates (λEdd ≃ 0.5).

1 INTRODUCTION

Understanding how galaxies form and evolve in the first billion years of cosmic evolution is one of the driving science objectives of the recently launched JWST. With its unprecedented sensitivity, JWST has already started to revolutionize the field. JWST/NIRCam observations have enabled to detecting galaxies at z ≳ 10 (Castellano et al. 2022; Labbe et al. 2023a; Tacchella et al. 2023a; Naidu et al. 2022b; Adams et al. 2023; Robertson et al. 2023; Tacchella et al. 2023b), providing the first observational constraints on the UV luminosity function and star formation history in the first 500 Myr of the universe (Harikane et al. 2022; Donnan et al. 2023), showing a surprising lack of evolution in the number density of bright galaxies at z ≳ 10. While some of these sources must be interpreted with caution, given their possible confusion with dusty, star-forming galaxies at z ∼ 5 (Naidu et al. 2022a; Arrabal Haro et al. 2023; Zavala et al. 2023), JWST is showing that galaxies in the early universe are characterized by a diversity of properties (Barrufet et al. 2023; Whitler et al. 2023), including unexpectedly strong dust obscuration at 8 < z < 13 (Rodighiero et al. 2023), massive quiescent galaxies at 3 < z < 5 (Carnall et al. 2023), and even low-mass quenched galaxies out to z ≃ 7 (Looser et al. 2023). JWST/NIRSpec observations have spectroscopically confirmed several tens of these systems (see Nakajima et al. 2023 and references therein), including galaxies at 10 ≤ z ≤ 13.2 (Curtis-Lake et al. 2023; Arrabal Haro et al. 2023; Bunker et al. 2023; Robertson et al. 2023). Most of these galaxies show strong nebular lines, indicative of extreme excitation conditions, with high-ionization parameters and low-metallicities (Cameron et al. 2023; Curti et al. 2023; Sanders et al. 2023a,b).

Interestingly, the spectra of some of these galaxies reveal the presence of broad permitted lines (H α, H β, Hγ) without a counterpart in the forbidden lines, hence consistent with being emitted from the Broad Line Region (BLR) of weak Active Galactic Nuclei (AGNs). Two systems with estimated BH masses of MBH ∼ 107 M have been reported by Kocevski et al. (2023): CEERS-1670 at z = 5.2, and CEERS-3210 at z = 5.6, the latter being a heavily obscured AGN caught in a transition phase between a dust-obscured starburst and an unobscured quasar, similar to the transitioning quasar GNz7q at z ∼ 7.19 reported by Fujimoto et al. (2022). Übler et al. (2023) identified an AGN in the compact galaxy GS-3073 at z = 5.6, estimating a BH mass of ∼108 M. Larson et al. (2023) reported a broad H β emission line in CEERS-1019, a galaxy at z = 8.7, suggesting a BH mass of MBH ∼ 107 M, although their detection has a lower significance (2.3σ). At even higher redshift, Maiolino et al. (2023a) argued that the z = 10.6 galaxy GNz-11 (Bunker et al. 2023; Tacchella et al. 2023b) hosts a nuclear BH with mass MBH ∼ 2 × 106 M accreting at super-Eddington (SE) rates. The new extensive analysis of JWST/NIRSpec observations of GNz11 presented in Maiolino et al. (2023a) suggests extremely high gas densities, typical of the BLR of AGNs, together with other AGN-like features (ionization lines typical of AGNs, ionization cones, and a Ly α halo similar to that observed in quasars, see Scholtz et al. 2023; Maiolino et al. 2023b), while Wolf–Rayet (WR) stars alone (Senchyna et al. 2023) can not account for many of the observed spectral properties.

In addition to GNz11, Maiolino et al. (in preparation) have identified 10 AGNs at 4 < z < 8 in three tiers of the JADES survey through the detection of broad Balmer emission lines, with also possible candidate merging BH systems, as traced by two broad H α components. The inferred BH masses in these AGN sample are in the range between 5 × 107 M and down to 3 × 105 M, probing the regime expected for heavy BH seeds formed by the collapse of Super Massive Stars (SMSs) (Ferrara et al. 2014).

A statistical census of broad H α emission line systems among 185 galaxies at 3.8 ≤ z ≤ 8.9 confirmed with NIRSpec has led to the discovery of 10 type-1 AGNs (corresponding to |$\sim 5~{{\ \rm per\ cent}}$| of the sample) with masses MBH ∼ 106 − 107 M at z = 4 − 6.9 (Harikane et al. 2023), two of which are the systems reported by Kocevski et al. (2023). It is important to note that broad emission lines can also be powered by mergers or outflows, but in these cases broad lines are seen both in permitted and forbidden lines, unlike these JWST sources, whose spectra show the typical narrow–forbidden/broad–permitted duality of AGNs.

Additional high redshift faint AGN candidates have been identified by means of high-resolution imaging and spectral energy distribution (SED) modelling (Furtak et al. 2023; Endsley et al. 2023; Labbe et al. 2023b; Onoue et al. 2023), including a system reported by Bogdan et al. (2023); this is a luminous X-ray source in the gravitationally lensed galaxy UHZ1, which has a solid photometric redshift of z ≃ 10.3 (Castellano et al. 2022). They interpret the source as a 107 − 108 M accreting BH, with a bolometric luminosity of Lbol = [2 − 8] × 1045 erg s–1. Interestingly, the estimated stellar mass of the galaxy is |$M_{\rm star} = 0.4^{+1.9}_{-0.2} \times 10^8 M_\odot$| (Castellano et al. 2022).1, suggesting that system maybe experiencing the Outsize Black Hole Galaxy (OBG) stage predicted by theoretical models for high redshift galaxies which are the birth places of heavy BH seeds (Agarwal et al. 2013; Valiante et al. 2016; Natarajan et al. 2017a).

These candidate AGNs, are often selected with red colours at long wavelengths, interpreted as reddened AGN emission, and a bluer component a short wavelength, possibly associated with scattered light from the AGN or with faint unobscured star formation in the host (Furtak et al. 2022; Labbe et al. 2023b). However, the possibility of compact dusty star formation can not be ruled out, and it would be essential to unambiguously identify the presence of a faint AGN through broad wing Balmer emission lines.

Recently, Matthee et al. (2023) have conducted a systematic search for broad line H α emitters over |$\simeq 230 \, {\rm arcmin}^2$| covered by the EIGER (Kashino et al. 2023) and FRESCO (Oesch et al. 2023) JWST programmes, corresponding to a volume of 6 × 105 Mpc3 over 4 ≤ z ≤ 6. They reported the identification of 20 systems at z = 4.2 − 5.5 that have broad H α components with line widths from 1200 − 3700 km s–1, that they interpreted as being powered by accretion onto SMBHs with masses ∼107 − 108 M, measuring a number densities of ≃ 10−5 cMpc−3 in the UV luminosity range −21 ≤ MUV ≤ −18.

Given the high rate of faint AGNs discovered with JWST, two questions naturally arise:

are we surprised to find these hidden little monsters2 in the hearth of galaxies at z ∼ 5 − 11? And can we use their estimated properties to constrain the nature and growth of the first BH seeds?

2 GREAT EXPECTATIONS AND AMAZING CONFIRMATIONS

The recent discovery of more than 40 new faint AGNs (with luminosities <1046 erg s–1) at z > 4 in ongoing surveys with JWST, including two tentative systems at z = 8.7 and 10.3, and the currently most distant AGN at z = 10.6, is an amazing confirmation of the rich BH landscape at cosmic dawn predicted by theoretical models.

The existence of hundreds of quasars at 6 ≤ z ≤ 7.64 (see Fan, Banados & Simcoe 2022 for a recent review, and Bañados et al. 2018; Yang et al. 2020; and Wang et al. 2021 for the three systems that mark the current quasars redshift frontier at z = 7.54, 7.52, and 7.64, respectively) show that 109 − 1010 M super massive black holes (SMBHs) are able to form in ∼700 Myr, by efficiently growing mass onto smaller BH seeds (see Volonteri, Habouzit & Colpi 2021 and Inayoshi, Visbal & Haiman 2020 for recent reviews).

The nature of BH seeds is still largely unconstrained. Astrophysical BH seeds have been proposed to form in different flavours: light BH seeds with masses up to a few |$\sim 100 \, M_\odot$| may be the remnants of massive Population III (Pop III) stars; heavy BH seeds with masses of ∼104 − 105 M may form by the direct collapse of a SMS; medium weight BH seeds, with masses ∼103 M may result from runaway stellar or BH collisions in dense stellar clusters. The physical conditions enabling the formation of BH seeds are tightly linked to gas cooling, heating, and accretion in metal-poor environments, and their occurrence goes hand-in-hand with the formation of early cosmic structures.

Theoretical models which attempt to incorporate a self-consistent BH seeding prescription find that the genealogy of SMBHs with masses >109 M powering quasars at z ∼ 6 is characterized by a rich variety of BH seed progenitors (see e.g. Volonteri & Natarajan 2009; Valiante et al. 2016; Ricarte & Natarajan 2018; Spinoso et al. 2023). Only a small fraction (<10 − 20 per cent) of BH seeds efficiently grow and contribute to the mass of the most massive SMBHs at z > 6 (Sassano et al. 2021), while the vast majority experience a more quiet evolution and populate the low end of the mass and luminosity function, which potentially encode information on their growth mode (Trinca et al. 2022).

2.1 Expected number of JWST detectable AGNs and challenges to their detection

Within a scenario where SMBHs at z ∼ 6 originate from the Eddington-limited (EL) growth of light and heavy BH seeds, Valiante et al. (2018) find that efficiently growing heavy BH seed progenitors are within the sensitivity limit of NIRCam3 out to z ∼ 16. The number of detectable seeds depends on the survey area and on the degree of obscuration. Trinca et al. (2023b) estimate that, in a survey like CEERS, between 48 to 175 AGNs with BH masses in the range MBH = [106 − 108] M would be observable at 5 ≤ z < 7, and between 8 and 21 at 7 ≤ z ≤ 10. JADES-Deep is expected to find between 12 and 63 AGNs with with BH masses in the range MBH = [104 − 106] M at 7 ≤ z ≤ 10, and between 5 and 32 at z ≳ 10. These numbers are predicted to be even larger in models accounting for phases of SE accretion (Lupi et al. 2014; Madau, Haardt & Dotti 2014; Volonteri, Silk & Dubus 2015; Pezzulli, Valiante & Schneider 2016), although in this regime BH growth is very intermittent, as suggested by high-resolution hydrodynamical simulations (Regan et al. 2019; Massonneau et al. 2023; Sassano et al. 2023).

Although these numbers are still significantly larger than the confirmed AGNs (Harikane et al. 2023; Larson et al. 2023; Matthee et al. 2023; Maiolino et al. 2023a), spectroscopy in these fields is still incomplete and does not yet match the depth of imaging. In addition, the theoretical predictions reported earlier have been estimated only on the basis of the flux limits and area coverage of ongoing JWST surveys, without estimating the contamination of the host galaxy and how this hampers the AGN identification. A common feature of faint AGNs discovered by JWST is that their identification relied on the detection of broad hydrogen (and in some cases helium) permitted lines, while their colours are hardly distinguishable from a metal-poor star-forming galaxies, and often dominated by the stellar continuum of their host galaxies. Prior to the spectroscopic confirmation by Kocevski et al. (2023), CEERS-1670 was suggested to be either a strong broad-line emitter or a SE accreting BH by Onoue et al. (2023) (the source was named CEERS-AGN-z5-1 in this study), on the basis of spectral fitting of deep NIRCam photometry (complemented by additional Canada–France–Hawaii Telescope, Hubble Space Telescope (HST), and Spitzer Space Telescope photometric data) and its point-source morphology; however, the possibility of an extremely young galaxy with moderate dust attenuation was not completely ruled out. In addition, Harikane et al. (2023) report that high-resolution HST and JWST/NIRcam images of 7 out of 10 faint AGNs identified as broad H α emission line systems show extended morphologies, indicating a significant contribution to the total light from their host galaxies. Hence, commonly adopted colour–colour selection techniques that have been developed (Volonteri et al. 2017; Natarajan et al. 2017b; Barrow, Aykutalp & Wise 2018; Valiante et al. 2018; Goulding & Greene 2022; Trussler et al. 2023) would not be able to identify these systems, unless the BH is overmassive and the AGN provides a dominant contribution to the total emission of the galaxy (Volonteri, Habouzit & Colpi 2023). In addition, Kocevski et al. (2023); Larson et al. (2023); Übler et al. (2023), and Harikane et al. (2023) show that these JWST discovered AGNs would not be classified as such using standard spectral diagnostic techniques based on narrow line ratios that are commonly adopted to discriminate galaxies dominated by the AGN emission or by stellar emission, such as the BPT diagram (Baldwin, Phillips & Terlevich 1981) or the OHNO diagram (Backhaus et al. 2022). In fact, in the BPT diagram they occupy the region with high [OIII]/H β and low [NII]/H α, similar to other high-z galaxies without broad line emission (Nakajima et al. 2023; Sanders et al. 2023a) and within the z ∼ 0 classical star-forming regions. Similarly, in the OHNO diagram, the systems lie in the AGN region of the diagram but overlap with other high-z star-forming galaxies (Trump et al. 2023). The observed line ratios of all these systems are consistent with what predicted by photoionization models with low-metallicity and high-ionization parameters, and alternative diagnostic diagrams have been proposed to separate stellar emission from AGN emission in these early galaxies, such as the EW(HeIIλ4686) versus HeIIλ4686/H β diagram (Nakajima & Maiolino 2022), where the source GS-3073 is compatible with an AGN and clearly separated from star-forming galaxies (Übler et al. 2023).

What are the physical properties of JWST detectable BH seeds, and can we identify – among the simulated samples – systems that could represent synthetic counterparts of the currently most distant AGN, GNz11 at z = 10.6 (Maiolino et al. 2023a), and of CEERS-1019 at z = 8.7, a system that may potentially host an AGN (Larson et al. 2023)?

2.2 Expected properties of JWST detectable AGNs

We base our investigation on the Cosmic Archaeology Tool (CAT) model presented in Trinca et al. (2022) (hereafter T22). CAT is a semi-analytical model that allows us to follow the co-evolution between nuclear BHs and their host galaxies all the way from the formation of the first stars and their BH remnants in the star-forming molecular-cooling halos at z ≳ 20 − 30, down to the assembly of systems hosting bright quasars and less luminous AGNs at 4 ≤ z ≲ 7.5. The model runs on a sample of dark matter (DM) halo merger histories at 4 ≤ z ≤ 30 generated through the galaxy formation model galform (Cole et al. 2000; Parkinson, Cole & Helly 2008) based on the Extended Press Schechter theory (EPS)4 We adopt a redshift dependent mass resolution that corresponds to a DM halo virial temperature of 1200 K, and we simulate the formation histories of |$[10^9-10^{14}]\, M_\odot$| DM halos at z = 4 by sampling their mass interval into 11 logarithmically spaced bins and by generating – for each bin – 10 independent merger trees.

In each halo, we follow gas infall, cooling, and star formation following the same implementation described in T22. Pop III stars form in metal poor star-forming regions with metallicity Z < Zcr = 10−3.8Z, with a mass that is stochastically sampled from a Larson initial mass function (IMF) with a characteristic mass of mch = 20 M in the mass range 10 Mm* ≤ 300 M. When ZZcr, Pop II/I stars form according to a Larson IMF with a characteristic mass of mch = 0.35 M in the mass range 0.1 Mm* ≤ 100 M. Stars are evolved according to their characteristic evolutionary time-scales, and we account for their radiative, chemical, and mechanical feedback by computing their metallicity- and age-dependent radiative output (in the hydrogen ionizing and Lyman–Werner bands), their mass- and metallicity-dependent metal/dust yields, and the energy deposited by core-collapse and pair-instability supernova explosions. This allows us to follow the gas, metal, and dust content in each galaxy, and to estimate the effects of photoionization and photo-dissociation feedback, which regulate the star formation efficiency in molecular-cooling halos.

The adopted BH seeding prescriptions are the same as in T22: light BH seeds form as remnants of Pop III stars, and their masses ranges from 10s to 100s M, depending on the Pop III star formation efficiency. Following each Pop III star formation episode, the most massive among BH remnants is assumed to migrate to the centre of the host halo, becoming its nuclear light BH seed. Heavy BH seeds are assumed to form with a constant mass of 105 M from the collapse of SMSs in atomic-cooling halos (Tvir ≥ 104 K), when metal and dust cooling is still inefficient (ZZcr), and when molecular cooling is suppressed by a strong illuminating LW flux.5 Hence, in CAT heavy and light seeds co-exist and either of the two types of seed forms when the physical conditions allow. This means that the BH population predicted at each given redshift descends from both light or heavy BH seeds, with the exact match at late times determined by the ability or inability of the seeds to grow, as described in the following.

Once formed, BH seeds can grow by gas accretion and mergers. In particular, following galaxy major mergers (defined as those with interacting halos mass ratio greater than 1/10), the nuclear BHs are assumed to merge with their host galaxies. In doing so we do not account for the (admittedly poorly constrained, see Bortolas et al. 2020) dynamical friction time-scale, as well as for the subsequent binary hardening phase, and kick velocity of the newly formed BHs. All these processes may have important effects on the BH occupation fraction and on the estimated BH binary merger rates (see the discussion in Volonteri et al. 2021). While an implementation of physically motivated prescriptions to account for these processes in CAT is currently underway (Trinca et al. in preparation), we do not expect these processes to significantly affect our predictions, given that nuclear BHs grow pre-dominantly by gas accretion rather than by mergers.

Here, we consider CAT predictions for two alternative BH accretion models: in the first model, the accretion onto light and heavy BH seeds is EL (the model referred to as the reference model in T22); in the second model, light and heavy BH seeds can accrete at SE rates in short phases triggered by galaxy major mergers (SE, the model referred to as merger driven model in T22). The maximum duration of the SE phases is 1/10 of the host halo dynamical time-scale (Δtburst ≃ 2 − 7 Myr at z ≃ 10 − 20, see Section 2.4.2 in T22), but SE accretion is generally sustained for ≲ 2 − 5 Myr, consistent with predictions based on high-resolution hydrodynamical simulations (Regan et al. 2020; Massonneau et al. 2023; Sassano et al. 2023). When computing the BH bolometric luminosity, in the EL model we take a fixed radiative efficiency of ϵr = 0.1, while in the SE model we adopt a slim disc solution with the parametric form proposed by Madau et al. (2014), and we assume a constant BH spin of 0.572, so that when the BHs are not accreting at SE rates their ϵr → 0.1.

CAT has been used to explore observational imprints of light and heavy BH seeds and their growth mode on the redshift evolution of the BH mass and luminosity functions (see T22), their detectability within ongoing JWST surveys (Trinca et al. 2023b), the nature of UV-bright galaxies detected by JWST at z ≥ 10 (Trinca et al. 2023a), and the expected imprints of the first stellar populations and accreting BH seeds on the global 21 cm signal at cosmic dawn (Ventura et al. 2023).

BH masses and luminosities of AGNs detectable by JADES and CEERS at 5 ≲ z ≲ 11 are shown in Fig. 1. In each panel, the BH population grown from light and heavy BH seeds at the given redshift is shown by yellow points, while systems with luminosity brighter than the limiting magnitude of JADES and CEERS surveys are circled in violet and turquoise, respectively. Hence, turquoise points mark systems that can be observed by both surveys, while violet points mark systems that can only be detected by JADES. The upper and lower panels show, respectively, the systems predicted by CAT in the EL and SE models.

BH masses and bolometric luminosities predicted by CAT (T22) at 5 ≲ z ≲ 11 assuming the Eddington limited model (EL, top row) and the SE model (bottom row), where BHs can accrete at SE rates in short phases triggered by galaxy major mergers. In each panel, light yellow points represent all the systems at the corresponding redshift, which can descend from light and/or heavy BH seeds. Systems with a luminosity above the limiting sensitivity of JADES and CEERS are circled in violet and turquoise, respectively. The diagonal-dashed lines show reference values of λEdd = Lbol/LEdd = 0.01, 0.1, 1, and 10, and the vertical grey lines mark the initial heavy BH seeds mass adopted in CAT (105 M⊙). We also show the luminosities and BH masses identified in JWST surveys by Kocevski et al. (2023) (yellow pentagons at z ∼ 5), Übler et al. (2023) (blue star at z ∼ 5), Harikane et al. (2023) (green hollow stars at z ∼ 5 and ∼7); the candidate AGN in CEERS-1019 identified by Larson et al. (2023) is marked with a black star at z ∼ 9, the AGN in GNz11 reported by Maiolino et al. (2023a) is shown as a filled red star at z ∼ 11, and AGNs at lower redshift identified by Maiolino et al. (2023c) are shown by hollow red stars. In the rightmost panels we also show as a black star with errorbars the candidate AGN at z ≃ 10.3 reported by Bogdan et al. (2023). The properties of some of the observed systems at z ∼ 5 are consistent with predictions of the EL model, but the luminosity and BH masses of the remaining systems at z ∼ 5 and of most of the systems at z ≳ 7 are better reproduced in the SE model.
Figure 1.

BH masses and bolometric luminosities predicted by CAT (T22) at 5 ≲ z ≲ 11 assuming the Eddington limited model (EL, top row) and the SE model (bottom row), where BHs can accrete at SE rates in short phases triggered by galaxy major mergers. In each panel, light yellow points represent all the systems at the corresponding redshift, which can descend from light and/or heavy BH seeds. Systems with a luminosity above the limiting sensitivity of JADES and CEERS are circled in violet and turquoise, respectively. The diagonal-dashed lines show reference values of λEdd = Lbol/LEdd = 0.01, 0.1, 1, and 10, and the vertical grey lines mark the initial heavy BH seeds mass adopted in CAT (105 M). We also show the luminosities and BH masses identified in JWST surveys by Kocevski et al. (2023) (yellow pentagons at z ∼ 5), Übler et al. (2023) (blue star at z ∼ 5), Harikane et al. (2023) (green hollow stars at z ∼ 5 and ∼7); the candidate AGN in CEERS-1019 identified by Larson et al. (2023) is marked with a black star at z ∼ 9, the AGN in GNz11 reported by Maiolino et al. (2023a) is shown as a filled red star at z ∼ 11, and AGNs at lower redshift identified by Maiolino et al. (2023c) are shown by hollow red stars. In the rightmost panels we also show as a black star with errorbars the candidate AGN at z ≃ 10.3 reported by Bogdan et al. (2023). The properties of some of the observed systems at z ∼ 5 are consistent with predictions of the EL model, but the luminosity and BH masses of the remaining systems at z ∼ 5 and of most of the systems at z ≳ 7 are better reproduced in the SE model.

At all redshifts, the systems detectable by JADES extend to lower BH masses and luminosities, as expected by the fainter limiting magnitude of the survey compared to CEERS (Trinca et al. 2023b). In both models, the bolometric luminosity increases with BH mass, but with some discontinuities: in the EL model, BHs descending from light seeds grow very inefficiently, and their masses and luminosities remain ≲ 104 M and |${\rm Log} (L_{\rm bol}) [\rm erg/s] \lesssim 39.5$| at all redshifts. Above this mass and luminosity, there appears to be a mass gap, as – by construction – the minimum heavy BH seed mass adopted in CAT is MBH = 105 M (grey vertical lines). While this gap could be filled if a more continuous distribution of initial BH seeds masses were adopted (see e.g. Ferrara et al. 2014; Sassano et al. 2023)6, in the present version of the EL model all the systems detectable with JADES and CEERS appear to be descendants of heavy BH seeds (see T22 for a discussion on the consequences of the adopted accretion model on the evolution of the BH mass and luminosity functions).

In the SE model, even light BH seeds grow more efficiently and the mass gap disappears7 However, there is a luminosity gap, with two disjoint, parallel diagonal tracks, reflecting the short duration of the SE phase, which shift part of the systems on the upper track. The comparison with the diagonal dashed lines shows that systems on the upper track are characterized by Eddington ratios in the range 1 ≲ λEdd = Lbol/LEdd ≲ 10, while systems on the lower track have λEdd ≲ 0.01 − 0.1, smaller than in the EL model. The upper track becomes progressively less populated with cosmic time, due to the decreasing number of galaxy major mergers and their decreasing gas content (see also Pezzulli et al. 2016 for a similar finding).

As expected, in the SE model the number of detectable systems is larger than in the EL model, particularly at the highest redshift, where most of the systems are detectable when accreting at SE rates. Since BH growth is faster in the SE model, CEERS and JADES are predicted to observe systems with a broader range of BH masses, at all redshifts, with JADES reaching the limit of the initial heavy BH seed mass at z ≲ 11, with a few systems even crossing this limit at z ∼ 7 − 9, hence potentially probing BHs originated from lighter BH seeds. At lower redshifts, SE accreting BHs with masses ∼105 M become progressively rarer, as these systems have already grown, loosing memory of their initial BH seed mass (Valiante et al. 2018). Among the sample of AGNs detected by JWST, the BH mass and luminosity of a systems at z ∼ 5 are broadly consistent with predictions of the EL model, but most of the sources at z ≥ 7 are better reproduced by the SE model.

The distribution of simulated systems in the BH mass – stellar mass plane for the EL and the SE models is shown in Fig. 2, using the same colour-coding adopted in Fig. 1. In the z ∼ 5 panels, we also show the sample of local broad-line AGNs (light purple data points) from Reines & Volonteri (2015), and their best fit relation (grey shaded region), together with the updated fit to the larger sample from Greene, Strader & Ho (2020) (azure shaded region). We stress that – while pre-JWST observations based on dynamical mass measurements of z ≃ 6 quasars suggested a strong evolution of the relation with redshift (see Fan et al. 2022 for a recent review) – recent NIRCam/NIRSpec observations of two quasar at z ∼ 6 − 7 show that their location in the BH mass – stellar mass plane is consistent with the distribution at low redshift (Ding et al. 2023), suggesting that the most massive BHs inhabit the most massive galaxies.8

Similar to Fig. 1, but showing the BH mass – stellar mass relation predicted by CAT (T22) at 5 ≲ z ≲ 11 for the Eddington limited model (EL, top row) and the SE model (bottom row). The diagonal-dashed lines show reference values of MBH/Mstar from 10−4 to 1, and the horizontal grey lines mark the initial heavy BH seeds mass adopted in CAT (105 M⊙). The grey and azure shaded areas in the z ∼ 5 panels show, respectively, the fit of from Reines & Volonteri (2015) to their local sample of broad-line AGNs (light purple points), and the updated fit to the scaling relation from Greene et al. (2020). We also show the estimated BH and stellar masses (or upper limits) of systems identified in JWST surveys by Kocevski et al. (2023) (yellow pentagons at z ∼ 5), Übler et al. (2023) (blue star at z ∼ 5), Harikane et al. (2023) (green hollow stars at z ∼ 5 and ∼7); the candidate AGN in CEERS-1019 identified by Larson et al. (2023) is marked with a black star at z ∼ 9, the AGN in GNz11 reported by Maiolino et al. (2023a) is shown as a filled red star at z ∼ 11, and AGNs at lower redshift identified by Maiolino et al. (2023c) are shown by hollow red stars. In the rightmost panels we also show as a black star with errorbars the candidate AGN at z ≃ 10.3 reported by Bogdan et al. (2023). At all redshifts, most of the observable systems have MBH/Mstar ≲ 0.01 in the EL model, and 0.01 ≲ MBH/Mstar ≲ 0.1 in the SE model.
Figure 2.

Similar to Fig. 1, but showing the BH mass – stellar mass relation predicted by CAT (T22) at 5 ≲ z ≲ 11 for the Eddington limited model (EL, top row) and the SE model (bottom row). The diagonal-dashed lines show reference values of MBH/Mstar from 10−4 to 1, and the horizontal grey lines mark the initial heavy BH seeds mass adopted in CAT (105 M). The grey and azure shaded areas in the z ∼ 5 panels show, respectively, the fit of from Reines & Volonteri (2015) to their local sample of broad-line AGNs (light purple points), and the updated fit to the scaling relation from Greene et al. (2020). We also show the estimated BH and stellar masses (or upper limits) of systems identified in JWST surveys by Kocevski et al. (2023) (yellow pentagons at z ∼ 5), Übler et al. (2023) (blue star at z ∼ 5), Harikane et al. (2023) (green hollow stars at z ∼ 5 and ∼7); the candidate AGN in CEERS-1019 identified by Larson et al. (2023) is marked with a black star at z ∼ 9, the AGN in GNz11 reported by Maiolino et al. (2023a) is shown as a filled red star at z ∼ 11, and AGNs at lower redshift identified by Maiolino et al. (2023c) are shown by hollow red stars. In the rightmost panels we also show as a black star with errorbars the candidate AGN at z ≃ 10.3 reported by Bogdan et al. (2023). At all redshifts, most of the observable systems have MBH/Mstar ≲ 0.01 in the EL model, and 0.01 ≲ MBH/Mstar ≲ 0.1 in the SE model.

Fig. 2 shows that the BH mass gap is apparent in the EL model, where light BH seeds descendants with masses <104 M are undermassive with respect to their host galaxies. Even for the more efficiently growing heavy BH seeds, their masses remain ≲ 106 M until the host galaxy has grown to a mass of at least Mstar > 108 M, above which the growth of the BHs adjust to MBH/Mstar ≲ 0.001, with a large scatter (see also Habouzit, Volonteri & Dubois 2017 and Anglés-Alcázar et al. 2017). Interestingly, a small fraction of systems observable with JADES at z ≳ 7 is made by heavy BH seeds with ≲ 106 M which are over-massive with respect to their host galaxies, as expected if these systems keep memory of their peculiar birth conditions (Agarwal et al. 2013; Natarajan et al. 2017b; Valiante et al. 2018; Visbal & Haiman 2018).

In the SE model, instead, the vast majority of detectable BHs grow at a pace almost similar to that of their hosts, with MBH/Mstar ≃ 0.01 − 0.1; while this result might reflect the relatively simplified treatment of the SE accretion in CAT (see the discussion in T22 and in Section 3), it certainly shows how the assumed BH growth translates into a different evolution of the BH mass – stellar mass relation, even among models which rely on the same BH seeding conditions (Volonteri 2012).

At all redshifts, most of the observable systems have 0.001 ≲ MBH/Mstar ≲ 0.01 in the EL model, and 0.01 ≲ MBH/Mstar ≲ 0.1 in the SE model. The comparison with JWST detected AGNs shows that their luminosities are better reproduced by the SE model, while their host galaxy masses are better matched by the EL model, particularly at z ≳ 7. We also note that if the estimated photometric redshift, BH mass and host galaxy stellar mass of UHZ1 at z = 10.3 were to be confirmed by future observations (Castellano et al. 2022; Bogdan et al. 2023), this system – shown by the small black star with errorbars in the rightmost panel of the figure – would be more overmassive (in terms of BH to stellar mass ratio) than predicted by both models.

We nevertheless identify the subsample of detectable sources with BH masses that match the measured BH masses estimated for GNz11 (Maiolino et al. 2023a) and CEERS-1019 (Larson et al. 2023), and whose stellar mass (for sources in the EL model) or BH luminosity (for sources in the SE model) is consistent with the observed values. For these subsample of systems, which provide our best synthetic analogues for GNz-11 and CEERS-1019, we compute further the emission spectra and analyse their evolutionary tracks.

2.3 Emission spectra of the first accreting BHs and their hosts

The best synthetic candidates for GNz11 at z = 10.6 and CEERS-1019 at z = 8.68 are extracted from the source catalogs predicted by CAT in both the EL and SE models. We first select systems with a BH mass comparable to that estimated for the two observed sources (Larson et al. 2023; Maiolino et al. 2023a). Then, we also match two or more additional properties such as the estimated AGN bolometric luminosity, Eddington ratio (λEdd), host galaxy stellar mass, and metallicity (Larson et al. 2023; Tacchella et al. 2023b).

For CEERS-1019 we also explore the possibility that the system does not host an accreting nuclear BH, since the identification of the BLR through the broad H β emission line has a low signal-to-noise ratio (2.3σ). In this case, we select our synthetic candidates to have stellar masses and star formation rates within the observed range, with no additional condition on the BH mass and luminosity.

For each of these systems, we compute an incident continuum spectrum from the accreting BH and from star formation in the host galaxy, based on the properties (redshift, BH mass, and bolometric luminosity, stellar mass, age, and metallicity) predicted by CAT. We then use these spectra as input to cloudy photoionization models (we adopt version 22.00 of cloudy; Ferland et al. 2017) to compute the transmitted fluxes, adopting the interstellar medium (ISM) properties predicted by CAT for each system (gas metallicity, dust-to-gas, and dust-to-metals mass ratios), similarly to what done in our previous study (Valiante et al. 2018).

A systematic investigation of the SEDs for all the systems predicted by CAT to be detectable by JWST (and their position in spectral diagnostic diagrams) will be presented in Valiante et al. (in preparation). Here we show the SEDs for candidates which minimize the difference with the photometric data of CEERS-1019 and GNz11, according to our cloudy models.

2.3.1 Construction of the SEDs

The adopted models for the intrinsic emission spectra of accreting BHs and their hosts (used to define the incident radiation in cloudy), and the main set-up conditions adopted to run cloudy models are briefly described further.

Following Valiante et al. (2018), the primary optical/UV emission from the BH accretion disc is described as the sum of multicolour blackbody spectra, normalized to the BH bolometric luminosity. The X-ray emission from the hot corona is modelled as a power law with an exponential cut-off, |$L_\nu \propto \nu ^{-\Gamma +1} {\rm e}^{-h\nu /E_{\rm c}}$|⁠, with Ec = 300 keV and the photon index Γ in the hard X-ray energy band is assumed to depend on the Eddington ration λEdd as Γ2-10keV = 0.32LogλEdd + 2.27. We also include the metallicity-dependent X-ray reflection component using the reflection-only solution of the PEXRAV model in the xspec code (see Valiante et al. 2018 for more details).

The emission from the stellar component is computed with the spectral evolution code pegase v2.09 (Fioc & Rocca-Volmerange 2019), using as an input the ages and metallicities of the stellar populations formed in the host galaxy of the accreting nuclear BH predicted by CAT. We note here that when computing the stellar emission, we do not consider the full star formation history of the host galaxy, but we simply assume that all the stars formed in a burst and evolve over a times-scale of Δt = 2 (1) Myr, corresponding to the typical redshift-dependent time interval in CAT at z = 8 − 9 (10 − 11). We compute their cumulative emission considering their metallicity to be the same as the metallicity of their parent star forming cloud and a stellar age equal to Δt.

In cloudy models, we use the ISM properties (metallicity, dust-to-gas, dust-to-metal ratios) predicted by CAT to characterize the enrichment of the gas cloud. For systems with metallicity and dust-to-metal mass ratio close to local interstellar values (that we assume to be Z = Z = 0.01524, and ξd = ξ = 0.36, respectively), we use the interstellar abundances and depletion factors from Gutkin, Charlot & Bruzual (2016)10, except for Nitrogen, for which we adopt the depletion factor from Dopita et al. (2006, see Nakajima et al. 2018; Nakajima & Maiolino 2022). For sources with sub-solar metallicity, the corresponding abundance of elements lighter than Zinc are computed by rescaling the solar values by the metallicity predicted by CAT. In addition, following Gutkin et al. (2016), the depletion factors associated with different dust-to-metal mass ratios (ξd < ξ) are obtained with a linear interpolation of the depletion factors associated with ξd = 1, 0.36, and 1. To account for dust physics, we rescale the default ISM grains abundance of cloudy by the dust-to-gas mass ratio predicted by CAT. Grain size distribution and optical properties are from Mathis, Rumpl & Nordsieck (1977) and Martin & Rouleau (1991), already included in cloudy.

The clumping properties of the emission-line gas are defined by the filling factor (or clumping factor), f. For this study, we investigated both a uniform distribution model, f = 1, and clumpy gas models with f = 0.05, 0.1, and 0.5. For a given value of the dust-to-gas mass ratio (provided by CAT for each source), the adopted value of f translates into a different value of the visual extinction AV.

Following Nakajima & Maiolino (2022), in all models presented here, we assumed a constant Hydrogen density |$\rm log (n_H/cm^{-3})=3.0$| and we stop cloudy photoionization runs when the electron temperature, Te, drops below 100 K (in stellar models), or when a neutral Hydrogen column density of 1021 cm−2 is reached (in accreting BH models).

The resulting SEDs of systems that have been selected to be representative of the best synthetic candidates for CEERS-1019 at z = 8.68 and GNz11 at z = 10.6 are shown in Figs 3 and 4. In particular, we show the SEDs that minimize the difference with the photometric data points reported by Larson et al. (2023) for CEERS-1019 and by Tacchella et al. (2023b) for the point source of GNz11, which is attributed to an accreting BH by Maiolino et al. (2023a). It is important to emphasize that while the spectra of the observed sources show broad lines attributed to the BLR of the BH, cloudy does not calculate the line broadening caused by the BH potential.

Flux density predicted with cloudy for CEERS-1019 candidates selected from the EL model (upper panel) and from the SE (lower panel) model catalogs at z ∼ 9 predicted by CAT (see Section 2.2). The solid azure and dashed black lines in both panels show the SED of the accreting BH and stellar component, respectively, assuming that all the stars are formed in a single burst and a stellar age of 2 Myr. NIRcam and MIRI photometric data points from Larson et al. (2023) are shown in different colours, with bars indicating the filter band width. It is important to emphasize that while the observed spectrum shows a broad line attributed to the BLR of the BH, cloudy does not calculate the line broadening caused by the BH potential (see text).
Figure 3.

Flux density predicted with cloudy for CEERS-1019 candidates selected from the EL model (upper panel) and from the SE (lower panel) model catalogs at z ∼ 9 predicted by CAT (see Section 2.2). The solid azure and dashed black lines in both panels show the SED of the accreting BH and stellar component, respectively, assuming that all the stars are formed in a single burst and a stellar age of 2 Myr. NIRcam and MIRI photometric data points from Larson et al. (2023) are shown in different colours, with bars indicating the filter band width. It is important to emphasize that while the observed spectrum shows a broad line attributed to the BLR of the BH, cloudy does not calculate the line broadening caused by the BH potential (see text).

Flux density predicted with cloudy for GNz11 candidates selected from the EL and SE model catalogs at z ∼ 11 predicted by CAT (see Section 2.2). In the upper panel we show the emission from accreting BHs and compare these to the photometric data points reported by Tacchella et al. (2023b) for the point source of GNz11, which is attributed to an accreting BH by Maiolino et al. (2023a). The spectra are shown for a SE ($\rm \lambda _{Edd}=2.4$) accreting BH candidate with (azure solid line, AV = 0.4) and without (azure dotted line,AV = 0) dust attenuation. The sub-Eddington ($\rm \lambda _{Edd}=0.4$) candidate found in the z = 11 EL model catalogue is represented by the black dot-dashed line. For this latter system we predict the same metallicity as for the SE model, 0.06Z⊙ (from CAT), but AV = 1.0 (from the adopted cloudy model). For all these systems the BH mass is $1.5\times 10^{6} \, \rm M_\odot$. Using the same colour coding, in the bottom panel we show the emitted spectra of the corresponding host galaxies, whose stellar mass, metallicity, and dust attenuation are reported in the legend. For all these models, we have assumed all the stars to form in a single burst and a stellar age of 1 Myr. The grey shaded region shows the flux density predicted for a log(Mstar/M⊙) = 9 galaxy adopting the star formation history inferred by Tacchella et al. (2023b for the extended component, assuming a fixed metallicity of Z = 10−1.4Z⊙, and a dust attenuation in the range AV = 0 − 2. The model predictions are compared to the photometric data points reported by Tacchella et al. (2023b) for the extended component of GNz11. It is important to emphasize that while the observed spectrum shows a broad line attributed to the BLR of the BH, cloudy does not calculate the line broadening caused by the BH potential (see text).
Figure 4.

Flux density predicted with cloudy for GNz11 candidates selected from the EL and SE model catalogs at z ∼ 11 predicted by CAT (see Section 2.2). In the upper panel we show the emission from accreting BHs and compare these to the photometric data points reported by Tacchella et al. (2023b) for the point source of GNz11, which is attributed to an accreting BH by Maiolino et al. (2023a). The spectra are shown for a SE (⁠|$\rm \lambda _{Edd}=2.4$|⁠) accreting BH candidate with (azure solid line, AV = 0.4) and without (azure dotted line,AV = 0) dust attenuation. The sub-Eddington (⁠|$\rm \lambda _{Edd}=0.4$|⁠) candidate found in the z = 11 EL model catalogue is represented by the black dot-dashed line. For this latter system we predict the same metallicity as for the SE model, 0.06Z (from CAT), but AV = 1.0 (from the adopted cloudy model). For all these systems the BH mass is |$1.5\times 10^{6} \, \rm M_\odot$|⁠. Using the same colour coding, in the bottom panel we show the emitted spectra of the corresponding host galaxies, whose stellar mass, metallicity, and dust attenuation are reported in the legend. For all these models, we have assumed all the stars to form in a single burst and a stellar age of 1 Myr. The grey shaded region shows the flux density predicted for a log(Mstar/M) = 9 galaxy adopting the star formation history inferred by Tacchella et al. (2023b for the extended component, assuming a fixed metallicity of Z = 10−1.4Z, and a dust attenuation in the range AV = 0 − 2. The model predictions are compared to the photometric data points reported by Tacchella et al. (2023b) for the extended component of GNz11. It is important to emphasize that while the observed spectrum shows a broad line attributed to the BLR of the BH, cloudy does not calculate the line broadening caused by the BH potential (see text).

2.3.2 CEERS-1019

For CEERS-1019, we show the predicted emission spectra for the accreting BHs (azure solid lines) and host stellar populations (dashed black lines) of two systems, selected from the EL (upper panel) and SE models (lower panel). The two candidates have comparable BH mass and were selected among the closest CEERS-1019 synthetic counterparts (in terms of BH mass, luminosity, and stellar mass) predicted by CAT at z ∼ 9.

In both panels, we report their main physical properties and we compare the emission spectra with NIRcam and MIRI photometric data points from Larson et al. (2023), that are shown in different colours, with bars indicating the filter band width.

In the EL model, the continuum is dominated by the emission from the host galaxy, which has a total stellar mass of log(Mstar/M) = 9.4. Here we have assumed that all stars formed in a single burst and a stellar age of 2 Myr, that is the typical time interval at z ∼ 8 − 9 in CAT. Despite this extreme assumption, some of the emitted lines are heavily contaminated by the accreting BH, which has a mass of log(MBH/M) = 6.9 and Eddington ratio of λEdd = 0.45.

In the SE model, instead, the continuum and emitted lines are dominated by the BH emission, which has a mass of log(MBH/M) = 6.8 and Eddington ratio of λEdd = 1.8. Due to the faster growth of the BH, the host galaxy has a stellar mass of only log(Mstar/M) = 8.0.

The metallicities of the two systems clearly reflect their different degree of evolution, being ≃ 0.08Z in the EL model, and a factor ∼1.6 smaller in the SE model. cloudy models with f = 0.05 provide AV = 0.3 − 0.6, in good agreement with the value inferred for CEERS-1019 (0.4 ± 0.2; Larson et al. 2023). The only exception is the host galaxy in the EL model (dashed line in the upper panel), for which we show the SED with attenuation AV = 1.0, obtained assuming f = 0.1. For a filling factor f < 0.1 (AV < 1.0) the continuum emission of this stellar component exceeds NIRCam observations at λ < 3.5 μm.

We find that the SE model, where the flux density is dominated by the accreting BH, provides a poor fit to the observed photometry of CEERS-1019. This conclusion generally applies to SE candidates with comparable BH mass and luminosity, and supports the interpretation provided by Larson et al. (2023), that is, that – if CEERS-1019 hosts an AGN – the emission from the host galaxy comprises between |$\sim 80 - 99.5~{{\ \rm per\ cent}}$| of the UV/optical flux, depending on the adopted AGN template.

We also explore the possibility that CEERS-1019 does not host an AGN, that is, that the 2.3σ detection of the H β broad line component is not confirmed by future observations. In this case, we select from CAT simulated galaxies at z ≃ 9 systems with stellar masses and star formation rates consistent with the values inferred by Larson et al. (2023), without imposing any additional condition on the nuclear BH mass and luminosity. In the SE model, all the selected galaxies host nuclear BHs with masses log(MBH/M) = 7.4 − 8.5 and λEdd = 0.009 − 0.03, that is, more massive and less active than inferred by Larson et al. (2023). In the EL model, we find that |$\simeq 10~{{\ \rm per\ cent}}$| of the selected galaxies host nuclear BHs that have properties similar to the system shown in the top panel of Fig. 3, while the rest of the sample host smaller (log(MBH/M) = 5.25 − 6.5), and less luminous BHs (λEdd = 0.03 − 0.3).

2.3.3 GNz11

For GNz11, Fig. 4 compares the predicted emission spectra powered by accreting BHs (upper panel) and their host galaxies (lower panel) with the photometric data reported by Tacchella et al. (2023b) for the central point source (which Maiolino et al. 2023a attribute to the AGN emission) and the extended component, respectively.

In the upper panel, we show the flux density predicted for two accreting BHs with mass |$M_{\rm BH}\sim 1.5 \times 10^6 \, \rm M_\odot$|⁠, selected from the SE model (λEdd = 2.4, azure lines) and from the EL model (λEdd = 0.4, black dashed-dotted line). Using the same colour-coding, in the bottom panel we show the flux density predicted for their host galaxies. The two accreting BHs are hosted by galaxies with different stellar mass (log(Mstar/M) = 8.6 in the EL model, and log(Mstar/M) = 8 in the SE model), but with the same ISM metallicity (Z = 0.06Z).

To best match the photometry of the point-source and extended component, we performed different cloudy runs assuming different filling factors, and found that the best match is provided by f = 0.05 for the SE model, corresponding to AV = 0.4 (azure solid line), and f = 0.5 for the EL model, corresponding to AV = 1 (the latter choice provides the best fit to the photometry of the extended component). For the SE model, we also show a model with no dust obscuration (AV = 0, azure dotted line). We find that the SED of the SE model candidate is in good agreement with the observed photometry of the point source, especially if dust extinction is neglected, while the EL model candidate is too faint. By assuming f = 0.05, we can reduce the visual extinction to AV = 0.3 (more consistent with that inferred by Tacchella et al. 2023b), and the emission from the EL accreting BH increases (by a factor ≈3), but it is still too small to reproduce the observations. With the same assumption, however, the emission of the stellar component is too bright to match the photometry of the extended component.

We recall that the stellar emission is computed assuming that the stars are formed in a single burst and, in this case, have an age of 1 Myr (approximately the CAT time interval at z ∼ 10 − 11). To appreciate the impact of considering the contribution of older stellar populations, we also show the flux density predicted for a galaxy with total stellar mass of log(Mstar/M) = 9 adopting the star formation history inferred by the SED fitting of the extended component of GNz11 by Tacchella et al. (2023b), with the grey shaded region corresponding to models with different degrees of dust obscuration, in the range AV = 0 − 2. This shows that a better modelling of the star formation history of the host galaxy is required to properly account for the observations, and that the galaxies hosting the SE accreting BHs may be too faint to account for the observed flux density of the extended component, particularly at 2–4 μm.

Our analysis suggests that the emission of the central point source in GNz11 can be reproduced by SE accreting BHs, with mass and luminosity comparable to the values inferred by Maiolino et al. (2023a); however, these systems are hosted in galaxies which are too faint to match the observed emission of the extended component.

2.4 Tracking the origin and fate of z ≃ 9 − 11 SMBHs

The analysis of the emission properties of the best candidates of GNz11 and CEERS-1019 (among which we show in Figs 3 and 4 some representative cases) supports the interpretation that the central point source of GNz11 could be powered by a |$1.5 \times 10^{6} \, \rm M_\odot$| SE accreting BH (λEdd ≃ 2 − 3), while CEERS-1019 may harbour a more massive BH, with MBH ≳ 107 M, accreting sub-Eddington (λEdd ≲ 0.5), with a dominant emission from the host galaxy.

Here we go one step further, and we investigate the origin and fate of these systems. In Fig. 5 we show their BH mass as a function of redshift (lower horizontal axis) and time (upper horizontal axis). Evolutionary tracks for the best candidates of CEERS-1019 at z = 8.7 (filled black star) in the EL model, and of GNz11 at z = 10.6 (filled red star) in the SE model are shown with blue and orange solid lines in the upper and lower panels, respectively. We also show with a smaller filled black star with errorbars the candidate system in UHZ1 at z = 10.3. The dark blue and red coloured points along the tracks mark the occurrence of BH mergers, while grey points indicate the mass of the initial BH seed and the mass at the observed redshift of GNz11 and CEERS-1019. We also show the BH masses and redshifts of AGNs identified by JWST at z < 7, using the same symbols and colours adopted in Figs 1 and 2. The lavender data points indicate the BHs measured in quasars at z ∼ 6 − 7.5 (Inayoshi et al. 2020). The grey horizontal bands identify the mass ranges corresponding to light and heavy BH seeds in CAT. In both panels, the dotted lines with triangles indicate the maximum BH mass that is predicted by CAT at 4 ≤ z ≤ 20 (yellow), and the minimum BH mass that would be detectable in a survey like CEERS (turquoise) and JADES (violet). We also mark with vertical dotted lines the maximum redshift at which BHs would be observable in these two surveys, that correspond to the redshift at which the minimum detectable BH mass equals the maximum BH mass predicted by each model.

BHs evolutionary tracks predicted by CAT for the Eddington limited model (EL, top panel) and the SE model (bottom panel). In top (bottom) panel, the tracks of the best candidates for CEERS-1019 and GNz11 are shown with blue (orange) solid lines, with darkblue (red) data points indicating the occurrence of BH mergers, and the grey data points indicating the BH mass at seed formation. The horizontal grey shaded regions mark the mass ranges adopted in CAT for light and heavy BH seeds formation. The candidate AGN in CEERS-1019 identified by Larson et al. (2023) is marked with a black star at z = 8.7 (tH = 0.570 Gyr), the candidate AGN in UHZ1 at z ≃ 10.3 (tH = 0.46 Gyr) reported by Bogdan et al. (2023) is shown by a smaller black star with errorbars. Finally, the AGN in GNz11 reported by Maiolino et al. (2023a) is shown as a red star at z = 10.6, (tH = 0.44 Gyr), and AGNs at lower redshift identified by Maiolino et al. (2023c) are shown by hollow red stars. We also show the BH mass and redshift of systems identified in JWST surveys by Kocevski et al. (2023) (yellow pentagons), Übler et al. (2023) (blue star), Harikane et al. (2023) (green hollow stars), and the BH masses measured in quasars at z ∼ 6 − 7.5 (lavender data points, Inayoshi et al. 2020). Finally, in both panels, triangles connected with dotted lines indicate the maximum BH mass predicted in each model (yellow), and the minimum BH mass detectable in a survey like CEERS (turquoise) and JADES (violet). The vertical dotted lines show the corresponding maximum observable redshift in the two surveys. According to the EL model, systems like GNz11 and CEERS-1019 are among the most massive BHs at their corresponding redshift, and must originate from heavy BH seeds formed at 14 ≤ z ≤ 16. In the SE model, they can originate from rapidly growing light seeds formed at z ≥ 20 and/or from heavy BH seeds formed at 14 ≤ z ≤ 16, and do not represent the most extreme objects. In both models, we find that GNz11 may be a progenitor of CEERS-1019, if the latter hosts an AGN (see text), and that these BHs continue to grow, reaching a BH mass comparable to those measured in quasars at 6 ≤ z ≤ 7.5. The estimated mass of the candidate AGN in UHZ1 exceeds the maximum BH mass predicted by EL model at z = 10.3, but is smaller than the maximum BH mass predicted by the SE.
Figure 5.

BHs evolutionary tracks predicted by CAT for the Eddington limited model (EL, top panel) and the SE model (bottom panel). In top (bottom) panel, the tracks of the best candidates for CEERS-1019 and GNz11 are shown with blue (orange) solid lines, with darkblue (red) data points indicating the occurrence of BH mergers, and the grey data points indicating the BH mass at seed formation. The horizontal grey shaded regions mark the mass ranges adopted in CAT for light and heavy BH seeds formation. The candidate AGN in CEERS-1019 identified by Larson et al. (2023) is marked with a black star at z = 8.7 (tH = 0.570 Gyr), the candidate AGN in UHZ1 at z ≃ 10.3 (tH = 0.46 Gyr) reported by Bogdan et al. (2023) is shown by a smaller black star with errorbars. Finally, the AGN in GNz11 reported by Maiolino et al. (2023a) is shown as a red star at z = 10.6, (tH = 0.44 Gyr), and AGNs at lower redshift identified by Maiolino et al. (2023c) are shown by hollow red stars. We also show the BH mass and redshift of systems identified in JWST surveys by Kocevski et al. (2023) (yellow pentagons), Übler et al. (2023) (blue star), Harikane et al. (2023) (green hollow stars), and the BH masses measured in quasars at z ∼ 6 − 7.5 (lavender data points, Inayoshi et al. 2020). Finally, in both panels, triangles connected with dotted lines indicate the maximum BH mass predicted in each model (yellow), and the minimum BH mass detectable in a survey like CEERS (turquoise) and JADES (violet). The vertical dotted lines show the corresponding maximum observable redshift in the two surveys. According to the EL model, systems like GNz11 and CEERS-1019 are among the most massive BHs at their corresponding redshift, and must originate from heavy BH seeds formed at 14 ≤ z ≤ 16. In the SE model, they can originate from rapidly growing light seeds formed at z ≥ 20 and/or from heavy BH seeds formed at 14 ≤ z ≤ 16, and do not represent the most extreme objects. In both models, we find that GNz11 may be a progenitor of CEERS-1019, if the latter hosts an AGN (see text), and that these BHs continue to grow, reaching a BH mass comparable to those measured in quasars at 6 ≤ z ≤ 7.5. The estimated mass of the candidate AGN in UHZ1 exceeds the maximum BH mass predicted by EL model at z = 10.3, but is smaller than the maximum BH mass predicted by the SE.

Our results suggest that there might be an evolutionary link between systems with a BH mass similar to GNz11 at z = 10.6, and systems with a mass similar to CEERS-1019 at z = 8.7. However, their formation history is very different in the two models: in the EL model, systems like GNz11 and CEERS-1019 originate from heavy BH seeds formed at 14 ≤ z ≤ 16, and are among the most massive BHs at their corresponding redshift. In the SE model, instead, all the tracks start at z ≳ 20 from light BH seeds, although mergers can occur with heavy BH seeds as they cross the heavy BH seed mass, at z ≤ 16, indicating that the genealogy of systems like GNz11 and CEERS-1019 could be characterized by both light and heavy BH seeds ancestors. This is consistent with the results presented in T22 (see in particular their Fig. 7), which show that less than |$\simeq 20~{{\ \rm per\ cent}}$| of BHs with masses ≃ 106 M − 107 M at 8 ≤ z ≤ 10 descend from a heavy BH seed. We also note that in the SE model, the evolutionary tracks are always well below the maximum BH mass predicted by the model, meaning that systems similar to GNz11 and CEERS-1019 are not the most massive among the BH population at their corresponding redshift. Conversely, the candidate AGN at z = 10.3 reported by Bogdan et al. (2023), if confirmed, would be among the most extreme BHs, in terms of BH mass, predicted by the SE models, and too massive to be formed in the EL model. Interestingly, at z < 9 − 11 the selected systems continue to grow in mass: in both models they could be the progenitors of SMBHs powering quasars at 6 ≤ z ≤ 7, and of the most massive among the AGNs detected with JWST at z < 7. Conversely, JWST observed AGNs with mass ≲ 3 × 107 M at z < 7 must originate from less efficiently growing seeds (see for instance the scenario proposed by Sassano et al. 2023), and their evolution is therefore unrelated to GNz11-like or CEERS-1019-like systems.

If we compare the maximum BH mass predicted in each model with the minimum BH mass that would be observable by CEERS and JADES, we find that CEERS is not expected to observe accreting BHs beyond z ≃ 9 in the EL model, and CEERS-1019 represent the most massive BH at z = 8.7 and it is ∼0.6 dex larger than the minimum detectable BH mass. JADES has the sensitivity to observe accreting BHs out to z ≃ 12, and GNz11 is ∼0.4 dex smaller than the maximum BH mass at z = 10.6 and ∼0.6 dex larger than the minimum detectable BH mass. In the SE model these figures change considerably: we predict that CEERS would be able to detect BHs with masses ≳ 3 × 105 M − 106 M out to z ≃ 15, and JADES would be able to detect BHs with masses ≳ 105 M − 3 × 105 M out to z ∼ 18.

While BHs with masses comparable to those of GNz11 and CEERS-1019 can be explained in both models, the predicted discovery space is very different in the two scenarios: we can expect that JWST will continue to discover small but luminous AGNs at z ≲ 9 − 12. However, if future surveys as sensitive as CEERS and JADES were to discover accreting BHs at z > 12, this would be a strong indication that early BH seeds must experience phases of SE growth.

3 DISCUSSION AND CONCLUSIONS

The recent discovery of more than 40 new AGNs at z > 4 in ongoing surveys with JWST, including one system at z = 10.6, with estimated BH masses ≤107 M, offers the unprecedented opportunity to constrain theories for the formation and evolution of SMBHs with masses >109 M powering quasars at z > 6. Given the high rate of AGNs discovered with JWST, two questions naturally arise: are we surprised to find these accreting BHs in the nuclei of galaxies at z ∼ 5 − 11? And can we use their estimated properties to constrain the nature and growth of the first BH seeds?

Prompted by these questions, we use the CAT (T22) semi-analytical model to explore the population of BHs and their host galaxies at 5 ≤ z ≤ 11, and compare their predicted properties to those estimated for the AGNs and host galaxies recently detected by JWST.

Following T22, we assume that light BH seeds form as remnants of Pop III stars and their masses range from 10s to 100s M, while heavy BH seeds form with a constant mass of 105 M from the collapse of supermassive stars. Starting from these two BH seeds flavours, we consider two alternative BH accretion models: in the first model, BH accretion is EL, while in the second model, light and heavy BH seeds can accrete at SE rates in short phases triggered by galaxy major mergers (SE). We note here that our description of SE accretion via a slim disc (Lupi et al. 2014; Madau et al. 2014; Volonteri et al. 2015; Pezzulli et al. 2016) is simplified (see the discussion in T22). High-resolution hydrodynamical simulations show that in this regime BH growth is very intermittent, due to radiative and mechanical feedback from the growing BH (Regan et al. 2019; Massonneau et al. 2023), and that efficient, SE growth may be possible onto heavy BH seeds (Inayoshi et al. 2022), and less likely onto seeds with masses ≤103 M (Smith et al. 2018; Sassano et al. 2023). With this caveat in mind, in the present study we have considered both the EL and the SE models to explore how the assumed BH growth compares with current observational constraints provided by the recently JWST detected AGNs.

We first consider the predicted distribution of nuclear BHs and their host galaxies in the BH mass – luminosity and BH – stellar mass planes at 5 ≤ z ≤ 11, and investigate the properties of BHs that are luminous enough to be detectable by JWST in surveys like JADES and CEERS (see also Trinca et al. 2023b). We find that

  • in the EL model, BHs descending from light seeds grow very inefficiently, their masses and luminosities remain ≲ 104 M and |${\rm Log} (L_{\rm bol}) [\rm erg/s] \lesssim 39.5$|⁠, and are undermassive with respect to their host galaxies. All the systems detectable with JADES and CEERS appear to be descendants of heavy BH seeds. Even for these more efficiently growing seeds, their masses remain ≲ 106 M until their host galaxy has grown to a mass of at least Mstar > 108 M, above which the growth of the BHs adjust to MBH/Mstar ≲ 0.001, with a large scatter. At all redshifts, most of the observable systems have 0.001 ≲ MBH/Mstar ≲ 0.01. Interestingly, a small fraction of systems observable with JADES at z ≳ 7 is made by heavy BH seeds with ≲ 106 M which are over-massive with respect to their host galaxies, as expected if these systems keep memory of their peculiar birth conditions (Agarwal et al. 2013; Natarajan et al. 2017b; Valiante et al. 2018; Visbal & Haiman 2018).

  • In the SE model, even light BH seeds grow more efficiently. However, there is a luminosity gap in the BH mass – luminosity plane, with two disjoint, parallel diagonal tracks, corresponding to 1 ≲ λEdd ≲ 10 and λEdd ≲ 0.01 − 0.1. This reflects the short duration of the SE phase, during which systems are shifted on the upper track, which becomes progressively less populated with decreasing redshift. The number of detectable systems is larger than in the EL model, particularly at z ≳ 12, where most of the systems are detectable when accreting at SE rates. This prediction could be test observationally with surveys like JADES and CEERS. At all redshifts, detectable BHs evolve with MBH/Mstar ≃ 0.01 − 0.1, and are overmassive with respect to the canonical value MBH/Mstar ≃ 0.001 (Reines & Volonteri 2015; Greene et al. 2020).

  • The comparison with JWST detected AGNs at 4 < z < 10.6 shows that their estimated BH luminosities are better reproduced by the SE model, while their host galaxy masses are better matched by the EL model. In particular, we identify the subsample of detectable sources with BH masses that match the estimated values for GNz11 (Maiolino et al. 2023a) and CEERS-1019, if this system hosts an AGN (Larson et al. 2023), and whose stellar mass (for sources in the EL model) or luminosity (for sources in the SE model) is consistent with the observed values. For these subsample of systems, which provide our best synthetic AGN analogues for GNz-11 and CEERS-1019, we compute their SEDs by coupling CAT predictions with the cloudy photoionization code, and compare them with the observed JWST photometry. Our analysis supports the interpretation that the central point source of GNz11 could be powered by a |$1.5 \times 10^{6} \, \rm M_\odot$| SE accreting BH (λEdd ≃ 2 − 3), while the observed emission from CEERS-1019 is dominated by the host galaxy: if CEERS-1019 harbours an accreting BH, we find it to have a mass of MBH ≃ 107M, and to be accreting at sub-Eddington rates (λEdd ≃ 0.5). Yet, these conditions are met in only 10  per cent of the systems predicted by CAT EL model at z ≃ 9, when these are selected to have stellar masses and star formation rates consistent with the values inferred for CEERS-1019 (Larson et al. 2023). In the rest of the EL sample, our synthetic counterparts of CEERS-1019 host smaller (log(MBH/M) = 5.25 − 6.5), and less active BHs (λEdd = 0.03 − 0.3), or – in the SE model – more massive (log(MBH/M) = 7.4 − 8.5) and less active BHs (λEdd = 0.009 − 0.03).

  • We also analyse the evolutionary tracks of CEERS-1019 and GNz11 synthetic analogues in the EL and in the SE models to investigate the origin and fate of these systems. We find that there might be an evolutionary link between systems with a BH mass similar to GNz11 at z = 10.6 and systems with a BH mass at z = 8.7 similar to CEERS-1019, if this is an AGN. However, their formation history is very different in the two models: in the EL model, systems like GNz11 and CEERS-1019 originate from heavy BH seeds formed at 14 ≤ z ≤ 16, and are among the most massive BHs at their corresponding redshift. In the SE model, instead, systems similar to GNz11 and CEERS-1019 are not the most massive among the BH population at their corresponding redshift; their evolutionary tracks start at z ≳ 20 from light BH seeds, although mergers can occur with heavy BH seeds at z ≤ 16, indicating that they could descend from both light and heavy BH seeds ancestors. In both models GNz11-like or CEERS-1019-like systems continue to grow at z < 9 − 11, reaching masses comparable to the SMBHs powering quasars at 6 ≤ z ≤ 7, and of the most massive among the AGNs detected with JWST at z < 7. Conversely, the sample of JWST AGNs with mass ≲ 3 × 107 M at z < 7 must originate from less efficiently growing seeds.

Our results suggest that the AGNs discovered in ongoing JWST surveys are providing an amazing confirmation of the rich BH landscape at cosmic dawn predicted by theoretical models. AGNs with properties similar to those estimated for JWST detected AGNs, including the most distant source, GNz11 and CEERS-1019, if this hosts an AGN, can be explained in scenarios where BHs originate from EL gas accretion onto heavy BH seeds, or if light and heavy BH seeds can experience short phases of SE accretion. We expect that surveys like CEERS and JADES will continue to discover small but luminous AGNs at z ≲ 9 − 12 (Trinca et al. 2023b). However, if JWST surveys as sensitive as CEERS and JADES were to discover accreting BHs at z > 12, and/or the BH mass estimated for the AGN candidate hosted in UHZ1 with a photometric redshift of z = 10.3 were to be confirmed with future observations, we would have a strong indication that conditions enabling SE accretion, such as a high-optical depth and high-gas supply rates from large scales (Mayer 2019), should be naturally met and that feedback in the SE phase does not curtail significantly its duration.

ACKNOWLEDGEMENTS

We are very grateful to the Referee, Mike Strauss, for his very constructive and insightful comments. RS and LG acknowledge support from the Amaldi Research Centre funded by the MIUR programme ‘Dipartimento di Eccellenza’ (CUP:B81I18001170001) and from the INFN TEONGRAV specific initiative. MV thanks the Physics Department of Sapienza University of Rome and the Amaldi Research Center for the kind hospitality and the vibrant scientific environment. RM acknowledges support by the Science and Technology Facilities Council (STFC), by the ERC through Advanced Grant 695671 “QUENCH”, and by the UKRI Frontier Research grant RISEandFALL. RM also acknowledges funding from a research professorship from the Royal Society.

DATA AVAILABILITY

The simulated data underlying this article will be shared on reasonable request to the corresponding author.

Footnotes

1

We note here that the stellar mass of the system has been inferred using galaxy SED templates, that with no AGN contribution. However, the contribution of the AGN to the observed photometry should be subdominant, given the large gas column density (⁠|$N_{\rm H} \simeq 10^{24} \, {\rm cm}^{-2}$|⁠) inferred by the analysis of the X-ray data (Bogdan et al. 2023).

2

We use here a definition taken from the title of the paper by Kocevski et al. (2023), which reads Hidden Little Monsters: Spectroscopic Identification of Low-Mass, Broad-Line AGN at z > 5 with CEERS.

3

This estimate was based on the photometric performance requirements of the instrument, namely the point source faintest fluxes that can be obtained at a signal-to-noise ratio SNR  = 10 in a 10 ks integration (Valiante et al. 2018).

4

We assume a Lambda cold dark matter (ΛCDM) cosmological model with the following parameters: |$\Omega _\Lambda = 0.685$|⁠, ΩM = 0.315, h = 0.674, Ωb = 0.05 (Planck Collaboration et al. 2020).

5

The latter condition is usually expressed as JLWJcr = 300, where JLW is the cumulative flux into the Lyman–Werner energy band in units of 10−21 erg s−1 cm−2 Hz−1 sr−1 (see T22).

6

Medium weight BH seeds with masses ∼103 M are predicted to form by runaway collisions in the dense stellar clusters (Omukai, Schneider & Haiman 2008; Devecchi & Volonteri 2009; Devecchi et al. 2012), and a more continuous distribution between medium weight and heavy BH seeds is found within the super-competitive accretion model recently explored by Chon & Omukai (2020).

7

We note here that in our merger-driven SE accretion model light BH seeds grow more efficiently than predicted by other studies (see e.g. Smith et al. 2018). We will return to this point in Section 3.

8

We note here that the inferred relation at z ≃ 6 is strongly affected by potential selection biases and by using gas tracers (Volonteri & Stark 2011; Lupi et al. 2019), and that massive quasar hosts were indeed expected by evolutionary scenarios that attempt to reproduce their inferred large dust masses (Valiante et al. 2014).

10

The abundances computed by Gutkin et al. (2016) also account for both primary and secondary nucleosynthesis production of Nitrogen and Carbon.

References

Adams
N. J.
et al. ,
2023
,
MNRAS
,
518
,
4755

Agarwal
B.
,
Davis
A. J.
,
Khochfar
S.
,
Natarajan
P.
,
Dunlop
J. S.
,
2013
,
MNRAS
,
432
,
3438

Anglés-Alcázar
D.
,
Faucher-Giguère
C.-A.
,
Quataert
E.
,
Hopkins
P. F.
,
Feldmann
R.
,
Torrey
P.
,
Wetzel
A.
,
Kereš
D.
,
2017
,
MNRAS
,
472
,
109

Arrabal Haro
P.
et al. ,
2023
,
Nature
.
in press, preprint
()

Backhaus
B. E.
et al. ,
2022
,
ApJ
,
926
,
L161

Baldwin
J. A.
,
Phillips
M. M.
,
Terlevich
R.
,
1981
,
PASP
,
93
,
5

Bañados
E.
et al. ,
2018
,
Nature
,
553
,
473

Barrow
K. S. S.
,
Aykutalp
A.
,
Wise
J. H.
,
2018
,
Nature Astron.
,
2
,
987

Barrufet
L.
et al. ,
2023
,
MNRAS
,
522
,
449

Bogdan
A.
et al. ,
2023
,
preprint
()

Bortolas
E.
,
Capelo
P. R.
,
Zana
T.
,
Mayer
L.
,
Bonetti
M.
,
Dotti
M.
,
Davies
M. B.
,
Madau
P.
,
2020
,
MNRAS
,
498
,
3601

Bunker
A. J.
et al. ,
2023
,
A&A
.
in press, preprint
()

Cameron
A. J.
et al. ,
2023
,
preprint
()

Carnall
A. C.
et al. ,
2023
,
MNRAS
,
520
,
3974

Castellano
M.
et al. ,
2022
,
ApJ
,
938
,
L15

Chon
S.
,
Omukai
K.
,
2020
,
MNRAS
,
494
,
2851

Cole
S.
,
Lacey
C. G.
,
Baugh
C. M.
,
Frenk
C. S.
,
2000
,
MNRAS
,
319
,
168

Curti
M.
et al. ,
2023
,
MNRAS
,
518
,
425

Curtis-Lake
E.
et al. ,
2023
,
Nature Astronomy
,
7
,
622

Devecchi
B.
,
Volonteri
M.
,
2009
,
ApJ
,
694
,
L302

Devecchi
B.
,
Volonteri
M.
,
Rossi
E. M.
,
Colpi
M.
,
Portegies Zwart
S.
,
2012
,
MNRAS
,
421
,
1465

Ding
X.
et al. ,
2023
,
Nature
.
621
,
51

Donnan
C. T.
et al. ,
2023
,
MNRAS
,
518
,
6011

Dopita
M. A.
et al. ,
2006
,
ApJS
,
167
,
177

Endsley
R.
,
Stark
D. P.
,
Whitler
L.
,
Topping
M. W.
,
Chen
Z.
,
Plat
A.
,
Chisholm
J.
,
Charlot
S.
,
2023
,
MNRAS
.

Fan
X.
,
Banados
E.
,
Simcoe
R. A.
,
2022
,
ARA&A
,
61
,
373

Ferland
G. J.
et al. ,
2017
,
RMxAA
,
53
,
385

Ferrara
A.
,
Salvadori
S.
,
Yue
B.
,
Schleicher
D.
,
2014
,
MNRAS
,
443
,
2410

Fioc
M.
,
Rocca-Volmerange
B.
,
2019
,
A&A
,
623
,
143

Fujimoto
S.
et al. ,
2022
,
Nature
,
604
,
261

Furtak
L. J.
et al. ,
2023
,
 ApJ
,
952
,
142

Goulding
A. D.
,
Greene
J. E.
,
2022
,
ApJL
,
938
,
L9

Greene
J. E.
,
Strader
J.
,
Ho
L. C.
,
2020
,
ARA&A
,
58
,
257

Gutkin
J.
,
Charlot
S.
,
Bruzual
G.
,
2016
,
MNRAS
,
462
,
1757

Habouzit
M.
,
Volonteri
M.
,
Dubois
Y.
,
2017
,
MNRAS
,
468
,
3935

Harikane
Y.
et al. ,
2022
,
ApJS
,
259
,
20

Harikane
Y.
et al. ,
2023
,
preprint
()

Inayoshi
K.
,
Visbal
E.
,
Haiman
Z.
,
2020
,
ARA&A
,
58
,
27

Inayoshi
K.
,
Onoue
M.
,
Sugahara
Y.
,
Inoue
A. K.
,
Ho
L. C.
,
2022
,
ApJ
,
931
,
L25

Kashino
D.
,
Lilly
S. J.
,
Matthee
J.
,
Eilers
A.-C.
,
Mackenzie
R.
,
Bordoloi
R.
,
Simcoe
R. A.
,
2023
,
ApJ
,
950
,
L66

Kocevski
D. D.
et al. ,
2023
,
ApJ
,
954
,
L4

Labbe
I.
et al. ,
2023a
,
Nature
,
616
,
266

Labbe
I.
et al. ,
2023b
,
preprint
()

Larson
R. L.
et al. ,
2023
,
ApJ
,
953
,
L29

Looser
T. J.
et al. ,
2023
,
preprint
()

Lupi
A.
,
Colpi
M.
,
Devecchi
B.
,
Galanti
G.
,
Volonteri
M.
,
2014
,
MNRAS
,
442
,
3616

Lupi
A.
,
Volonteri
M.
,
Decarli
R.
,
Bovino
S.
,
Silk
J.
,
Bergeron
J.
,
2019
,
MNRAS
,
488
,
4004

Madau
P.
,
Haardt
F.
,
Dotti
M.
,
2014
,
ApJ
,
784
,
L38

Maiolino
R.
et al. ,
2023a
,
preprint
()

Maiolino
R.
et al. ,
2023b
,
preprint
()

Martin
P. G.
,
Rouleau
F.
,
1991
, in
Malina
R. F.
,
Bowyer
S.
, eds,
First Berkeley Colloquium on Extreme Ultraviolet Astronomy
.
Pergamon Press
,
New York
, p.
341

Massonneau
W.
,
Dubois
Y.
,
Volonteri
M.
,
Beckmann
R. S.
,
2023
,
A&A
,
669
,
143

Mathis
J. S.
,
Rumpl
W.
,
Nordsieck
K. H.
,
1977
,
ApJ
,
217
,
L425

Matthee
J.
et al. ,
2023
,
preprint
()

Mayer
L.
,
2019
, in
Latif
M.
,
Schleicher
D.
eds,
Formation of the First Black Holes
.
World Scientific Publishing
, p.
195

Naidu
R. P.
et al. ,
2022a
,
preprint
()

Naidu
R. P.
et al. ,
2022b
,
ApJ
,
940
,
L14

Nakajima
K.
,
Maiolino
R.
,
2022
,
MNRAS
,
513
,
5134

Nakajima
K.
et al. ,
2018
,
A&A
,
612
,
94

Nakajima
K.
,
Ouchi
M.
,
Isobe
Y.
,
Harikane
Y.
,
Zhang
Y.
,
Ono
Y.
,
Umeda
H.
,
Oguri
M.
,
2023
,
ApJ Supplement Series
.
in press, preprint
()

Natarajan
P.
,
Pacucci
F.
,
Ferrara
A.
,
Agarwal
B.
,
Ricarte
A.
,
Zackrisson
E.
,
Cappelluti
N.
,
2017a
,
ApJ
,
838
,
L117

Natarajan
P.
,
Pacucci
F.
,
Ferrara
A.
,
Agarwal
B.
,
Ricarte
A.
,
Zackrisson
E.
,
Cappelluti
N.
,
2017b
,
ApJ
,
838
,
L117

Oesch
P. A.
et al. ,
2023
,
MNRAS
,
525
,
2864

Omukai
K.
,
Schneider
R.
,
Haiman
Z.
,
2008
,
ApJ
,
686
,
L801

Onoue
M.
et al. ,
2023
,
ApJ
,
942
,
L17

Parkinson
H.
,
Cole
S.
,
Helly
J.
,
2008
,
MNRAS
,
383
,
557

Pezzulli
E.
,
Valiante
R.
,
Schneider
R.
,
2016
,
MNRAS
,
458
,
3047

Planck Collaboration
et al. .,
2020
,
A&A
,
641
,
A6

Regan
J. A.
,
Downes
T. P.
,
Volonteri
M.
,
Beckmann
R.
,
Lupi
A.
,
Trebitsch
M.
,
Dubois
Y.
,
2019
,
MNRAS
,
486
,
3892

Regan
J. A.
,
Haiman
Z.
,
Wise
J. H.
,
O’Shea
B. W.
,
Norman
M. L.
,
2020
,
The Open Journal of Astrophysics
,
3
,
1

Reines
A. E.
,
Volonteri
M.
,
2015
,
ApJ
,
813
,
L82

Ricarte
A.
,
Natarajan
P.
,
2018
,
MNRAS
,
481
,
3278

Robertson
B. E.
et al. ,
2023
,
ApJ
,
942
,
L42

Rodighiero
G.
,
Bisigello
L.
,
Iani
E.
,
Marasco
A.
,
Grazian
A.
,
Sinigaglia
F.
,
Cassata
P.
,
Gruppioni
C.
,
2023
,
MNRAS
,
518
,
19

Sanders
R. L.
,
Shapley
A. E.
,
Topping
M. W.
,
Reddy
N. A.
,
Brammer
G. B.
,
2023a
,
preprint
()

Sanders
R. L.
,
Shapley
A. E.
,
Topping
M. W.
,
Reddy
N. A.
,
Brammer
G. B.
,
2023b
,
preprint
()

Sassano
F.
,
Schneider
R.
,
Valiante
R.
,
Inayoshi
K.
,
Chon
S.
,
Omukai
K.
,
Mayer
L.
,
Capelo
P. R.
,
2021
,
MNRAS
,
506
,
613

Sassano
F.
,
Capelo
P. R.
,
Mayer
L.
,
Schneider
R.
,
Valiante
R.
,
2023
,
MNRAS
,
519
,
1837

Scholtz
J.
et al. ,
2023
,
preprint
()

Senchyna
P.
,
Plat
A.
,
Stark
D. P.
,
Rudie
G. C.
,
2023
,
preprint
()

Smith
B. D.
,
Regan
J. A.
,
Downes
T. P.
,
Norman
M. L.
,
O’Shea
B. W.
,
Wise
J. H.
,
2018
,
MNRAS
,
480
,
3762

Spinoso
D.
,
Bonoli
S.
,
Valiante
R.
,
Schneider
R.
,
Izquierdo-Villalba
D.
,
2023
,
MNRAS
,
518
,
4672

Tacchella
S.
et al. ,
2023a
,
 MNRAS
,
522
,
6236

Tacchella
S.
et al. ,
2023b
,
ApJ
,
952
,
74

Trinca
A.
,
Schneider
R.
,
Valiante
R.
,
Graziani
L.
,
Zappacosta
L.
,
Shankar
F.
,
2022
,
MNRAS
,
511
,
616

Trinca
A.
,
Schneider
R.
,
Valiante
R.
,
Graziani
L.
,
Ferrotti
A.
,
Omukai
K.
,
Chon
S.
,
2023a
,
preprint
()

Trinca
A.
,
Schneider
R.
,
Maiolino
R.
,
Valiante
R.
,
Graziani
L.
,
Volonteri
M.
,
2023b
,
MNRAS
,
519
,
4753

Trump
J. R.
et al. ,
2023
,
ApJ
,
945
,
L35

Trussler
J. A. A.
,
Conselice
C. J.
,
Adams
N. J.
,
Maiolino
R.
,
Nakajima
K.
,
Zackrisson
E.
,
Ferreira
L.
,
2023
,
MNRAS
.
in press, preprint
()

Übler
H.
et al. ,
2023
,
A&A
.
in press, preprint
()

Valiante
R.
,
Schneider
R.
,
Salvadori
S.
,
Gallerani
S.
,
2014
,
MNRAS
,
444
,
2442

Valiante
R.
,
Schneider
R.
,
Volonteri
M.
,
Omukai
K.
,
2016
,
MNRAS
,
457
,
3356

Valiante
R.
,
Schneider
R.
,
Zappacosta
L.
,
Graziani
L.
,
Pezzulli
E.
,
Volonteri
M.
,
2018
,
MNRAS
,
476
,
407

Ventura
E. M.
,
Trinca
A.
,
Schneider
R.
,
Graziani
L.
,
Valiante
R.
,
Wyithe
J. S. B.
,
2023
,
MNRAS
.

Visbal
E.
,
Haiman
Z.
,
2018
,
ApJ
,
865
,
L9

Volonteri
M.
,
2012
,
Science
,
337
,
544

Volonteri
M.
,
Natarajan
P.
,
2009
,
MNRAS
,
400
,
1911

Volonteri
M.
,
Stark
D. P.
,
2011
,
MNRAS
,
417
,
2085

Volonteri
M.
,
Silk
J.
,
Dubus
G.
,
2015
,
ApJ
,
804
,
L148

Volonteri
M.
,
Reines
A. E.
,
Atek
H.
,
Stark
D. P.
,
Trebitsch
M.
,
2017
,
ApJ
,
849
,
L155

Volonteri
M.
,
Habouzit
M.
,
Colpi
M.
,
2021
,
Nat. Rev. Phys.
,
3
,
732

Volonteri
M.
,
Habouzit
M.
,
Colpi
M.
,
2023
,
MNRAS
,
521
,
241

Wang
F.
et al. ,
2021
,
ApJ
,
907
,
L1

Whitler
L.
,
Stark
D. P.
,
Endsley
R.
,
Leja
J.
,
Charlot
S.
,
Chevallard
J.
,
2023
,
MNRAS
,
519
,
5859

Yang
J.
et al. ,
2020
,
ApJ
,
897
,
L14

Zavala
J. A.
et al. ,
2023
,
ApJ
,
943
,
L9

This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://dbpia.nl.go.kr/journals/pages/open_access/funder_policies/chorus/standard_publication_model)