<figure><img src="image_25.png" /></figure>

Figure 24: Illustration of the determination of the HOD from spatial clustering measurements. (a) Schematic of an autocorrelation measurement for radio AGN; the large circles represent dark matter halos, the ovals galaxies, and those with orange “jets” indicate radio AGN. The correlation function includes separations between AGN in separate halos (the two-halo term) and those within the same halo (the one-halo term). (b) Measurement of the angular autocorrelation function $\omega(\theta)$ for radio AGN at $z \sim 0.7$ from Petter et al. (2024). The theoretical autocorrelation function for dark matter is shown by the lower dotted curve; the upper dotted curve is a scaling (by the square of the linear bias) of the dark matter $\omega(\theta)$. The orange shaded region shows a fit to a full HOD model, including the one-halo term on small scales. (c) Halo occupation distributions for radio AGN (orange and red shaded regions) illustrating their evolution with redshift, and compared to the HOD for optical quasars (Prada et al., 2023). Each HOD consists of a term for central galaxies plus a power-law function rising to higher halo mass, representing satellites. The HODs clearly demonstrate the preference for higher halo masses for radio AGN compared to optically-selected quasars. Source: adapted from Figs 4 & 9 of Petter et al. (2024).

range of satellite fractions (from <10% up to ~40%) for AGN in typical group-scale haloes (e.g., Richardson et al., 2012; Alam et al., 2021; Comparat et al., 2023; Petter et al., 2024).

A related result is the measurement of clustering as a function of SMBH mass as determined through virial SMBH mass estimators (see §5.1). These measurements suggest that while there is a weak relation between clustering and $L_{AGN}$ as seen in many previous studies, there is a strong, roughly linear correlation between $M_{halo}$ and $M_{BH}$ (e.g., Krumpe et al., 2015, 2023), with no clear dependence on Eddington ratio. This indicates that the primary relationship of AGN to their host halos is via SMBH mass, as may be expected given the relations between $M_{BH}$ and galaxy mass, and galaxy mass and $M_{halo}$. This interpretation is consistent with the picture presented in §4.2 in which AGN are short-lived and largely stochastic events in the lifecycle of their hosts.

An interesting but somewhat puzzling set of measurements relate the strength of clustering with AGN obscuration. In the simplest “unified” model of AGN obscuration, obscuration depends only on small (torus) scale orientation which should introduce no connection between obscuration and large-scale environment. However, measurements using mid-IR selected AGN (initially with *Spitzer*-IRAC, then with larger *WISE* samples) provided evidence that obscured AGN reside in more massive halos than their unobscured counterparts. These studies have used both spatial autocorrelations as well as CMB lensing measurements, with increasingly higher-significance results with larger and more precise data sets (e.g., Hickox et al., 2011; Donoso et al., 2014; DiPompeo et al., 2017a; Mitra et al., 2018; Petter et al., 2023). Some studies of X-ray selected AGN in the low-redshift ($z < 0.5$) have found correspondingly stronger clustering for obscured sources (Powell et al., 2018; Krumpe et al., 2018). These results immediately suggest that processes on scales larger than the unified model “torus” are linked to AGN obscuration, as predicted by some theoretical models of

AGN fuelling (e.g., Blecha et al., 2018). Empirical modelling of various physical scenarios and selection effects suggest that an evolutionary scenario is needed to reproduce the clustering results (DiPompeo et al., 2017b; Whalen et al., 2020); see §5.3.

However, other recent studies using X-ray, and optical, and mid-IR selected AGN have found no significant difference between obscured and unobscured AGN (e.g., Viitanen et al., 2023) or even stronger clustering for the unobscured sources (e.g., Allevato et al., 2014; Córdova Rosado et al., 2024). This indicates a range of potential physical processes and/or substantial selection effects that are, as yet, not well understood, and opens up the possibility of new insights with the next generation of large spectroscopic galaxy and AGN surveys.

### 4.5. The complex connection (or lack thereof) between mergers and AGN activity

A great deal of research has focused on the potential for galaxy mergers to provide the drivers for rapid SMBH growth. This is largely motivated by early observational results from the local Universe suggesting that the most powerful starbursts and AGN are associated with major galaxy mergers, and may be part of a sequence in which rapid obscured star formation is followed by obscured and then unobscured quasar activity (e.g., Sanders et al., 1988). This merger-AGN connection was then extended to a broader theoretical framework to explain the co-evolution of SMBHs and galaxy bulges (e.g., Hopkins et al., 2008) and provided compelling reasons to search for connections between these processes.

#### 4.5.1. Techniques for measuring the merger-AGN connection

In light of the long relevant timescales and complexities of the observational signatures involved, a range of different techniques have been utilized to test the AGN-merger connection:

(a) Perhaps the most widely used and well-known method is to determine the fraction of AGN host galaxies that are
in mergers and compare to “control” galaxies that do not show AGN signatures but are otherwise matched in properties such as redshift, stellar mass, and SFR (e.g., Cisternas et al., 2011b; Villforth et al., 2017). These analyses involve a number of different challenges, including the detectability of tidal features and other merger signatures at low surface brightness (particularly in the presence of a bright nuclear source), difficulties in creating appropriately matched galaxy samples, and statistical constraints, especially for the AGN samples which are often limited in size. Identification of merger signatures can use visual classification, often dividing into merger subclasses such as major and minor mergers and interacting systems. Quantitative techniques include measurements of asymmetry or clumpiness in the galaxy images (e.g., Lotz et al., 2004), or more recently the use of advanced machine learning algorithms (e.g., Pearson et al., 2019).

(b) An alternative method is to consider all galaxies, separated into mergers and non-mergers (or subclasses therein) and to determine the fraction of each class that show AGN signatures or the corresponding average AGN luminosity (e.g., Goulding et al., 2018; Comerford et al., 2024). This approach uses similar techniques to overcome the challenges of identifying mergers and control galaxies, but has the advantage of better statistical precision in utilizing the full galaxy population.

(c) Finally, several studies have expanded on technique (b) by focusing on the AGN fraction in interacting galaxies as a function of their separation distance, thus tracing the incidence of AGN over the course of the merger process (e.g., Ellison et al., 2008; Li et al., 2008; Scott and Kaviraj, 2014; Ellison et al., 2019). These can use spectroscopically identified kinematic pairs, or simply close separations in photometric samples. These studies are sometimes augmented by samples of late-stage mergers that can only be identified on the basis of disturbed morphologies rather than close pairs of distinct galaxies (e.g., Ellison et al., 2013). Recently, the post-merger evolution (i.e., after coalescence) has been explored by using dating of stellar populations in post-starburst galaxies, finding that the AGN activity declines after a peak around coalescence (Ellison et al., 2024).

### 4.5.2. Overview of results for merger fuelling

Using the technique (a) above, a broad range of studies collectively find no clear difference in merger fractions between AGN and “inactive” galaxies that are matched in key properties, leading to the conclusion that the bulk of AGN activity is not associated with mergers (e.g., Cisternas et al., 2011b; Villforth, 2023). There is evidence that mergers become more prevalent at the highest luminosities or for dust-reddened systems (e.g., Treister et al., 2012; Glikman et al., 2015) although these results have been contested (e.g., Villforth, 2023).

Alternatively, techniques (b) and particularly (c) have yielded clear positive associations of mergers with AGN activity. These results are especially pronounced for close separations and/or advanced stages of mergers (e.g., Ellison et al., 2013, 2019).

This is especially true for obscured AGN (identified for example in the mid-IR or X-rays), for which the enhancement seems increasingly clear at all redshifts (e.g., Satyapal et al., 2014; Kocevski et al., 2015; Goulding et al., 2018; Ricci et al., 2021); see also §5.3. The particular connection of obscured AGN activity with mergers has been predicted by a number of theoretical works (e.g., Hopkins et al., 2008; Blecha et al., 2018).

The ultimate conclusion from these works is that some SMBH growth is connected to galaxy mergers and interactions, and it may be especially important for fuelling of the most powerful and rapid growth (e.g., when the AGN luminosity is driven above the “knee” of the luminosity function for a substantial period of time). However, mergers are probably not associated with the bulk of the AGN population and corresponding SMBH growth, which is likely driven by secular processes (see §3). Overall, any physical process that drives gas within the potential of the SMBH is likely to lead to AGN activity.

### 4.5.3. Dual AGN

An exciting new area of investigation into the connections between mergers and AGN activity over the last decade is the identification of dual AGN. While dual AGN (or binaries) likely provide direct insights into SMBHs in post-merger systems, current observational facilities primarily enable the detection of dual AGNs in separate galaxies (especially at high redshift), prior to any potential galaxy merger. Recent discoveries by JWST have confirmed the presence of approximately 30 candidate dual and triple AGNs with separations ranging from a few to tens of kpc, extending out to $z \sim 5$ (e.g., Perna et al., 2023b,a; Li et al., 2024b); these observational constraints can constrain the subgrid modeling of large-scale cosmological simulations (Rosas-Guevara et al., 2019; Volonteri et al., 2022; Chen et al., 2023; Puerto-Sánchez et al., 2025). This is a growing field of research and new techniques such as varstrometry with Gaia data or the Gaia multi-peak method combined with follow-up observations using HST or MUSE adaptive-optics spectroscopy have been successful to confirm multiple sources (Chen et al., 2022b; Mannucci et al., 2022, 2023; Ciurlo et al., 2023; Scialpi et al., 2024). While these systems are crucial for understanding the hierarchical growth of SMBHs through mergers, dual AGNs with separations of less than a kpc are even more interesting for SMBH fuelling. Even more difficult to confirm, observed dual systems with separations of hundreds of parsecs have been detected (e.g., Müller-Sánchez et al., 2015; Koss et al., 2023). Searching for pc or even sub-pc scale binary systems is feasible through observations with very long baseline interferometric (VLBI) facilities, such as the Very Long Baseline Array (e.g., Gabányi et al., 2016) and European VLBI Network (e.g., Venturi et al., 2020). These arrays can achieve angular resolution of $< 0.5$ mas at $\sim$cm wavelengths. This resolution corresponds to $< 1$ pc at $z = 0.1$, allowing VLBI to separate jet or accretion flow emission from two AGN even in extremely close binary orbits, and more modest orbits have already been observed (e.g., Deane et al., 2014). In the future, higher-frequency observations with the Next Generation Event Horizon Telescope (Ayzenberg et al., 2023) could achieve even higher resolution of $\sim 0.15$ µas, corresponding to spatial resolution of $< 13$ pc across
all redshifts.

In addition to directly resolving binary AGN, it is possible to identify them via the analysis of specific emission lines, periodic light curves, shifts in the radial velocity of broad emission lines, or variability (e.g., Comerford et al., 2013; Graham et al., 2015). Capturing the gas dynamics, inflows, and outflows at small scales in these systems, both theoretically and observationally, will be key to understanding AGN fuelling and the role of the galactic environment and mergers (De Rosa et al., 2019, for a review).

Finally, it is worth noting that a brand new probe of dual SMBHs has emerged in the past few years with the detection of the stochastic gravitational wave background from pulsar timing observations (EPTA Collaboration et al., 2023; Agazie et al., 2023; Reardon et al., 2023; Xu et al., 2023). These results provide a fascinating *integrated* constraint on the population of massive binary SMBHs over cosmic volumes; the large amplitude of the stochastic background signal points toward a tension with the most straightforward models of SMBH-galaxy relationships and merger rates (Sesana et al., 2018; EPTA Collaboration et al., 2024).

## 4.6. Black holes and AGN in dwarf galaxies

Alongside the work showing that SMBH growth is prevalent across a wide range of host galaxy properties and large-scale environments, discussed in §4.2–4.5 above, the last decade has seen concerted efforts to find evidence of lower-mass BHs in low-mass galaxies, with stellar masses of $M_{\star} \le 10^{9.5} M_{\odot}$ (hereafter “dwarf galaxies”). Dwarf galaxies are likely to have experienced a less dramatic growth and merger history compared to more massive galaxies. As a result, they may offer indirect constraints on BH formation through “archeology”. Specifically, the mass of their BHs may have remained close to their initial mass at birth, and thus the number of dwarf galaxies hosting a BH (i.e., the occupation fraction) serves as a crucial constraint on the formation efficiency of BHs. Furthermore, probing lower mass BHs (reaching the elusive “intermediate-mass” BH regime of $10–10^5 M_{\odot}$; Greene et al., 2020) and constraining the extent of their growth is vital to understand the assembly of the SMBHs that appear ubiquitous in the higher mass galaxy population. Given their importance—and the availability of powerful new facilities to push into this low-mass regime—the last decade has seen studies of BHs and AGN in dwarf galaxies emerging as a distinct research field.

Currently, there is an ongoing systematic search for the presence of BHs in dwarf galaxies in the local Universe (Greene et al., 2008; Reines and Volonteri, 2015; Reines and Comastri, 2016; Mezcua et al., 2016; Mezcua, 2017; Chilingarian et al., 2018; Greene et al., 2020; Reines et al., 2020). The least massive BHs were found in RGG118, with a mass of $5 \times 10^4 M_{\odot}$ (Baldassare et al., 2015), and in NGC4395 with a mass estimate of $9100–4 \times 10^5 M_{\odot}$ (Woo et al., 2019). The identification and confirmation of these low-mass BHs rely on multi-wavelength observations including optical spectroscopy, observations in X-ray, radio, and IR (e.g., Greene and Ho, 2004; Greene et al., 2008; Dong et al., 2012; Chilingarian et al., 2018; Liu et al.,

2018b; Mezcua et al., 2018, 2019; Kaviraj et al., 2019), and remain challenging due to potential contamination from various sources, such as young starburst regions in low-metallicity environment (reassembling AGN), winds from both supernovae (SN) and massive stars, and emission from X-ray binaries (e.g., Gallo et al., 2008, 2010; Reines et al., 2011; Mezcua et al., 2016; Pardo et al., 2016; She et al., 2017; Mezcua et al., 2018; Fornasini et al., 2018; Schirra et al., 2020). Recent observations revealed that, in the low-mass regime, the presence of BHs is not confined to the centres of galaxies; e.g., Reines et al. (2020) discovered that half a dozen dwarf galaxies were hosting offset radio sources that could be AGN, and faint off-centre AGN signatures were detected in dwarf galaxies through MaNGA IFU observations (Mezcua and Domínguez Sánchez, 2020).

Most constraints indicate an AGN occupation fraction in dwarf galaxies, across all luminosity ranges, of $\le 0.03$, with significant variations spanning more than an order of magnitude (Shi et al., 2008; Aird et al., 2012; Schramm and Silverman, 2013; Reines et al., 2013; Miller et al., 2015; Lemons et al., 2015; Pardo et al., 2016; Aird et al., 2018; Mezcua et al., 2018; Kaviraj et al., 2019; Mezcua and Domínguez Sánchez, 2020; Birchall et al., 2020, 2022). Several studies have derived constraints on the BH occupation fraction, requiring corrections for inactive BHs missed in observational samples. These studies suggest that over half of galaxies with stellar masses of $\le 10^{10} M_{\odot}$ may harbor a central BH (Miller et al., 2015; Trump et al., 2015; Nguyen et al., 2019). In large-scale cosmological simulations, the BH occupation fraction is of almost unity in low-redshift low-mass galaxies. However, results for the AGN occupation fraction vary significantly, spanning several orders of magnitude due to differences in the subgrid modeling of BH seeding, growth, feedback and galaxy processes among simulations (Haidar et al., 2022).

## 4.7. Summary of the key drivers of progress

The results highlighted above demonstrate the exciting progress that has been made in the past decade on the connection between SMBH growth and their host galaxies and environments. Interestingly, the bulk of these advances have been driven by new conceptual insights as well as novel data analysis and theoretical frameworks, rather than the availability of new data. In the *previous* decade, the explosion of extragalactic multiwavelength surveys provided observational foundations but also raised significant puzzles that were only resolved in subsequent years. The progress in this past decade was enabled by a deeper understanding of AGN and galaxy evolution timescales and their impact on observables, along with a more sophisticated handling of selection effects via forward modelling and advanced statistical techniques. Theoretical insights have come from simple “pen-and-paper” analytical thinking as well as adapting these ideas to numerical simulations. Much of this progress happened collectively, with advances catalysed by collaborations, conferences, and any discussions. The upshot is a clearer and more comprehensive SMBH and galaxy evolution framework that has been widely adopted by the field.

Of course, as always with astronomy, new facilities have played a significant role in this progress; new sensitive wide-
field sky surveys (for example with *WISE* and *LOFAR* for obscured and jet-emitting AGN, respectively) have enabled greater statistics for clustering and host galaxy measurements, while *JWST* has provided a number of leaps forward including studies of mergers and dwarf AGN.

Looking forward, exciting progress will be driven by an-
other revolution in data volume, through new sensitive wide-
area survey instruments. In the optical and near-infrared, facil-
ities such as Rubin LSST (Ivezić et al., 2019), *Roman Space
Telescope* (Akeson et al., 2019), *SphereX* (Crill et al., 2020),
and the culmination or initiation of a number of ground-based
spectroscopic surveys (e.g. DESI: DESI Collaboration et al.
2024; SDSS-V: Kollmeier et al. 2017; WEAVE: Jin et al. 2024,
4MOST: de Jong et al. 2019; MOONS: Taylor et al. 2018), will
increase the number of AGN and host galaxy measurements by
an order of magnitude or more. In addition, large sky surveys at
other wavelengths, for example with planned or proposed facil-
ities in the radio (Square Kilometer Array: Braun et al., 2015),
IR (PRIMA: Moullet et al., 2023), and X-rays (NewAthena:
Cruise et al. 2024; AXIS: Reynolds et al. 2023) will further add
to multiwavelength AGN samples and provide vital constraints
on their properties.

Taking full advantage of these huge new data volumes will
require new techniques in the analysis of observations. There
has already been significant progress in machine learning tech-
niques to identify AGN from multiwavelength imaging, spec-
troscopy, and timing (e.g., Geach, 2012; Faisst et al., 2019; Car-
vajal et al., 2023; Savić et al., 2023; Hviding et al., 2024) and
these will become especially valuable as sample sizes increase.
However, the last decade has shown us the power of concep-
tual advances that are driven not just by new larger datasets,
but communities of scientists from different viewpoints coming
together to assess the findings and draw a new consensus.

In addition to large sky surveys, new smaller-field instru-
ments will provide deeper and high-resolution views of AGN
and their hosts. Ongoing studies with *JWST* continue to push
the boundaries in studying the host galaxies and environments
of AGN (e.g., Saxena et al., 2024; Bonaventura et al., 2025),
while the E-ELT (Gilmozzi and Spyrromilio, 2007) and other
future extremely large ground-based telescopes will produce
unprecedented sensitivity and spatial resolution for separating
and characterizing the nuclear and stellar components of galax-
ies. Of particular interest are IFU observations that will yield
spatially-resolved maps of velocities, ionized gas properties,
and stellar populations to increasingly high redshifts.

Finally, our understanding of SMBH populations, and in par-
ticular mergers, will be revolutionized with the advent of space-
based gravitational wave astronomy with *LISA* (Amaro-Seoane
et al., 2023; Colpi et al., 2024). *LISA* will build on the integrated
constraints from PTAs (which probe relatively widely separate
SMBH binaries) by directly observing the *coalescence* of indi-
vidual merging-SMBH systems. This will yield corresponding
insights into AGN fuelling via the rates of galaxy and SMBH
mergers measured across all of cosmic time.

5. What fuels the rapid growth of the most massive (and also the first) black holes?

The previous section focused on the demographics and prop-
erties of the overall AGN population, building on the prior sec-
tion which looked at how gas is accreted onto the SMBH. In
this section we take a more detailed look at the most rapidly
growing black holes, systems typically referred to as quasars,
where the accretion physics and host-galaxy/SMBH environ-
ments may be different to the more general AGN population.

Quasars were originally identified as systems where the AGN strongly outshines the host galaxy at optical wavelengths, although more recently quasars are often defined as an AGN with a luminosity above a set threshold, usually the knee in the luminosity function (e.g., a bolometric luminosity of > 10⁴⁵ erg s⁻¹, which broadly corresponds to Lₓ > 10⁴⁴ erg s⁻¹ and an absolute magnitude of M_B <-23). While adopting a specific quasar definition is essential for individual studies, since there is not a single universal definition of a quasar, it is more problematic to stick rigidly to a single definition in a review article without leaving out many important studies. Therefore, in this review we adopt a pragmatic quasar definition: those systems in which the authors referred to them as quasars.

The majority of quasar research has focused on the optically
bright unobscured subset of the population, where at least a
component of the accretion disk is seen in the optical wave-
band, allowing for direct constraints on the mass-accretion and
SMBH properties: we will mostly focus on this population
here. Over the last decade spectroscopic samples of unobscured
quasars have expanded by an order of magnitude from ≈ 100 k
(Schneider et al., 2010) to of-order a million (e.g., Lyke et al.,
2020; Flesch, 2023), providing sufficient source statistics for
more-nuanced SMBH accretion investigations: over the next 2
years the DESI quasar survey alone will increase the quasar de-
mographics to over ≈ 3 M, pushing down to lower-luminosity
systems (e.g., Alexander et al., 2023; Chaussidon et al., 2023).
Large advances have also been made in the identification of the
highest-redshift quasars: a decade ago the first *z* ≈ 7 quasar
had just been identified (Mortlock et al., 2011) while dedicated
deep searches over the last decade have revealed an order of
magnitude more *z* > 7 quasars (Fan et al., 2023). The recently
launched *JWST* telescope has now identified lower luminosity,
but less secure, candidate AGN out to *z* ≈ 10 (e.g., Adamo et al.,
2024a). When possible we will also refer to the (theoretically)
equally populous obscured quasar population but note that, due
to challenges in spectroscopically identifying obscured quasars,
their samples number in the thousands.<fnref n="25" />

<footnote n="25"><sup>25</sup>There are two key challenges in spectroscopically identifying obscured quasars (1) the difficulty in reliably selecting a complete sample of obscured quasars since the optical-near-IR emission will (typically) be dominated by the host galaxy), and (2) their optical faintness in the optical waveband compared to unobscured quasars. We expect this inequity to be addressed over the coming decade with sensitive large-scale multi-object spectrographs (DESI; WEAVE; 4MOST) with surveys targeting X-ray, infrared, and radio AGN (e.g., Smith et al., 2016; Merloni et al., 2019; Duncan et al., 2023; Andonie et al., 2025) to construct a more complete census of the AGN and quasar population; see §4.1 for the merits of different multi-wavelength AGN selection approaches.</footnote>
<figure><img src="image_26.png" /></figure>

Figure 25: Bolometric luminosity (top) and black-hole mass (bottom) versus spectroscopic redshift for quasars from SDSS DR14 (grey dots) and dedicated $z > 4$ quasar identification studies (dark blue circles), compared to candidate high-redshift AGN identified by JWST-NIRSpec with either broad Balmer emission lines (magenta squares) or just narrow emission lines (green triangles); the three highest-redshift JWST AGN plotted in Fig. 26 are labelled. The bolometric luminosities and black-hole masses of three nearby moderately luminous to luminous AGN (NGC 1068, Mkn 231, PDS 456) are plotted to provide context; PDS 456 is plotted for two different black-hole mass measurements, from virial scaling relationships (Reeves et al., 2009) and from GRAVITY+ spatially resolved measurements (GRAVITY Collaboration et al., 2024). Primary sources: Glikman et al. (2010), Willott et al. (2010), McGreer et al. (2013, 2018), Giallongo et al. (2015), Yang et al. (2017b, 2019b, 2021), Matsuoka et al. (2018b, 2022), Onoue et al. (2019), Shen et al. (2019), Rakshit et al. (2020), Fan et al. (2023), Goulding et al. (2023), Greene et al. (2023), Harikane et al. (2023), Kocevski et al. (2023, 2024), Larson et al. (2023), Maiolino et al. (2023, 2024b), Scholtz et al. (2023). Bogdán et al. (2024), Furtak et al. (2024), Juodžbalis et al. (2024), Kokorev et al. (2024), Lai et al. (2024), Lin et al. (2024), Matthee et al. (2024), Mazzolari et al. (2024b), Napolitano et al. (2024), Perger et al. (2024), Takahashi et al. (2024), Treiber et al. (2024), Tripodi et al. (2024), Juodžbalis et al. (2025), Taylor et al. (2025).

Fig. 25 shows the bolometric luminosity–redshift and SMBH mass–redshift plane for spectroscopically identified quasars and high-redshift JWST-identified AGN. Dedicated selection and follow-up studies of $z > 4$ quasars has allowed for the identification of “faint” quasars out to $z \approx 7$ with properties which overlap the parameter space explored by JWST. However, thanks to the high spectroscopic sensitivity of NIRSpec, JWST is now starting to identify the majority of the AGN population (lower luminosities; less-massive SMBHs) out to $z > 8$.

In the following sub sections we first discuss approaches and progress in measuring the masses of SMBHs in distant quasars, a key parameter to determine the relative SMBH growth rate (§5.1) before investigating the physics of accretion and the surprising overall uniformity in the central-engine properties of unobscured quasars out to $z \approx 7$ (§5.2). We then review the evidence for different quasar populations representing distinct phases in the overall lifetime of quasars, investigating the well-known but still largely unconfirmed quasar-evolutionary model (§5.3). We finally look at the advances in identifying and measuring the properties of the highest-redshift quasars and AGN ($z \approx 4–10$) over the last decade, systems which can place constraints on the origins of the first seed black holes (§5.4), including the latest remarkable discoveries made by JWST.

## 5.1. Black-hole mass measurements of distant quasars

A key advantage in the identification of unobscured quasars is that some of the physical properties of the accretion process can be inferred due to the direct detection of emission from

the accretion disk (see §4 of Padovani et al., 2017, for a brief review). We highlight inferred because the accretion-disk emission from quasars typically peaks in the far-UV ($\approx 10–100$ nm; e.g., Jin et al., 2012; Mitchell et al., 2023) where absorption from gas within our own Galaxy precludes direct measurements of the peak of the accretion-disk emission. Consequently, as discussed in §4.1 of AH12, it is necessary to apply uncertain bolometric corrections (e.g., Duras et al., 2020; Shen et al., 2020) to the observed quasar emission to estimate the overall accretion-disk luminosity and infer the mass-accretion rate, which will be dependent on the accretion-disk structure and Eddington ratio (e.g., Hagen et al., 2024; Kang et al., 2024, see §3.4). We note, in the absence of obscuration and host-galaxy dilution affects, the most accurate bolometric corrections will be at wavelengths closest to the accretion-disk emission peak.

For the gas in the vicinity of the SMBH ($< 10$ pc for a SMBH of $10^8 M_{\odot}$; see Footnote 1 and Fig. 1), the kinematics are dominated by the SMBH rather than the host galaxy, providing the potential to estimate the mass of the SMBH. For nearby systems where the spatial resolution is high and the targets are bright, it is possible to directly measure the gas kinematics as a function of distance from the SMBH, allowing for an accurate model of the gravitational potential and hence the SMBH mass (e.g., Kormendy and Richstone, 1995; Kormendy and Ho, 2013, for reviews). However, for the majority of the more distant quasars, we need to resort to mass-scaling relations, which provide an estimate of the SMBH mass based on the widths of the broad emission lines and an estimate of the distance to the broad-line
region (BLR; see Shen, 2013, for a review). These mass-scaling relations are referred to as virial SMBH mass estimators as they assume the BLR gas is in virial equilibrium with the SMBH (see Eqn. 8 from AH12). As the BLR radius is typically too small to be directly measured from imaging data, a technique called “reverberation mapping” is used to measure the time lag between accretion disk and broad-line luminosity changes (e.g., Blandford and McKee, 1982; Peterson, 1993; Shen, 2013), which provides a direct measurement of the accretion disk–BLR distance.

The foundations of the virial SMBH mass technique were established in the decade prior to AH12 (e.g., Peterson et al., 2004; Bentz et al., 2009) from reverberation campaigns of tens of local AGN, where a tight relationship was found between the luminosity of the AGN and the radius of the BLR.<fnref n="26" /> This important discovery showed that SMBH masses can be estimated from single-epoch optical spectroscopic observations with individual mass uncertainties of ≈ 0.5 dex (Vestergaard and Peterson, 2006; Shen, 2013), an aspect exploited to estimate SMBH masses for 100,000s of distant quasars (e.g., Shen et al., 2011; Rakshit et al., 2020). Prior virial SMBH mass estimators had been determined for local AGN, rather than distant luminous quasars, and are typically determined for the Hβ emission line, which becomes inaccessible for optical spectroscopy at z > 0.7 and, instead, the more-poorly calibrated Mg II and/or C IV broad emission lines typically need to be adopted for distant quasars (e.g., Shen and Liu, 2012; Coatman et al., 2016).<fnref n="27" />

Reverberation mapping (RM) of distant quasars is more challenging than for local AGN because their larger SMBH masses mean the BLR response time to accretion disk changes will be of order years rather than the days/month timescale for lower-luminosity AGN, and their fainter optical magnitudes necessitate longer exposure times. So the previous approach of undertaking RM on a source by source basis, “sequentially”, is too inefficient. The catalyst to allow for significant advances has been the advent of large field of view multi-object spectrographs over the last decade which have allowed for RM for a significant number of distant quasars to be observed in “parallel” as part of the larger SDSS and DES survey programmes (e.g., SDSS-RM and OzDES-RM; Shen et al., 2015; Yuan et al., 2015; Malik et al., 2023; Shen et al., 2024). These RM campaigns have found flatter luminosity-radius relationships for distant quasars when tracing broad Mg II than that found for local AGN using broad Hβ (Homayouni et al., 2020; Yu et al., 2023), but they remain broadly consistent with a similar overall scatter. Microlensing of the BLR of distant gravitationally lensed quasars provides a complementary approach to measuring the BLR size. While the sample sizes are small, studies have found results broadly consistent with those obtained from RM, once biases have been taken into account (Guerras et al., 2013; Hutsemékers et al., 2024). These key measurements provide reassurance on the reliability of prior SMBH mass estimations.

The vast majority of distant quasar SMBH mass estimates have relied on the width of broad Mg-II as the proxy for the gas kinematics, which is a less reliable mass tracer than the Balmer broad lines (Shen et al., 2019). However, over the last decade the increase in sensitive near-IR spectroscopic facilities (e.g., Gemini-GNIRS; VLT-XShooter; JWST-NIRSpec) have allowed for SMBH masses of distant quasars based on the more accurately calibrated Hβ emission line and the empirical relation between the BLR size and rest-frame 5100 Å continuum luminosity (Vestergaard and Peterson, 2006); i.e., the same tracers employed for low-redshift Universe studies. Initial findings indicate that Hβ-based mass estimates are, on average, consistent with, or slightly lower than, Mg-II-based masses (Ding et al., 2023; Yang et al., 2023c; Bertemes et al., 2024; Yue et al., 2024c).

While it is currently impossible to directly measure the BLR radius for the majority of distant quasars due to the need for µas resolution, the recent GRAVITY+ instrument on the VLT-I interferometer (GRAVITY+ Collaboration et al., 2022) provides the potential to measure the spatial offset between the blue and red side of the broad emission lines for the brightest quasars hosting the most massive SMBHs. For example, Abuter et al. (2024) reported the dynamical measurement of the SMBH mass for a *z* ≈ 2 quasar with GRAVITY+ finding it was 0.7 dex and 0.4 dex lower than the estimate based on single-epoch Hβ and Hα, respectively; see also Fig. 25 for virial mass–GRAVITY+ differences for the PDS 456 quasar. The GRAVITY+ wide project (Drescher et al., 2022) promises to extend this study to hundreds of *z* ≈ 2 quasars, providing important calibrations of the virial SMBH estimators for distant quasars. This will also allow any potential cross calibration issues due to differences in the radius measurement approaches (i.e., responsivity weighted versus emissivity weighted) to be addressed.

## 5.2. Accretion environment of distant quasars

The accretion disk is the origin of the high luminosities of quasars and, while neutral hydrogen within our own galaxy prevents us from directly detecting the peak of the emission at far-UV wavelengths, the near-UV-optical waveband traces the cooler Rayleigh-Jeans tail of the accretion-disk emission and provides a plethora of emission lines which yield indirect constraints on the accretion-disk properties (e.g., Mathews and Ferland, 1987; Osterbrock and Ferland, 2006). The emission from quasars at wavelengths shorter and longer than the broad far-UV-optical waveband is typically due to processes associated with the broader accretion-disk environment: thermal IR emission from the dusty obscuring torus heated by the accretion disk, synchrotron radio emission from an AGN-driven jet/wind/shocks, X-ray emission due to inverse Compton scattering from electrons in a high-energy accretion-disk “corona”, in addition to host-galaxy emission (stellar emission; star formation); see Figs. 1 and 2. We may therefore also expect a clear connection between this broader accretion disk environment and the accretion disk itself.

Remarkably, the properties of the central engine of optically selected quasars appear to have changed little with cosmic time

<footnote n="26">²⁶It is important to note that these relationships are defined for populations and the BLR properties (e.g., geometry; orientation; ionisation; metallicity) are captured as a single “fudge factor” (*f*) which can differ significantly between individual quasars.</footnote>

<footnote n="27">²⁷Note that C IV-derived SMBH masses require a blueshift-dependent correction to account for the outflow component (Coatman et al., 2017); see §5.2.</footnote>
<figure><img src="image_27.png" /></figure>

Figure 26: (left): rest-frame UV spectral composite of $z = 6.3–7.6$ quasars (violet) from Yang et al. (2021) compared to the rest-frame UV-optical spectral composite of $z \approx 1.5$ normal quasars (blue) from Fawcett et al. (2022) and JWST-NIRSpec spectra of three $z \approx 10$ AGN candidates: GNZ-11 (orange; Maiolino et al., 2024b), GHZ9 (green; Napolitano et al., 2024) UHZ-1 (red; Goulding et al., 2023; Bogdán et al., 2024); (right): rest-frame UV-optical spectral composite of normal quasars (blue) and red quasars (salmon; red) at $z \approx 1.5$ from Fawcett et al. (2022, 2023) compared to the spectral composite of JWST-identified “little red dots” (violet) at $z = 2–6$ from Setton et al. (2024). Key emission and absorption lines are indicated and labelled.

out to at least $z \approx 7$, despite the huge changes in cosmic environment over this broad time frame.<fnref n="28" /> The composite rest-frame UV-optical spectrum of $z \approx 6–7$ quasars is largely indistinguishable from that of lower-redshift quasars (see Fig.26), at least when matched in luminosity, with essentially identical continuum emission slopes and emission-line properties (i.e., emission line strengths and widths) indicating no significant changes in the overall accretion-disk emission (e.g., Shen et al., 2019; Yang et al., 2021; D’Odorico et al., 2023). This result indicates that typical quasars already have “mature” metal-rich central regions even out to $z \approx 7$ (Lai et al., 2022). The average X-ray properties and the relationship between the soft X-ray and near-UV emission (referred to as the $\alpha_{OX}$ relationship, the equivalent power-law slope between the rest-frame UV emission at 2500 Å and 2 keV in the X-rays), which traces the relationship between the accretion disk and the “corona”, and the 2500 Å luminosity of quasars, also appears to be largely unchanged with redshift (e.g., Salvestrini et al., 2019; Vito et al., 2022; Wang et al., 2021a; Rankine et al., 2024).<fnref n="29" /> The overall UV-mid-IR SED of quasars also lacks strong redshift evolution, indicating no broad changes in the relationship between the accretion-disk emission and the hot dust at the inner edge of the obscuring “torus” (e.g., Leipski et al., 2014; Lyu and Rieke, 2022; Auge et al., 2023). In terms of host-galaxy properties, careful subtraction of the quasar emission in rest-frame optical-near-IR imaging (greatly facilitated by JWST observations), has shown that quasars out to $z \approx 4$ are typically found in massive galaxies ($M_{stellar} > 10^{10} M_{\odot}$), and the consequent SMBH-host galaxy mass relationship shows no strong evolution with redshift (e.g., Suh et al., 2020; Mountrichas, 2023; Tanaka et al., 2024; Sun et al., 2024), in contrast to the JWST-identified high-redshift AGN (see §5.4.2). The average far-IR emission from quasars is found to increase with redshift, tracking the large-scale changes in dust-obscured star formation with cosmic time for typical massive star-forming galaxies (e.g., Schreiber et al., 2015; Stanley et al., 2017; Schulze et al., 2019; Calistro Rivera et al., 2021) and implying that the majority of this emission is from the host galaxy rather than from the central engine; see §4.2 and §5.4.2 for further host galaxy and environmental constraints.

The essence of these results had been established in the decade prior to AH12, although the increase in data over the last decade have strengthened these conclusions. The fractions of different quasar sub populations also appear to be consistent, with broadly similar fractions of radio-loud (see §4.1.4) and broad-absorption-line quasars (BALQSOs)<fnref n="30" /> identified out to $z \approx 6$, at least for optically selected quasars (e.g., Bañados et al., 2018; Liu et al., 2021; Yang et al., 2021; Bischetti et al., 2023); we note there is evidence for an increase in the fraction of BALQSOs at $z > 6$, not clearly driven by quasar luminosity or accretion rate, suggesting an increase in accretion-driven winds at early epochs.

To first order, the lack of change in the overall properties of quasars with cosmic time gives the impression of little diversity in the central engine. However, this is not the case and significant variation is seen between individual systems: indeed, a major focus of quasar research over the last few decades has

<footnote n="28">²⁸We note that most studies have used SDSS-selected quasars which have a specific selection towards blue quasars and so a lack of differences may be partially a selection effect. For example, see §5.4 for discussion of the significant differences in the properties of the high-redshift candidate AGN uniquely discovered by JWST.</footnote>

<footnote n="29">²⁹The lack of a strong redshift dependence in the $\alpha_{OX}-2500$Å luminosity relationship allows carefully selected quasar sub samples to be used as standard candles out to at least $z \approx 3$, far exceeding the most distant Type 1a supernovae (e.g., Risaliti and Lusso, 2019; Lusso et al., 2020).</footnote>

<footnote n="30">³⁰BALQSOs are quasars where BALs blueward of the broad-emission lines are detected in the rest-frame UV (e.g., Weymann et al., 1981, 1991). The BALs occur when high-velocity outflowing gas passes across the line of sight between the quasar accretion disk and the observer, causing broad and complex absorption features related to the velocities of the individual absorbing clouds. In the majority of BALQSOs, the BALs are only associated with high-ionisation emission lines such as C IV, and the system is typically called a HiBAL, but a smaller fraction of BALQSOs have additional BALs associated with low-ionisation lines such as Mg II, in which case the system is often referred to as a LoBAL (e.g., Voit et al., 1993; Hall et al., 2002).</footnote>
<figure><img src="image_28.png" /><figcaption>Figure 27: Eddington-ratio–black-hole mass plane for quasars at $z \approx 2$ from the SDSS survey DR16, colour coded to indicate from left to right: C IV blueshift, He II EW, $\alpha_{OX}$, and relative excess of hot-dust at 2 µm. A key region of the left-hand panel highlights the largest C IV blueshifts, which is repeated in the other panels and shows the connection between C IV outflowing gas, the accretion-disk SED (as traced by He II EW and $\alpha_{OX}$), and the connection to excess hot dust. Source: adapted from (Temple et al., 2023) but recalculated and reprojected onto the $L/L_{Edd} - M_{BH}$ plane. The right-hand panel calculates the hot-dust excess following Temple et al. (2021) but using the (Temple et al., 2023) sample and projecting the results onto the $L/L_{Edd} - M_{BH}$ plane (M. Temple, priv. communication).</figcaption></figure>

been identifying the driver(s) of the diversity seen in the spectra of quasars. While several significant empirical correlations have been identified between quasar emission lines and continua (e.g., Baldwin, 1977; Boroson and Green, 1992; Sulentic et al., 2000; Richards et al., 2011; Shen and Ho, 2014; Marziani et al., 2018; Panda, 2024), the greatest progress in understanding the physical origin of these relationships has come in the last decade thanks to the order of magnitude increase in the number of sensitive optical spectra for quasars, allied to developments in data-analysis techniques, broader overall sensitive wavelength coverage, and accretion-disk models.

Foundational early quasar research had previously identified the C IVλ1549 broad emission line as a key spectral diagnostic, given the anti correlation of the C IV equivalent width with the quasar UV continuum luminosity (the famous “Baldwin” effect; Baldwin, 1977). More recently, a connection has been found between the C IV equivalent width (EW) and the blueshift of C IV (e.g., Richards et al., 2011; Rivera et al., 2022), indicating that a significant wind component is present in broad C IV that is not seen in the lower-ionisation broad emission lines. It was found early on that Eddington ratio must be a driver (e.g., Baskin and Laor, 2005; Shen and Ho, 2014; Sulentic and Marziani, 2015) but that property alone could not explain the observed relationships. Greater clarity was made when considering the strength of He IIλ1640, a sensitive indicator of the unseen far-UV radiation field and, hence, the accretion-disk SED (Rankine et al., 2020).<fnref n="31" /> Investigating the He IIλ1640 and C IVλ1549 properties of quasars across the Eddington-ratio–$M_{BH}$ plane, Temple et al. (2023) showed that the C IV blueshift is tightly connected to both Eddington ratio and black-hole mass with a distinct region in the Eddington

ratio–$M_{BH}$ plane where the strongest outflows are found and He II is weak. This key result is shown in the first two panels of Fig. 27. The weakness of He II indicates a relatively “soft” accretion-disk SED which prevents the gas from becoming over ionised, allowing for the launch of a line-driven wind through radiation pressure (e.g., Murray et al., 1995; Proga et al., 2000; Giustini and Proga, 2019) and producing blue-shifted C IV emission (Yong et al., 2017; Matthews et al., 2023): quasars with a lower SMBH mass or lower Eddington ratio either overionise the gas or produce insufficient photons to radiatively drive the wind. As can be seen in Fig. 27, the same region in the Eddington ratio–$M_{BH}$ plane where strong C IV blueshifts and weak He II is found also corresponds to relatively weak X-ray emission (as traced by $\alpha_{OX}$), providing further verification of a “soft” accretion-disk SED; indeed, Timlin et al. (2021) recently discovered a tight relationship between $\alpha_{OX}$ and He II equivalent width which demonstrates that $\alpha_{OX}$ provides an indirect indicator of the accretion-disk SED peak wavelength. Overall, these results are in good qualitative agreement with physical accretion-disk models which predict the changes in the accretion-disk SED, X-ray corona, and the strength of line-driven winds with Eddington ratio and $M_{BH}$ mass (e.g., Kubota and Done, 2018; Giustini and Proga, 2019); see §3.4.

Approximately 10–40% of quasars are referred to as BALQSOs (e.g., Allen et al., 2011), quasars with powerful outflows observed in absorption (see Footnote 30). An area of considerable uncertainty since their original discovery has been the relationship between BALQSOs and normal “non-BAL” quasars. Thanks to a suite of complementary multi-wavelength and time-domain studies over the last decade there is growing evidence that most BALQSOs are normal quasars with “soft” accretion-disk SEDs (either intrinsically “soft” or obscuration absorbs the hard ionising photons) allowing for the launch of powerful line-driven winds from the accretion disk even as far out as the dusty “torus” (e.g., Czerny and Hryniewicz, 2011; Giustini and Proga,

<footnote n="31">³¹He IIλ1640 provides a sensitive indicator of the > 54 eV photon rate and, hence, the far-UV continuum. He II is produced from the recombination of He III to He II and so is driven by the ionisation rate of He II to He III (e.g., Timlin et al., 2021).</footnote>
2019). Orientation angle almost certainly plays a role given the connection between dust and BALQSO activity, particularly for LoBALs (e.g., Hall et al., 2002; Gibson et al., 2009; Urrutia et al., 2009; Matthews et al., 2016), and the fact that the BAL wind is unlikely to have a 100% covering factor. Key evidence includes:

(a) The fraction of BALQSOs within the quasar population increases with decreasing He II equivalent width, as does the strength of the BAL outflow and dust reddening, with LoBALs (BALQSOs with the most powerful and optically thick outflows) having the weakest He II emission (e.g., Baskin et al., 2013; Hamann et al., 2019). Indeed, the optical spectral properties for the majority subset (HiBALs) are consistent with those of normal quasars when the emission lines are reconstructed for the impact of the BAL features, with the exception that BALQSOs are not found with strong He II, large C IV EW, no C IV blueshift (e.g., Baskin et al., 2013; Rankine et al., 2020), as expected given Fig. 27.

(b) BALQSOs are often X-ray weak, producing less X-ray emission than expected from the $\alpha_{OX}-L_{2500}$ relationship, particularly for the most extreme LoBALs (e.g., Gallagher et al., 2006; Luo et al., 2014). The origin of the X-ray weakness often appears to be due to absorption, although some may be intrinsically X-ray weak where the “corona” is quenched or unable to form (e.g., Proga, 2005; Leighly et al., 2007; Luo et al., 2014; Liu et al., 2018a).<fnref n="32" /> These line-driven winds may also be a major contributor to the X-ray obscuration (Matthews et al., 2016, 2020).

(c) BAL features are found to vary on day-year timescales with the strongest BALs demonstrating the weakest variation, likely due to (predominantly) transverse motions of a clumpy absorbing gas (e.g., Capellupo et al., 2011; Filiz Ak et al., 2013). In extreme cases the BAL features can disappear altogether and the BALQSO becomes a normal quasar and vice versa (i.e., BAL features can appear in normal quasars; Rogerson et al., 2016) with even LoBALs transforming into HiBALs (e.g., Filiz Ak et al., 2013; McGraw et al., 2017; Yi and Timlin, 2021; Yi et al., 2022). The properties of a quasar following a BALQSO transition appear indistinguishable to those of normal quasars (Sameer et al., 2019), although the current samples are small. Based on the fraction of transiting quasars (≈ 2%) and the probed rest-frame timescales (≈ years) the typical BAL lifetime is calculated to be ≈ 100-1000 years.

We highlight that differences may occur for the rare FeLoBAL population, where detailed spectral analysis-physical modelling has shown that the launching radius is

sometimes on kpc scales rather than in the vicinity of the accretion disk, and where some FeLoBALs appear to accrete at low accretion rates (Choi et al., 2022a,b; Leighly et al., 2022).

A subset of the non-BAL quasar population with weak high-ionisation emission lines (C IV EW < 15Å; weak-line quasars: WLQs) but normal low-ionisation emission lines (e.g., Wu et al., 2011; Plotkin et al., 2015; Chen et al., 2024) are also often X-ray weak (e.g., Pu et al., 2020; Wang et al., 2022). The X-ray weak fraction (factor > 6 less X-rays than predicted by the $\alpha_{OX}-L_{2500}$ relationship) appears to increase with decreasing C IV EW, approaching ≈ 40% as compared to ≈ 6% for normal quasars (e.g., Tanimoto et al., 2020; Pu et al., 2020). The origin of the X-ray weakness appears to be similar to that speculated for X-ray weak BALQSOs (e.g., Luo et al., 2015; Liu et al., 2021; Ni et al., 2022a; Wang et al., 2022); see Footnote b.

Quasars with strong C IV blueshifts have also been found to exhibit excess near-IR emission, suggesting that systems with powerful winds have a larger covering factor of hot dust (Temple et al., 2021). In Fig. 27 (right panel) this result is extended to the full Eddington ratio-$M_{BH}$ plane of Temple et al. (2023) showing that the quasars with the strongest outflows have both the weakest He II and excess near-IR emission, which may be due to winds produced in the outer part of the accretion disk driving dust away from the equatorial regions and into the polar regions, such as that observed for some local AGN (see §3.3 and Fig. 10). Differences in the accretion-disk SED across the Eddington ratio-$M_{BH}$ plane may also help to explain the modest fraction (≈ 15%; Lyu et al., 2017) of quasars with relatively weak hot-dust emission, a population that appear to lack a polar-dust component (e.g., Lyu and Rieke, 2018, 2022).

Good advances in our understanding of the radio emission from quasars have been made over the last decade due to sensitive large field-of-view radio facilities such as LOFAR, allowing for the detection of the weaker “radio-quiet” emission from the majority of the quasar population. Using deep LOFAR observations, Rankine et al. (2021) has shown that radio-detected quasars reside in the same region of the C IV EW–blue shift plane as radio-undetected quasars. The radio-detection fraction of quasars increases with C IV blueshift but the radio-loud fraction decreases with increasing C IV blueshift such that the majority of radio-loud quasars have no C IV blueshift. The radio-loud quasars also have the strongest He II emission, indicating a strong ionising SED preventing a line-driven wind but allowing a powerful radio-emitting jet to form. BALQSOs and dust-reddened quasars follow this same basic trend, although the presence of a BAL or dust reddening leads to enhanced radio emission (Petley et al., 2024). Dust is a stronger driver of enhanced radio emission than either the C IV properties or BAL features, with a striking relationship found between the amount of dust reddening towards the quasar and the radio-detection fraction (Fawcett et al., 2023; Calistro Rivera et al., 2023; Petley et al., 2024). This radio–opacity connection (the dust is probably a proxy for general opacity, given dust–gas mixing in the ISM) is likely due to shocks caused by quasar-driven winds or low-power jets interacting with the circumnuclear/ISM environment (e.g., Zakamska and Greene, 2014; Nims et al., 2015). The identification of excess mid-IR emission from red quasars

<footnote n="32">³²The absorber is speculated to reside between the accretion disk and the BLR (possibly a geometrically and optically thick inner accretion disk, as expected at the highest Eddington ratios; e.g., Sadowski et al., 2014; Jiang et al., 2019.; see Fig. 11). This absorber further prevents the gas from becoming over ionised and facilitating even more powerful line-driven winds than those seen in normal quasars; see Higginbottom et al. (2024) for a recent radiative hydrodynamical simulation demonstrating the significant radiative challenges in not over ionising the gas</footnote>
suggests the dust may be circumnuclear (potentially polar dust; see Fig. 10) and entrained in a quasar-driven wind (Calistro Rivera et al., 2021; Zakamska and Alexandroff, 2023).

### 5.3. Evidence for a quasar evolutionary sequence

A key driver of quasar research over the last few decades has been testing if different quasar populations represent distinct phases within a broader evolutionary sequence. The catalyst in the now “classical” quasar evolutionary model is the merger of two gas-rich galaxies, which drives gas into the central regions and initiates intense dust-obscured star formation (e.g., Sanders et al., 1988; Canalizo and Stockton, 2001; Hopkins et al., 2008). The star formation serves to enhance the gas inflow towards the central SMBH, resulting in quasar activity. The quasar will be initially heavily obscured from both the star formation and the nuclear region. However, due to the interaction between quasar-driven winds (and potentially jets) with the gas and dust (see §6.3), the quasar will ultimately drive away the obscuring material, shut down the star formation and reveal an unobscured quasar. We note that the merger only serves to drive gas into the central ≈ kpc of the galaxy (e.g., Hopkins and Quataert, 2010), and other processes (e.g., gas torqued by stars; stellar winds; see §3.1) are required to drive the gas towards the nucleus. Consequently, given sufficient gas inflow, this evolutionary sequence could also be initiated by secular processes, either from within the galaxy or from interactions with a satellite galaxy. See Fig. 28 for a schematic representation of expected changes on the galactic and nuclear scale.

We have already discussed in the previous sub section and §3.5 that AGN and quasars are not static objects and can vary substantially during their mass-accretion “event” lifetimes, causing changes to their observed properties (e.g., spectral type; outflow signatures; hot-dust emission). Therefore, how can we distinguish between normal “expected” time-variable behaviour and distinct phases within a broader evolutionary sequence? A key difference is that the evolutionary model predicts clear trends over long timescales (e.g., a decrease in SFR with a decrease in nuclear obscuration) while the standard unified AGN model (essentially the “null hypothesis” for evolutionary model tests) predicts no significant long-term trends. The most reliable test of an evolutionary model will, therefore, be those based on studying different quasar (and galaxy) populations rather than individual sources. Inherent to any population studies, the selection approaches of the quasar populations must be self consistent, including a “control” sample selected in the same basic manner, to ensure any differences are not due to selection biases.³³³⁴ This can be challenging since

many quasar populations are rare and faint, particularly when obscured and/or distant. This has often led to research focused on small often heterogenous quasar populations, selected using an approach that differs to that adopted for the control sample, which can lead to apparent conflicting results in the literature on differences (or not) between quasar populations.

Over the last decade a growing number of studies, adopting a range of approaches but utilising well-defined control samples, have found convincing evidence for obscured AGN and quasars residing in statistically different host-galaxy environments to unobscured systems. Broadly, a larger fraction of AGN appear to be obscured in galaxy mergers than found for AGN in isolated galaxies (e.g., Kocevski et al., 2015; Del Moro et al., 2016; Donley et al., 2018; Goulding et al., 2018; Ellison et al., 2019); see also §4.5. A trend with increasing obscuration has also been found as a function of merger phase and companion galaxy separation. Using sensitive NuSTAR high-energy observations of nearby luminous and ultra-luminous infrared galaxies, Ricci et al. (2017a, 2021) found a significant increase in the fraction of obscured AGN in late merger phases (when the galaxies are close to coalescence; i.e., < 10 kpc), when compared to early merger phases.<fnref n="35" /> Adopting high-resolution optical/near-IR imaging and the high-energy near obscuration unbiased Swift-BAT sample, Koss et al. (2018) showed that the obscured AGN fraction increased yet further for the closest scales of < 3 kpc. Qualitatively similar results are also found for mid-IR selected AGN samples, again a near obscuration-unbiased selection approach, where the obscured AGN fraction has been found to systematically increase as a function of galaxy separation distance down to host-galaxy scales of at least ≈ 10 kpc (e.g., Satyapal et al., 2017; Weston et al., 2017; Barrows et al., 2023; Dougherty et al., 2024).

Numerical simulations predict that a large fraction of the obscuration towards luminous AGN and quasars is due to intense on-going dust-obscured star formation, which can be enhanced in gas-rich galaxy mergers (e.g., Hopkins et al., 2008; Blecha et al., 2018). Observationally, it has been shown that the host galaxy can contribute to a substantial fraction of the obscuration towards the AGN (e.g., Buchner and Bauer, 2017; Circosta et al., 2019; Gilli et al., 2022; Alonso-Tetilla et al., 2024). Focusing on mid-IR selected quasars, Chen et al. (2015) and Andonie et al. (2022, 2024) showed that the fraction of obscured quasars increases with SFR such that the majority of quasars are obscured at the highest SFRs (>$300 M_{\odot} yr^{-1}$). Using high-resolution ALMA observations, Andonie et al. (2024) found that the star formation at the highest SFRs is so compact (≈ 0.5–3 kpc radius) that dust within the host galaxy can obscure the quasar even up to Compton-thick column densities ($N_H > 1.5 \times 10^{24} \text{ cm}^{-2}$). There is also evidence that the dust-covering factor in the nuclear regions is larger in the most heavily obscured systems, potentially due to chaotic gas inflows and/or radiation-pressure driven dusty outflows (e.g., Ramos

<footnote n="33">³³This is analogous to the apples-to-oranges idiom: don't compare apples with oranges if you want to understand different varieties of apples!</footnote>

<footnote n="34">³⁴We note that, even in the evolutionary model, all quasars should host a dusty/molecular “torus”. Therefore, although the torus covering factor is expected to vary between the evolving populations, unless it is total then at least some unobscured quasars will reside in the earlier more-obscured phase (but aligned so we see the accretion-disk emission) and some obscured quasars will reside in the later more-unobscured phase. However, if different quasar populations predominantly reside in different evolutionary phases then we will see key differences in the average properties and distributions between the populations.</footnote>

<footnote n="35">³⁵While the major merger of gas-rich galaxies is not a requirement to driving significant gas inflow, it does provide an indicator (a “time stamp”) on what point a given galaxy is within the merger sequence (i.e., close separation galaxies are likely later in the sequence than far separation galaxies).</footnote>
<figure><img src="image_29.png" /><figcaption>Figure 28: Quasar evolutionary model schematic showing the gas-inflow trigger options (merger or secular), the key quasar phases with differing amounts of obscuration and star formation, and the end result of either a quiescent galaxy or an inactive galaxy where the AGN and star formation have temporarily ceased. The inset panel provides a zoom in on the nuclear regions of the key quasar phases, showing the dust obscuration and indicating gas outflow via winds and jets. Source: the inset panel is adapted from Fig. 19 of Lyu and Rieke (2018).</figcaption></figure>

<table><caption>Table 2: Selected extreme source populations identified over the last decade, along with brief details including a key reference that defines the population.</caption><thead><tr><td>Extreme source population</td><td>Brief selection description</td><td>Definition reference</td></tr></thead><tbody><tr><td>Extremely Red Quasar (ERQ)</td><td>Extreme mid-IR-optical slope selection: typically selects dust-reddened quasars</td><td>Ross et al. (2015); Hamann et al. (2017)</td></tr><tr><td>Hot Dust Obscured Galaxy (Hot DOG)</td><td>Extreme mid-IR-optical slope and extreme near-IR colour selection: typically selects heavily obscured quasars</td><td>Eisenhardt et al. (2012)</td></tr><tr><td>Little Red Dot (LRD)</td><td>Red optical-near-IR colours and compact near-IR morphology: mixture of object types including broad-line systems</td><td>Matthee et al. (2024)</td></tr></tbody></table>

Almeida et al., 2011; Yamada et al., 2021, 2024), and some obscured quasars may be completely buried so that even their AGN-powered narrow-emission lines are weak or undetected (e.g., Greenwell et al., 2022, 2024a). At the extreme end of the obscured quasar population are Hot Dust Obscured Galaxies (hot DOGs; see Table 2), which are amongst the most luminous AGN known and are heavily obscured with large dust covering fractors (e.g., Wu et al., 2012; Stern et al., 2014; Assef et al., 2015; Fan et al., 2016b). They appear to be growing close to their maximum rates but are not exclusively merger driven (Fan et al., 2016a; Wu et al., 2018; Díaz-Santos et al., 2021).

Red quasars, systems where the accretion-disk emission is seen in the optical waveband but is at least partially obscured by dust (up-to $A_V \approx 5$ mags), are often the lynch-pin in quasar-evolutionary scenarios. They are expected to represent a key transition population caught between the fully obscured and unobscured quasar phases where the quasar is still actively removing the gas and dust along the line of sight through quasar-driven outflows (i.e., winds and jets). Red quasars were first identified by selecting radio-detected quasars several decades ago (Webster et al., 1995) and appear to represent a substantial

fraction (> 10–20%) of the overall quasar population (e.g., Banerji et al., 2012; Glikman et al., 2012). The colours of red quasars makes them difficult to distinguish from other “contaminant” populations in the optical waveband, limiting the reliable selection of red quasars from most large optical surveys such as the SDSS to systems with modest amounts of dust obscuration ($A_V < 1$ mag; Richards et al., 2003; Fawcett et al., 2022), with the exception of Extremely Red Quasars (ERQs), selected using optical-mid-IR colours but identified in modest numbers in the SDSS survey (e.g., Ross et al., 2015; Hamann et al., 2017); see Table 2. Consequently, many red-quasar samples are constructed using specific selection criteria sensitive to larger levels of dust obscuration, often requiring extensive follow-up observational campaigns even to characterise relatively small samples and often lacking consistently selected control samples.

The majority of red quasar studies over the last decade have shown that they have the properties expected for a transition population (prevalence for galaxy mergers, accreting at high Eddington ratios, and driving powerful outflows; e.g., Urrutia et al., 2012; Glikman et al., 2012, 2015; Banerji et al., 2015; DiPompeo et al., 2018; Calistro Rivera et al., 2021; Stacey
et al., 2022; Kim et al., 2024a). ERQs appear to represent the extreme end of the red quasar distribution with the largest luminosities, highest-velocity outflows, and absorbing column densities (e.g., Zakamska et al., 2016; Perrotta et al., 2019; Ma et al., 2024; Vayner et al., 2024). However, depending on the selection of the “control” sample, normal blue quasars can also be identified with similar properties to those of red quasars (e.g., Zakamska et al., 2019; Villar Martín et al., 2020).

Recent research exploiting the excellent source statistics of the SDSS quasar survey and defining careful control samples (but confirmed by other quasar-selection approaches) made the surprising discovery that red quasars emit excess but comparatively weak radio emission over that found for blue quasars (e.g., Klindt et al., 2019; Fawcett et al., 2020, 2022; Calistro Rivera et al., 2021; Rosario et al., 2021; Glikman et al., 2022; Yue et al., 2024a). The radio emission is predominantly produced on scales smaller than the host galaxy (<2–10 kpc) with no significant differences between red and blue quasars with either extended radio emission or radio-loud emission. As mentioned in §5.2, this excess radio emission is connected to the dust extinction towards the quasar, suggesting it may be caused by shocks due to interactions between low power jets/quasar-driven winds with the gas/dust in the circumnuclear/ISM environment (e.g., Hwang et al., 2018; Calistro Rivera et al., 2021, 2023; Fawcett et al., 2023; Petley et al., 2024). This therefore appears to be AGN feedback “in action” where gas/dust is being removed along the line of sight, as expected for a “transition” population. There are surprisingly few other differences in the overall quasar SEDs other than tentative excess mid-IR emission (potentially “dust-entrained” outflows; Calistro Rivera et al., 2021) suggesting this “blow-out” phase is brief.

The large-scale environment can also help determine the relationship between different quasar populations (see §4.4.2 for an overview). Based on the two-point correlation function, Petter et al. (2022) found no significant clustering differences for SDSS-selected quasars at $0.8 < z < 2.2$ as a function of optical colour, with the reddest and bluest quasar populations found to have characteristic dark-matter halo masses of $M_{\text{halo}} \approx 3 \times 10^{12} h^{-1} M_{\odot}$; this result was confirmed by CMB lensing. Interpreted within an evolutionary sequence, these results indeed suggest any timescale differences between red and blue quasar phases must be relatively short. By contrast, as already highlighted in §4.4.2, most two-point clustering analyses of obscured and unobscured quasars find that obscured quasars reside in more massive dark-matter halos ($M_{\text{halo}} \approx 10^{13} h^{-1} M_{\odot}$; Donoso et al., 2014; DiPompeo et al., 2016; Powell et al., 2018; Petter et al., 2023), a result again backed up by CMB lensing (e.g., DiPompeo et al., 2016; Petter et al., 2023). These results cannot be explained by the standard orientation model and argue for a more complex evolutionary connection between obscured and unobscured quasars where the obscuration often occurs on larger scales than the dusty “torus” (see §4.4.2), in broad agreement with recent results (Andonie et al., 2022, 2024). However, large-scale optical spectroscopic confirmation is required to improve on the current photometric redshifts for the majority of the obscured quasars and to rule out a significant population of non-quasars “contaminating” the clustering

signal (see Footnote 25).

### 5.4. The high-redshift frontier of quasars

High-redshift ($z > 4$) quasars trace the early growth phase of SMBHs and provide key constraints on the formation of BH “seeds”. However, their identification is hindered by two main factors: their rarity, with a number density of $Gpc^{-3}$, and at the very highest redshifts, the significant absorption of their radiation in the observed optical waveband, by material in the foreground intergalactic medium (Fan et al., 2023, for a review). Identifying these objects has thus typically required sensitive surveys covering large volumes with wide-field imaging, such as the SDSS (Fan et al., 2001, 2003), the Canada-France-Hawaii Telescope Legacy Survey (CFHTLS) (Willott et al., 2007, 2010), and the Panoramic Survey Telescope and Rapid Response System 1 (Pan-STARRS 1; Bañados et al., 2016).

Fig. 29 charts the history of record-breaking quasar redshifts from the first identified quasar (Schmidt, 1963) to the present day (cf. Stark et al., 2025, for record-breaking galaxy redshifts). The $z \approx 5$ barrier was broken in 1999 (Fan et al., 1999) and the first $z > 7$ quasar was identified in 2011 (Mortlock et al., 2011), requiring near-IR photometry; the highest-redshift quasar lies at $z = 7.642$ (Wang et al., 2021b). JWST has pushed the identification of potential AGN to $z \approx 10$ with sensitive near-IR spectroscopy (Goulding et al., 2023; Bogdán et al., 2024; Maiolino et al., 2024b; Napolitano et al., 2024). See Fig. 26 for spectra of the three $z \approx 10$ systems, all of which have relatively narrow emission lines (FWHM<1000 km s⁻¹). The highest-redshift JWST AGN with clear broad Balmer emission lines (FWHM>1000 km s⁻¹) lie at $z = 8.500$ (Kokorev et al., 2023), $z = 8.632$ (Tripodi et al., 2024), $z = 8.679$ (Larson et al., 2023), and $z = 9.288$ (Taylor et al., 2025).

#### 5.4.1. Increasing the demographics of high-redshift AGN

Fig. 25 provides a comprehensive census of the bolometric luminosities and SMBH masses of high-redshift quasars. To date > 4000 and > 200 quasars have been identified at $z > 4$ and $z > 6$ (Perger et al., 2017; Fan et al., 2023), respectively.<fnref n="36" /> Substantially fewer high-redshift quasars have SMBH masses due to the need for near-IR spectroscopy to obtain reliable virial SMBH mass tracer emission lines (e.g. Mg II; see §5.1). The majority of the $z > 4$ quasars are extremely luminous ($L_{\text{bol}} \ge 10^{46}$ erg/s) and host massive black holes ($M_{\text{BH}} \approx 10^9 M_{\odot}$). We refer to these systems as **bright quasars**: they represent the tip of the active SMBH population, although in terms of their multi-wavelength and spectral properties they are similar to those of luminous quasars at lower redshifts; see §5.2. A major push over the last decade has been to identify and characterise fainter high-redshift quasars to gain a less biased view on the early assembly of more typical SMBHs, their environments, and co-evolution with their hosts.

At the time of writing, about a hundred faint high-redshift (non-JWST) quasars with bolometric luminosities ranging from

<footnote n="36">³⁶We note that the apparent paucity of quasars at $z \approx 5.3–5.7$ is due to the similarity of their observed optical colours with those of late-type stars, making their selection and identification particularly challenging (Yang et al., 2017b).</footnote>
<figure><img src="image_30.png" /></figure>

Figure 29: Record-breaking quasar redshifts as a function of their year of discovery, starting from the first “distant” quasar at z = 0.158 (Schmidt, 1963) to the highest-redshift non-*JWST*-identified quasar at z = 7.642 (Wang et al., 2021b) and the highest-redshift candidate AGN identified by *JWST* at z = 10.603 (Maiolino et al., 2024b); the record-breaking redshifts for the *JWST* systems are indicated with an open symbol (Kokorev et al., 2023; Goulding et al., 2023; Maiolino et al., 2024b; Bogdán et al., 2024). Key developments that have allowed for the identification of higher-redshift quasars are indicated on the top axis and the redshift where Ly α enters the near-IR waveband is plotted as a horizontal magenta line. Quasars at z ≥ 5 began to be uncovered about 20 years ago and are now being characterized by the latest telescopes: *JWST*, *Euclid*, and soon *Roman*. Looking ahead, the gravitational wave observatory LISA will enable the discovery of merging BHs at any redshift, complementing other new telescopes like *Athena* and *AXIS* that will detect the AGN X-ray emission.

$L_{\text{bol}} = 10^{45}$ to $10^{46}$ erg/s have been identified at $z = 6-7$ in dedicated deep quasar surveys, such as CFHQS (Canada-France High-redshift Quasar Survey; Willott et al., 2007, 2010) and SHELLQs (Subaru High-redshift Exploration of Low-Luminosity Quasars; Matsuoka et al., 2016, 2018a). We refer to these systems as **faint quasars**: the identification of these faint quasars has enabled the characterization of the quasar luminosity function at $z = 5$ (e.g., McGreer et al., 2018; Niida et al., 2020) and $z = 6$ (Matsuoka et al., 2018b), down to rest-frame UV $M_{1450} \sim -22.3$. These faint high-redshift quasars are powered by SMBHs with a wide range of masses ($M_{BH} \sim 10^{7.5} - 10^9 M_{\odot}$) and Eddington ratios ($f_{Edd} = 0.1-1$; Onoue et al., 2019; Takahashi et al., 2024); in broad terms they appear to be a mix of the lower Eddington-ratio tail of the more luminous $z \approx 6$ quasar population and lower mass, higher-Eddington ratio quasars. Comparisons to the extrapolation from observed X-ray luminosity functions, primarily constrained up to $z \sim 4$ (see e.g., Pouliasis et al., 2024, for new constraints at $z > 4$), indicate that optical-near-IR photometric selection approaches may overlook a large fraction of the quasar population at $z \ge 6$ due to obscuration (Vito et al., 2018).

*JWST* has now pushed the identification of high-redshift AGN populations down to bolometric luminosities of ~ $10^{43}$ – $10^{44}$ erg/s. We refer to these systems as **JWST AGN**: their spectral and multi-wavelength properties differ from those of comparable-redshift quasars, making them less-reliable candidate AGN; however, the presence of an AGN appears the most probable solution for the majority. Within less than 2 years, already > 100 UV-faint *JWST* AGN at $z \ge 4 - 8.5$

have been identified with weak broad (≈ 1000–4000 km s⁻¹) Balmer emission lines (e.g., Greene et al., 2023; Kocevski et al., 2023; Maiolino et al., 2023; Übler et al., 2023; Matthee et al., 2024; Taylor et al., 2024) and a similar number of narrow-line systems (Scholtz et al., 2023; Mazzolari et al., 2024b; Treiber et al., 2024), due to these systems being common even in small (<< deg²) *JWST* fields; see Fig. 25. Below we summarise the basic properties of the key detected *JWST* AGN populations.

**Broad-line JWST AGN:** ≈ 5–10% of the systems spectroscopically identified in high-redshift galaxy searches have broad Balmer emission lines (Harikane et al., 2023; Maiolino et al., 2023; Treiber et al., 2024). Spectroscopically these systems are not like normal quasars: they are metal poor and lack strong high-ionisation emission lines (such as C IVλ1549), properties more similar to those of AGN identified in local metal-weak dwarf galaxies (Trefoloni et al., 2024); see also §4.6. Assuming the broad emission lines are due to virialised gas around an accreting SMBH (see §5.1), they have masses of $M_{BH} \approx 10^6–10^8 M_{\odot}$ and Eddington ratios of ≈ 0.1–1; we note the absence of broadened forbidden lines makes a sole outflow origin for the broad lines unlikely. A narrow, typically blueshifted, absorption feature is seen in the Balmer lines for ≈ 10–20% of the systems, interpreted as a high column density of gas and potentially related to an outflow and/or the outer region of the BLR (Inayoshi and Maiolino, 2024; Juodžbalis et al., 2024; Kocevski et al., 2024; Lin et al., 2024; Matthee et al., 2024; Ji et al., 2025); note, < 0.1% of low-redshift broad-line AGN show this feature (Lin et al., 2024).

Very few systems are X-ray detected, either individually or
from deep stacking analyses, indicating they are systematically
X-ray weak by at least 1 order of magnitude (Kocevski et al.,
2024; Lambrides et al., 2024; Maiolino et al., 2024a); none are
radio detected but this is consistent with them being radio quiet
(Mazzolari et al., 2024a). The origin of the X-ray weakness
is unclear: it may be due to heavy Compton-thick absorption
(although given the high rest-frame energies probed at $z > 4$,
$> 2.5–40$ keV, the absorption would need to be extreme) or the
systems may be intrinsically X-ray weak. The most-promising
(and explored) option for the intrinsically X-ray weak scenario
is super-Eddington accretion which is expected to be faint due
to an intrinsically “soft” X-ray SED (like NLS1s) and absorp-
tion from the slim accretion disc (Madau and Haardt, 2024; In-
ayoshi et al., 2024; Pacucci and Narayan, 2024; King, 2025).
The lack of high-ionisation emission lines is consistent with a
super-Eddington scenario (Lambrides et al., 2024) and would
make these JWST systems similar to X-ray weak WLQs; see
§5.2. An implication of the super-Eddington accretion scenario
is the luminosity and/or the SMBH masses are overestimated by
a factor of at least a few (e.g., Maiolino et al., 2024a; Madau and
Haardt, 2024). About 10–30% of the broad-line systems are red
and have compact near-IR morphologies, with the fraction in-
creasing strongly with redshift (Hainline et al., 2024; Inayoshi,
2025; Ma et al., 2025): they are referred to as “little red dots”.

**Little Red Dots (LRDs):** significant attention has been focused on understanding the nature of LRDs; see Table 2. LRDs are selected to be compact in NIRCam imaging (<100–200 pc) and typically have steep red rest-frame optical colours with blue rest-frame UV colours, giving them an unusual “V” shaped UV-optical SED (Labbe et al., 2023).<fnref n="37" /> Depending on the exact selection, ≈ 20–80% host broad Balmer emission lines and have rest-frame UV-optical spectral, X-ray, and radio properties similar to the JWST broad-line systems (Greene et al., 2023; Akins et al., 2024; Ananna et al., 2024; Hainline et al., 2024; Kocevski et al., 2024; Matthee et al., 2024; Perger et al., 2024; Yue et al., 2024b). If interpreted as AGN then they appear to be dust obscured (up-to $A_V$ ≈ 5 mag) and the red optical colour is predominantly due to the reddened AGN while the blue UV colour is produced by scattered AGN emission or a young stellar population (Barro et al., 2023; Labbe et al., 2023; Greene et al., 2023; Kocevski et al., 2023; Killi et al., 2024); ≈ 30% of LRDs have extended rest-frame UV emission, as expected from host-galaxy emission, which is often asymmetric potentially due to galaxy interactions (Kocevski et al., 2024; Rinaldi et al., 2024). However, LRDs are weak at rest-frame near-IR wavelengths, suggesting a lack of obscuring dust on nuclear scales and therefore inconsistent with a standard obscured AGN scenario (Akins et al., 2024; Li et al., 2025; Pérez-González et al., 2024; Williams et al., 2024; Casey et al., 2024; Li et al., 2025). Few LRDs are variable in the rest-frame UV-optical casting some doubt on an AGN hypothesis (Kokubo and Harikane, 2024; Tee et al., 2024; Zhang et al., 2024). How-

ever, spectral monitoring of a $z \approx 7$ multiply lensed LRD over a uniquely long timescale (observed-frame 22 years; rest-frame 2.7 years) identified broad Balmer emission line variations consistent with those expected for a damped random walk model of an accretion disk around a $\approx 4 \times 10^7$ M$_\odot$ SMBH (Furtak et al., 2025).

Many LRDs have a spectral break around the wavelength of the Balmer break (Labbe et al., 2024; Ma et al., 2024; Wang et al., 2024). Analysing a large sample of LRD JWST spectra, Setton et al. (2024) found the inflexion in the “V” shaped SED corresponds to the Balmer break in ≈ 50% of sources; see Fig. 26 for the break in the composite LRD spectrum. They argued this break is unlikely to occur for dust-reddened AGN with a range of dust extinctions and host-galaxy strengths and, adopting “Occam’s razor”, proposed the properties of LRDs are from a single origin, potentially due to late-type stars. However, a growing number of studies have argued that an AGN buried in an extremely dense cocoon of gas (what is sometimes referred to as a “black-hole star” scenario) would also produce such feature and could simultaneously explain the Balmer absorption lines seen in the spectra of some JWST AGN (e.g., Inayoshi and Maiolino, 2024; de Graaff et al., 2025; Ji et al., 2025; Naidu et al., 2025). A potential consequence of the large gas density is scattering of the broad emission lines, producing an extended velocity component and artificially increasing SMBH masses by ≈ 1–2 orders of magnitude, when calculated using the virial SMBH estimator (Rusakov et al., 2025), although this scenario has been disputed for at least some JWST AGN (Juodžbalis et al., 2025). Overall, despite their unusual properties, it therefore appears likely that the majority of LRDs host AGN, potentially representing an early SMBH growth phase which would help explain their strong redshift dependence (Inayoshi, 2025) and overmassive SMBHs compared to their hosts; however, as for any broad-band photometric selection, it is quite likely there is some heterogeneity in the LRD population (see, for example Setton et al., 2024; Zwick et al., 2025, for some alternative non-SMBH and non-AGN scenarios).

**Narrow-line JWST AGN:** an abundant narrow-line AGN population is expected at high redshift with similar properties to the JWST-identified broad-line systems but lacking the broad Balmer emission lines. Their metal-poor nature means that typical AGN emission-line diagnostics (Baldwin, 1977; Veilleux and Osterbrock, 1987) are not reliable: indeed, the broad-line systems reside in the same region as low-mass star-forming galaxies. On the basis of new emission-line diagnostics developed to distinguish between high and low ionising continua in high-redshift low-mass galaxies, ≈ 20% of JWST-identified galaxies at $z = 2–9$ are found to have spectral properties suggesting the presence of AGN activity (Scholtz et al., 2023; Mazzolari et al., 2024c,b; Treiber et al., 2024). Their estimated bolometric luminosities (≈ $10^{42}–10^{45}$ erg s⁻¹) overlap with those of the broad-line systems but extend down to very low-luminosity systems. Few are X-ray detected, suggesting they are as X-ray weak as the broad-line systems (Maiolino et al., 2024a; Mazzolari et al., 2024b). The AGN signatures are clear in the NIRSpec spectra of some sources but they are weak over-

<footnote n="37">³⁷LRDs were first identified in high-redshift galaxy searches and, prior to spectroscopic confirmation, were interpreted as a massive (≈ $10^{10}–10^{11}$ M$_\odot$) and surprisingly ubiquitous “quiescent” galaxy population at $z \approx 6$ (Labbé et al., 2023), the so-called “universe breaker” galaxies</footnote>
all (e.g., see Fig. 1 of Mazzolari et al. (2024b)) and further confirmation of AGN in the weaker sources would be valuable.

**Very high-redshift JWST AGN:** to date, three $z > 10$ systems have been proposed as very-high redshift AGN: UHZ-1 at $z = 10.073$ (Goulding et al., 2023; Brogan et al., 2023), GHZ9 at $z = 10.145$ (Kovács et al., 2024; Napolitano et al., 2024), and GNz-11 at $z = 10.603$ (Maiolino et al., 2024b); see Fig. 26. The optical spectrum of UHZ-1 is consistent with co-eval star-forming galaxies with no evidence for an AGN but it is weakly X-ray detected at 2–7 keV, suggesting an optically obscured AGN. GHZ9 has high-ionisation lines, including C IV, consistent with an AGN and it is X-ray detected at 0.5–3 keV: the derived X-ray luminosity is also consistent with the 440 nm continuum luminosity, giving strong evidence for an AGN. GNz-11 is not X-ray detected but has several high-excitation and high-critical density lines (notably Nitrogen), including modestly broadened emission lines (FWHM ≈ 450 km s⁻¹), giving support for the presence of an AGN. Recent JWST MIRI spectroscopy of GNz-11 detected only narrow Hα emission without a clear broad component, casting some doubt on the broad-line AGN classification of GNz-11 (Álvarez-Márquez et al., 2024).

The rapid progress in identifying and studying the JWST AGN has come down to the ultra-sensitive near-IR–mid-IR spectroscopic and high-resolution imaging capabilities of the JWST facility. More work is required to understand why they differ from normal co-eval quasars and to rule out non-AGN solutions. However, at the time of writing, the presence of an AGN seems plausible for the majority. Although such systems were theoretically expected at high redshift from BH seed formation and growth models, their number density came as a surprise. Specifically, the number density of the JWST AGN is about $10^{-5}$ cMpc⁻³ in the UV magnitude range $M_{UV} = -21, -17$ at $z \sim 5$; the exact value depends on the magnitude range, selection, and volume covered by the observations (Matthee et al., 2024; Greene et al., 2023). The first constraints on the bolometric luminosity function indicate that the JWST broad-line AGN are as numerous, or even more numerous, than the AGN predicted by large-scale cosmological simulations when considering all simulated AGN are visible and none obscured (Matthee et al., 2024; Habouzit, 2024).

We illustrate this in Fig. 30 with the compilation of bolometric luminosity functions derived from different observational samples (explicitly selecting LRDs or broader AGN candidate samples) and the range probe by the AGN population of eight large-scale cosmological simulations (HORIZON-AGN, ILLUSTRIS, EAGLE, TNGs, SIMBA, ASTRID; Dubois et al., 2014a; Genel et al., 2014; Schaye et al., 2015; Pillepich et al., 2018; Nelson et al., 2019; Davé et al., 2019; Ni et al., 2022b). Including both broad-line AGN from various JWST surveys and narrow-line AGN identified in Scholtz et al. (2023) and the best fit of Grazian et al. (2024) yields the same conclusions when compared to semi-analytical models. The comparison between theoretical models (i.e. cosmological simulations and semi-analytical models) and observational constraints indicates that we either need to refine our theoretical models of BH physics to facilitate SMBH growth or to understand why some AGN can-

<figure><img src="image_31.png" /><figcaption>Figure 30: AGN bolometric luminosity function at $z \sim 5$. The shaded green region shows the combination of the luminosity functions, including Poisson uncertainties, from the HORIZON-AGN, ILLUSTRIS (TNG 50, TNG 100, TNG 300), EAGLE, SIMBA, and ASTRID large-scale cosmological simulations. Only AGN in galaxies with a total stellar mass of $M_{\star} \ge 10^{8.5}$ M$_{\odot}$ are considered. Symbols represent observational constraints from JWST for a sample of AGN candidates at $z \sim 5$ with a redshift range of $z \approx 4–6.5$ (Matthee et al., 2024; Greene et al., 2023; Maiolino et al., 2023; Harikane et al., 2023; Akins et al., 2024; Guo et al., 2024a). The grey curve is the pre-JWST luminosity function at $z \sim 5$ (Shen et al., 2020). The predicted luminosity function is lower than the observational constraints in some simulations, which in some cases only includes broad-line AGN (or a subset), highlighting the tension between theoretical models and JWST observations. Source: adapted from Fig. 4 of Habouzit (2024).</figcaption></figure>

didates may not actually be AGN or why their derived properties (SMBH mass and AGN luminosity) may be overestimated, for example. Theoretically, this includes processes such as BH formation, gas accretion, super-Eddington regime, SMBH dynamics and coalescence. Additionally, sustaining such efficient SMBH growth is likely closely linked to galaxy physics, such as feedback processes, morphology of the host, the gas reservoir and the rapid replenishment of cold gas.

### 5.4.2. Galactic and large-scale environments of high-redshift quasars and JWST AGN

Over the last decade since AH12, following the discoveries of high-redshift quasars (Fig. 29), the community has begun investigating their environments to understand their fuelling mechanisms and how they could have assembled so quickly. Significant progress has already been made in several key aspects, thanks to advances from new instruments such as ALMA and VLT-MUSE, as well as developments in theoretical studies. This includes the creation of dedicated semi-analytical models following the merger trees of the progenitors of quasars (Valiante et al., 2014) and cosmological simulations of quasar environments (Costa et al., 2014b), as described in §5.5 and in the Valiante et al. (2017) review, for example.

Quasars are thought to inhabit highly biased regions of the Universe, where SMBH growth is fuelled by an increased number of mergers and the convergence of filaments supplying
fresh gas. In the ΛCDM cosmology, this hypothesis would be supported by an observed increase in galaxy counts within the quasar's field of view. However, the observations of Lyman break galaxies (LBGs, probing a redshift uncertainty of Δz ~ 1) and more recently searches for lower-mass Lyman-alpha emitters (LAEs, Δz ~ 0.1) revealed a diversity of environments from overdense regions (e.g., Stiavelli et al., 2005; Willott et al., 2005; Zheng et al., 2006; Kim et al., 2009; Husband et al., 2013; Morselli et al., 2014; McGreer et al., 2014; Meyer et al., 2022) to environments consistent with blank fields (e.g., Bañados et al., 2013; Simpson et al., 2014; Mazzucchelli et al., 2017). Atypical fields, such as the Utsumi et al. (2010) discovery of an overdense ring-like structure at 3 Mpc from the quasar but no overdensity closer to the quasar, or the low density of LAEs within 3 Mpc (projected distance) found in the large-scale and LBG-overdense structure studied in Ota et al. (2018), could suggest negative feedback on the formation and regulation of galaxies at small physical scales (with potentially a larger effect on low-mass galaxies). Deep-and-wide observations are needed, such as from Subaru/Suprime-Cam or LBT/LBC, as used in Mignoli et al. (2020, with additional observation from MUSE for LAEs) for the spectroscopic confirmation of the overdensity presented with imaging in Balmaverde et al. (2017). [C II]λ158 µm observations on smaller scales with ALMA revealed that a significant fraction of quasars with $M_{BH} \ge 3 \times 10^8 M_{\odot}$ (≈ 16% of a sample of 25 quasars at $z > 5.9$) have a close companion within 100 kpc (see also Farina et al., 2017; Trakhtenbrot et al., 2017). ALMA allows for both the discovery of dusty companion galaxies that cannot be selected at any other wavelength, but also to determine accurate [C II] redshifts.

More recently, large observational campaigns have been designed to provide statistical constraints on the quasar environment by observing, for example, 25 quasars fields with JWST-NIRCam WFSS (Wang et al., 2023). While the long-standing absence of a consensus on the environments of high-redshift quasars was in part attributed to reduced observational sensitivity missing galaxies and small fields of view, the first JWST results demonstrate that the diversity of environments persist even with new deeper observations, including an overdense filamentary structure of 10 [O III]λ5007Å emitters at $z = 6.61$ (Wang et al., 2023), dense regions with close (potentially) merging galaxies with JWST NIRSpec IFU (Marshall et al., 2023, 2024), and fields consistent with average galaxy density (Eilers et al., 2024). Whether high-redshift quasars reside in a broad range of dark matter haloes and whether their strong feedback suppresses galaxy formation in their surroundings is still to be determined (Costa et al., 2014b; Costa, 2024; Habouzit et al., 2019, for a theoretical perspective from cosmological simulations).

At low redshift ($z \le 1$), quasar host galaxies exhibit diverse morphological properties (Ishino et al., 2020; Bettoni et al., 2015; Falomo et al., 2014; Bahcall et al., 1997). However, the brightness of quasars in the optical/near-IR has made it difficult to determine the host galaxies of luminous high-redshift quasars. Taking advantage of the relatively faint

millimetre-wavelength emission from the quasars, ALMA has placed constraints on the dynamical mass of the quasars hosts as well as their morphology. Recent [C II] observations resolving the hosts of $z \ge 6$ quasars with $\le 1.2$ kpc (Neeleman et al., 2021), revealed that about a third of high-redshift quasars have disturbed [C II] emission profiles characteristic of galaxy mergers, a third reside in extended disk galaxies (i.e., with smooth velocity gradients consistent with emission arising from rotating disks), and the last third show no clear velocity gradients and are consistent with dispersion-dominated galaxies with turbulent gas (i.e., compact bulge-like galaxies). The diversity of high-redshift quasar morphologies, similar to those at low redshifts suggests a diverse range of mechanisms for SMBH fuelling, not just galaxy mergers. Very high spatial resolution ALMA imaging down to 200–400 pc, similar to that achieved for nearby galaxies, has been achieved for a few high-redshift quasar hosts and, indeed, reveals complex velocity structure and nearby companions hinting towards gas funnelled to the galaxy centre by interactions, non-regular structures in which gas flow could be driven by bar motions, and surprisingly compact emission with regular disk emission and a lot of gas within the central 100 pc but with no clear signature of the presence of a massive SMBH on the gas dynamics (Walter et al., 2022; Venemans et al., 2020; Neeleman et al., 2019). In the near future, more statistical studies are needed to understand the relative role of the galactic morphologies and kinematics, whether they are intrinsically distinct or instead part of an evolutionary path of the quasars, and whether they are linked to distinct channels of SMBH fuelling. The new exquisite resolution achieved in observations allow comparison with recent cosmological hydrodynamical zoom-in simulations of quasars hosts following their gas and stellar kinematics with a resolution of 50 pc (Lupi et al., 2022). The simulations show that the host galaxies experience several mergers and a rapid evolution of the gas and stellar distributions, which can go from a perturbed disk, to a compact system, and later settle down into a more extended system in less than 200 Myr.

Over the years through the growth of observational samples of quasars with SMBH mass estimates, it has become possible to investigate the masses of the SMBHs powering the high-redshift quasars as a function of the dynamical masses of their hosts. Results have shown that bright quasars tend to be located above the $M_{BH} - M_{\star}$ mass scaling relation derived in the local Universe (e.g., Kormendy and Ho, 2013, assuming $M_{\star} = M_{dyn}$), while faint quasars are consistent with, and sometimes lie below, the local relation (Izumi et al., 2019). The stellar component of the hosts, however, remained a missing piece until the advent of JWST which prevented a direct $M_{BH} - M_{\star}$ comparison to local systems. Decomposition of the point-source quasar from its host galaxy requires a precise PSF, a method that was successfully employed with HST but limited to quasars at $z \le 2$. The exciting first detections of the stellar light from two quasar hosts with JWST were presented in Ding et al. (2023). This work was quickly followed by Stone et al. (2023); Yue et al. (2024c); Stone et al. (2024); Marshall et al. (2024), mostly for bright quasars. Although still with a
limited number of systems and large error bars on both $M_{BH}$ and $M_\star$, quasars reside mostly in galaxies within the range $M_\star \sim 10^{10} - 10^{1.5} M_\odot$, as shown in Fig. 15. Initial results are consistent with the ALMA constraints and suggest that faint quasars tend to align with the $M_{BH} - M_\star$ relation in the local Universe, while bright quasars have higher $M_{BH}/M_\star$ ratios, although not yet significantly different from the most massive SMBHs found in local large elliptical galaxies; see Fig. 15. These potential overmassive quasars reignite the debate on whether super-Eddington accretion can explain these systems. Although new high-resolution simulations of quasars suggest that prolonged super-Eddington phases can occur in gas-rich, high-redshift environments (Lupi et al., 2024a), even in the challenging scenario of relatively high spin and associated strong jet feedback, further investigation is required. This is because we are still unable to resolve the physical processes at very small scales in these simulations (e.g., accretion at the Bondi radius, early expansion of feedback-shocked gas).

In terms of their galaxy properties, JWST AGN are found in a wide range of hosts with stellar mass estimates from ~ $10^8$ to $10^{11}$ $M_\odot$ (the majority < $10^{10}$ $M_\odot$), for both broad-line and narrow-line systems (e.g. Maiolino et al., 2023; Harikane et al., 2023; Scholtz et al., 2023; Kocevski et al., 2024). However, it is important to caveat that contamination from the AGN leads to uncertainties in the stellar-mass estimate: the impact is less severe for the narrow-line candidates, which is dominated by the galaxy but the AGN can “contaminate” the Balmer, impacting SFR estimates. Indeed, uncertainties on both SMBH and stellar mass estimates can be as much as an order of magnitude (Li et al., 2024a; King, 2024; Lupi et al., 2024b, for discussions on using standard SMBH mass methods calibrated against local relations, selection biases, and $M_{BH} - M_\star$ measurement uncertainties). For example, SMBH mass measurements assume virialised gas, and could be over-estimated if the SMBH is in a super-Eddington regime, undergoing a TDE phase (see §3.2), or a significant fraction of the broadened emission is scattered (Rusakov et al., 2025). With these caveats in mind, a number of JWST AGN are found to have $M_{BH}/M_\star$ ratios similar to galaxies in the local Universe, although a non-negligible fraction have larger ratios making their SMBHs overmassive compared to their hosts, similar to that found for bright quasars but at lower masses; see Fig. 15. However, the large uncertainties on the derived quantities preclude strong conclusions on whether the JWST AGN host overmassive SMBHs and, consequently, their origins (e.g., massive BH seeds, BH seeds forming prior to the host galaxies, efficient growth, atypical host-galaxy evolution). The uncertain origin of the blue and red components of the SED in LRD systems further complicate the estimation of their $M_{BH}/M_\star$ mass ratios: if the continuum emission is primarily due to an AGN with absorption by dense neutral gas clumps on the shorter-wavelength side of the Balmer break then estimates of $M_\star$ based on a stellar origin for the continuum would be overestimated (Inayoshi and Maiolino, 2024). Finally, some of the outliers in the $M_{BH} - M_\star$ diagram (for relatively low stellar masses with $\le 10^9 M_\odot$) are found not to be outliers in $M_{BH} - \sigma$ or $M_{BH} - M_{dyn}$ diagrams with respect to what would

be expected from observations in the local Universe (Maiolino et al., 2023).

The majority of the AGN candidates detected by JWST were identified serendipitously in extensive spectroscopic surveys, in both average fields and quasar fields (with candidates in the foreground and background of the quasars). Despite the modest source statistics and cosmological volumes probed by the JWST surveys, it already seems likely that the JWST AGN reside in lower-mass haloes than the high-redshift quasars. Performing angular and projected cross-correlation analyses of 27 JWST AGN at $5 < z < 6$ with 679 photometrically selected galaxies, Arita et al. (2025) estimated a typical dark-matter halo mass for the JWST AGN of $M_{halo} \approx 3 \times 10^{11} h^{-1} M_\odot$, an order of magnitude lower than that of faint and bright high-redshift quasars (Arita et al., 2023); see also §4.4 & 5.3. Two-point cross correlation analyses of LRDs with co-eval H$\alpha$ emitting galaxies over $z \approx 4-6$ suggest that LRDs reside in similarly massive dark-matter halos of $M_{halo} \approx (1-2) \times 10^{11} h^{-1} M_\odot$ (Lin et al., 2025); complementary constraints are also placed by Pizzati et al. (2024) who, on the basis of the source density of LRDs at $z \approx 4-6$ and the relative number of coeval dark-matter haloes, argued that LRDs must reside in haloes of $M_{halo} < 5 \times 10^{11} h^{-1} M_\odot$ to avoid a halo occupancy of > 100%. These results indicate that the JWST AGN reside in dark-matter haloes an order of magnitude less massive than quasars (see §5.3), suggesting that they represent either an entirely different population or a much earlier SMBH growth phase; however, we note that at least one $z \approx 7$ LRD is found in a much denser large-scale environment (Schindler et al., 2024), suggesting some diversity in LRD environments.

On smaller spatial scales, there is evidence for a number of JWST dual AGN candidates at $z = 4-6$ (Maiolino et al., 2023) and $z = 3-3.5$ (Perna et al., 2023a,b); see §4.5 for dual AGN in regular quasars and AGN. This highlights the potential increased role of SMBH mergers in the growth of SMBHs at high redshift, which is theoretically expected to be much less significant than the dominant gas accretion channel, and the importance of galaxy interactions in triggering phases of enhanced gas accretion.

### 5.4.3. Origins of the seeds of massive black holes: models and observational constraints

Prior to JWST, observational facilities seeking evidence of the first BHs (i.e., the black hole seeds), could only probe the very tip of the BH distribution at $z \ge 5$ with the observations of extremely luminous and rare quasars with a number density of ~ $10^3 h^{-3}$, and powered by BHs of $M_{BH} \ge 10^8 M_\odot$. The existence of quasars at $z \ge 6$ more than a decade ago and the quasars discovered at $z = 7.5-7.6$ a few years ago suggest that their progenitors, referred to as BH seeds, must have formed with masses ranging from a few hundred solar masses (“light seeds”) up-to a million solar masses (“heavy seeds”) in the very early Universe. Some of the new JWST discoveries, such as the broad-line AGN at $z > 8.5$ with $M_{BH} \sim 10^8 M_\odot$ (Kokorev et al., 2023; Tripodi et al., 2024; Taylor et al., 2025), push even further the question of BH formation; see Fig. 31.
Contemporary theories on BH seed formation can be categorized into two main groups (Inayoshi et al., 2020; Volonteri et al., 2021, for recent reviews). First, seeds may have originated before even the formation of the first galaxies, arising from the collapse of high peak density fluctuations or cosmic string loops (e.g., Carr and Kühnel, 2020). These “cosmological” pathways predict seed formation as early as inflation and exhibit a broad range of initial masses ranging from 1 gram to several $10^5$ M$_{\odot}$, contingent on the precise time of formation. Alternatively, more widely studied pathways involve the formation of BH seeds alongside the first galaxies, possibly within either primeval or more evolved galaxies. The formation occurs through the gravitational collapse of a single stellar object or through the hierarchical merger of stars or stellar black holes. Each avenue for seed formation, including those described further below, carries its own set of advantages and drawbacks to explain the massive SMBHs sitting in the heart of almost all galaxies today. Importantly, these pathways are not theoretically mutually exclusive: the origins of the observed population of SMBHs may be diverse.

In the first theoretical endeavours, BHs were postulated to have formed from the remnants of the gravitational collapse of the first-generation stars, known as Pop III (Madau and Rees, 2001; Klessen and Glover, 2023). However, only the most massive Pop III stars would leave behind a seed with the capability to overcome the numerous challenges in evolving into a central SMBH. Furthermore, pair-instability supernovae occurring in massive Pop III stars with $140 - 260$ M$_{\odot}$ might not leave any remnants, leaving only stars with a mass of $\ge 260$, M$_{\odot}$ to form “light seeds” of a few $\sim 100$ M$_{\odot}$ at best. Understanding the potential contribution of this scenario to the overall BH population is further hindered by uncertainties regarding the initial mass function of Pop III stars, the shallow gravitational potential well of the host mini-halos, which could impede gas retention and, consequently, prevent efficient accretion onto the black holes. Additionally, challenges include the extended sinking time of low-mass BHs toward the galactic centre and the possible recoil of BHs from the central regions. The formation of more massive “heavy seeds” ($\sim 10^{4-5}$ M$_{\odot}$) would ease some of these important hurdles and make the existence of high-redshift quasars $< 1$ Gyr since the Big Bang ($z > 6$) less complex to elucidate (e.g., Bromm and Loeb, 2003; Begelman et al., 2006; Lodato and Natarajan, 2006). Under specific conditions, gas cooling and fragmentation can be delayed in primordial atomic cooling halos and instead collapse monolithically into one single protostar. With a sufficiently large accretion rate, a protostar can grow into a Super Massive Star (SMS) and collapse into a “heavy seed”. Conditions to avoid gas fragmentation have been increasingly studied over the last decade: intense photo-dissociating (Lyman-Werner) radiation from a neighboring galaxy (suppressing H$_2$ cooling, e.g., Visbal et al., 2014), dynamical heating from rapid halo mergers (increasing the heating rate, Wise et al., 2019), or large residual baryonic streaming motions from recombination (preventing gas infall and contraction into dark matter haloes), turbulence, and even mergers of evolved galaxies (Mayer and Bonoli, 2019; Mayer et al., 2024). Another important pathway is the

<figure><img src="image_32.png" /><figcaption>Figure 31: Black-hole mass versus redshift for several $z > 6$ quasars (grey circles) and $z > 8$ JWST AGN (grey triangles), focused on the $z = 8.632$ broad-line system from Tripodi et al. (2024) (yellow star) and the extreme $z = 6.327$ quasar from Mazzucchelli et al. (2023) (magenta circle) showing that such massive BHs could be related through Eddington-limited accretion from a direct collapse BH seed. Other potential sub Eddington and super Eddington scenarios are also shown from different seed masses as well as early SMBH growth simulations (Bennett et al., 2024) and semi-analytical models (Schneider et al., 2023; Cammelli et al., 2025). Source: adapted from Fig. 4 of Tripodi et al. (2024) and updated for recent $z > 8$ JWST AGN.</figcaption></figure>

formation of seeds with $\sim 10^{2-3}$ M$_{\odot}$ through hierarchical mergers of massive stars in compact and dense stellar clusters (e.g., Devecchi and Volonteri, 2009; Sakurai et al., 2017; Rantala and Naab, 2025). JWST probes the existence of such clusters as early as $z \sim 10$ (Adamo et al., 2024b; Vanzella et al., 2023; Mowla et al., 2024). Relatively massive BH seeds formed in stellar clusters would benefit from favourable dynamics and environment for rapid growth.

Although most of the theoretical pathways for BH seed formation have been around for quite some time$^{38}$, the last decade has been crucial to assess their detailed physical aspects (e.g., fragmentation of gas clouds or their collapse into protostars, accretion onto protostars and whether feedback processes inhibit growth) through cosmological simulations and semi-analytical models (e.g., Inayoshi and Tanaka, 2015; Chon et al., 2016; Regan et al., 2016; Wise et al., 2019; Latif et al., 2022).

What have we learned from recent observations with JWST? So far, JWST has provided three key insights. First, it has pushed the redshift record for the earliest observed SMBHs in the Universe to $z \approx 8.5-10.6$; see Figs. 25 & 29. The SMBH masses of these systems are highly uncertain but broadly appear to lie in the range of $\approx 10^7-10^8$ M$_{\odot}$, measured from the broad emission lines for the $z > 8.5$ broad-line systems (see §5.1) and estimated from the AGN luminosity and assuming an

<footnote>$^{38}$The paper of Martin Rees (1978) remains a solid and enlightening foundation of most contemporary theories (Rees, 1978).</footnote>
Eddington ratio for the $z \approx 10$ narrow-line systems. However, these SMBH masses are broadly consistent with that required to explain the previously known high-redshift quasars (Schneider et al., 2023). In other words, in terms of their growth history, these systems at $z \sim 8.5 - 10$ could be the progenitors of the high-redshift quasars (see discussion in Volonteri et al., 2023); see Fig. 31. The occupation fraction — the number of BHs in bins of galaxy properties such as stellar mass or velocity dispersion — can theoretically provide insights on the origin of BHs (Greene et al., 2020). However, more theoretical work is required to generate precise predictions for each potential BH formation channel, which can then be compared to observational data. Preliminary results from JWST suggest a potentially high number density of AGN (see §5.4.1), and as these systems and the host galaxies in which they reside become better understood, these observations will enable more accurate constraints on the AGN occupation fraction. One of the most striking discoveries by JWST is the potentially high $M_{BH}/M_\star$ mass ratios in a significant population of AGN across a wide redshift range ($z \ge 4$). If confirmed, these observations could suggest the presence of heavy seeds (Natarajan et al., 2024) and/or faster super Eddington or Eddington limited growth of SMBHs relative to the stellar populations of their host galaxies. However, the recent identification of a JWST broad-line AGN with $M_{BH} \approx 10^7 M_\odot$ in a $z = 7.04$ metal-poor galaxy ($< 1\%$ solar metallicity) argues against substantial super Eddington, or even Eddington limited accretion, for at least one system as the required gas inflow would drive significant star formation and hence metal enrichment (Maiolino et al., 2025). If a substantial fraction of the JWST AGN population are also found to reside in metal-poor galaxies then this would point towards more exotic BH formation scenarios such as primordial BHs to explain these systems (Maiolino et al., 2025); see §5.4.1 for more details. Therefore, at the time of writing, efficient growth of light BHs with phases of super-Eddington accretion, Eddington-limited growth of heavy seeds, and primordial origins of BHs, all remain viable options to explain the JWST AGN and it is quite likely that no one scenario explains the whole population (Schneider et al., 2023; Bennett et al., 2024; Dayal, 2024; Natarajan et al., 2024).

The seeds of BHs still fall below the sensitivity of JWST for which only AGN powered by $M_{BH} \ge 10^6 M_\odot$ are detectable and for which a SMBH mass estimate can be derived. However, by the mid-2030s the future LISA gravitational wave space antenna (Amaro-Seoane et al., 2023) will uncover all merging BHs in the mass range $\sim 10^4 - 10^7 M_\odot$ at any redshift (larger range at low redshift), while ground-based new-generation interferometers such as the Einstein Telescope (Branchesi et al., 2023) will complement the discovery space with BHs of $10 - 10^4 M_\odot$ (i.e., seed-mass BHs). These interferometers will open the door to new exciting discoveries, but will remain “silent” to non-merging BHs and will need to be completed by novel telescopes sensitive to electromagnetic emission from BH mass accretion. Next-generation X-ray telescopes such as NewAthena (Cruise et al., 2024) and AXIS (Reynolds et al., 2023), expected to be operational on a similar timescale as LISA and the Einstein Telescope, will hunt for

faint high-redshift AGN powered by low-mass BHs, potentially as low as $\sim 10^{4-5} M_\odot$, and will help confirm the nature of the high-redshift AGN candidates being uncovered by JWST.

### 5.4.4. Implication of recent discoveries on reionization

The epoch of reionization refers to the period during which the high-redshift IGM is converted from a fully neutral state to a fully ionized state. Prior to JWST, the prevailing wisdom was that reionization was primarily driven by star formation from the first dwarf galaxies, with AGN activity playing only a minor role. There are uncertainties about when reionization ended, although through the analysis of 70 quasar sightlines at $z > 5.5$, it appears to be around $z = 5.3$ (Bosman et al., 2022), and even greater uncertainty about when it exactly began (i.e., the redshift of the first Pop III stars, the first ionizing objects; Klessen and Glover, 2023). However, although numerous, the AGN candidates discovered by JWST are globally too faint to contribute significantly to reionization, with an upper limit of 30% to the ionizing budget (Dayal et al., 2024; Jiang et al., 2025).

### 5.5. Summary of the key drivers of progress

Many of the great advances in our understanding of the rapid and early growth of SMBHs over the last decade have been facilitated by significant increases in data and improvements in techniques, allied to theoretical model developments which have provided both the physical insight and clear observational predictions. The large order-of-magnitude increase in the number of identified quasars plus greater multi-wavelength datasets (from X-ray-radio wavelengths) has allowed for systematic comparisons of observed properties as a function of key physical parameter (e.g., AGN luminosity; SMBH mass; Eddington ratio) with sufficient source statistics to control for potential biases. It is now clear that the observed properties of quasars are tightly connected to the accretion-disk emission, which changes with SMBH mass and Eddington ratio (see Fig. 11), and the emission “felt” by the gas in the vicinity of the SMBH is modified depending upon opacity and orientation. This is unsurprising of course, given the accretion disk is the origin of the power of quasars, but connecting the observed properties of quasars to the underlying physics has required decades of observational and theoretical research.<fnref n="39" /> The greater source statistics have also allowed for a better characterisation of different quasar populations, critical for determining the role of environment (host and large scale) on the properties of quasars and whether quasars evolve through discrete phases, particularly when allied to careful population studies utilising control samples. These studies have mostly focused on the unobscured quasar population but, over the next decade, the operation of ambitious large-scale spectroscopic programmes targeting X-ray, infrared, and radio AGN will address this imbalance, providing obscured

<footnote n="39">³⁹It is perhaps sobering to note that if the gas in our galaxy did not absorb far-UV emission then these key discoveries may have been made decades earlier!</footnote>
quasar samples in the hundred thousands–millions (see Footnote 25).

Facilities have been a key driver of scientific progress for high-redshift AGN and quasars. Observations from JWST have expanded the redshift–luminosity discovery plane for AGN (and candidate AGN), allowing for rapid advances in our understanding of the early growth of SMBHs. The NIRSpec instrument has been crucial in driving these advances since the wide wavelength coverage, high spectral resolution, and high sensitivity allows for not only the identification of these high-redshift sources but also the determination of their physical properties. Larger field of view optical–near-IR imaging (e.g., Subaru Hyper Supreme Cam; PanSTARRS; DESI-legacy survey) and multi-object spectrographs, with increased spectroscopic multiplexing (e.g., SDSS-III, IV; DESI), have allowed for more efficient observing programmes and greater data volumes, providing unprecedented source statistics across the luminosity–redshift plane. New multi-wavelength facilities over the last decade (e.g., LOFAR in the radio; ALMA in the submm/mm; VLT-MUSE in the optical) have provided the broader multi-wavelength context to understand the physical properties of quasars and high-redshift AGN and their broader host galaxy and large-scale environments. An associated advance to these facility developments has also been greater community involvement by developing large programmes that exploit a large fraction of the available observing time but deliver data products to allow the wider astronomical community to exploit the data, driving efficient scientific progress.

Rapid advancements in computational power over the last decade have also enabled the execution of large simulations with 100 comoving Mpc box sizes, such as EAGLE, HORIZON-AGN, AND ILLUSTRIS. This has allowed for the evolution of thousands of galaxies and AGN to be captured across cosmic time (§4.1.5 & §4.2). A parallel advance has been an increase in the resolution for computationally expensive zoom-in simulations of high-redshift quasars and their environments (Costa et al., 2014b; Lupi et al., 2024a; Anglés-Alcázar et al., 2021), critical to accurately trace the physics on smaller spatial and temporal scales; see also §3.7.

**6. What impact do AGN have on the black-hole fuelling and star formation?**

The previous three sections have investigated the conditions for, and processes of, fuelling SMBH growth. As described, during episodes of SMBH growth (i.e., when the SMBH is identified as an AGN), a huge amount of energy, momentum, and mass can be released in the form of radiation, accretion disk winds, and relativistic jets of charged particles; see Figs. 1, 10, & 11. In this section we explore the current understanding of how AGN energy, momentum, and mass can produce a “feedback” effect on future SMBH growth and the properties of the host galaxies.

The co-existence of AGN and high-velocity gas in galactic nuclei has long been known (Crenshaw et al., 2003). Evidence

has sporadically pointed to large-scale outflows<sup>40</sup> (Tremonti et al., 2007; Morganti et al., 2007). However, it is only in the last ≈ 15 years that observations have shown that galaxy-wide impact driven by AGN activity, over scales ~ 100 pc – 10 kpc, is widespread. In quasars and local Seyferts, this impact takes the form of massive outflows, traced via P-Cygni-type profiles, broad emission lines and high-velocity absorption components (Fischer et al., 2010; Feruglio et al., 2010; Sturm et al., 2011; Rupke and Veilleux, 2011; Nesvadba et al., 2011; Aalto et al., 2012; Combes et al., 2013; Cicone et al., 2014).

For a SMBH lying on the local $M_{BH}-M_{\star}$ scaling relation, the cumulative energy released via accretion over the black hole’s lifetime exceeds the galactic binding energy by a factor ~ 100 (Fabian, 2012; King and Pounds, 2015). Indeed, even before the discovery of tight scaling relations linking SMBH masses and galaxy properties, energy release via AGN activity (commonly referred to as “AGN feedback”) had already been identified as a likely process to offset radiative cooling losses and expel gas from local and high-redshift, massive galaxies (Tabor and Binney, 1993; Haehnelt et al., 1998; Silk and Rees, 1998).

Two central questions underpin the ongoing debate surrounding AGN’s role in galaxy formation:

1. What fraction of the available energy is effectively transferred to the surrounding gas?

2. Which galactic component—the ISM or the CGM is most impacted?

This section reviews how our answers to these two central questions have evolved over the past decade. During this time, our understanding of AGN impact on galaxy evolution has been reshaped by parallel advances across several fields. These include a deeper (analytical) understanding of how winds and jets interact with galactic gas, the development of more sophisticated ‘sub-grid’ models for SMBH accretion and AGN feedback, higher-resolution simulations performed using new, more accurate (radiation-)hydrodynamic solvers. Additionally, the explosion of AGN-driven outflow detections—fuelled by the golden age of Integral Field Unit (IFU) spectroscopy and the maturity of ALMA—has further transformed our understanding.

Since feedback from AGN was first invoked on theoretical grounds (Tabor and Binney, 1993; Haehnelt et al., 1998; Silk and Rees, 1998), we start this section by reviewing the state-of-the-art theoretical treatment of AGN feedback in (large-scale) cosmological models of galaxy formation (§6.1). In §6.2, we ‘zoom-in’ on what has changed in our understanding of the small-scale processes driving AGN feedback. In §6.3, we assess what this progress has taught us about the impact of AGN on galaxy evolution, while in §6.4 we summarise the new observational insights into AGN feedback, reflecting on how these constrain the role of AGN in galaxy evolution.

<footnote n="40"><sup>40</sup>Here, we follow Harrison and Ramos Almeida (2024) and define an outflow as ISM/CGM material that has been swept up/entrained by some driving mechanism due to the AGN (such as an accretion disk wind, jet, or radiation pressure).</footnote>
<figure><img src="image_33.png" /></figure>

Figure 32: Despite the use of different hydrodynamic solvers, different models for star formation, stellar feedback, interstellar medium, virtually all cosmological models of galaxy evolution must invoke feedback from AGN in order to produce a realistic population of massive galaxies. According to cosmological hydrodynamic simulations, AGN feedback regulates the gas content of massive haloes and their central galaxies, modulating morphological transformations (left-hand panel), the stellar-to-halo mass ratio (central panel) and quiescent fraction of massive galaxies (right-hand panel), among a long catalogue of other galaxy properties (see text). On the left-hand panel (adapted from Fig. 2 of Dubois et al. 2016), the annotations give the halo mass (at the bottom) and the ratio of rotation to dispersion velocities for the stellar component. The central panel shows results from Wellons et al. (2023) (adapted from their Fig. 2), who compare simulations that do not include AGN feedback (coloured points), with simulations including various AGN feedback prescriptions (grey points) and with the semi-empirical relation (solid curve). The right-most plot shows the specific SFR (sSFR) for different runs the IllustrisTNG simulations (empty histograms; adapted from Fig. 1 of Terrazas et al. 2020) and observations (filled histogram). The fiducial simulation (black histogram) includes two modes of AGN feedback, which are required to produce the quiescent (low sSFR) galaxy population. The two smaller panels show example star-forming and quenched galaxies from the simulations (adapted from Fig. 3 of Terrazas et al. 2020).

## 6.1. AGN feedback in galaxy formation models

Every successful cosmological model of galaxy evolution now incorporates feedback from AGN. Despite different hydrodynamic solvers in numerical simulations, distinct models for star formation, the ISM and supernova feedback, there is broad consensus that AGN feedback is required to produce a realistic galaxy population in massive haloes with $M_{vir} \gtrsim 10^{12} M_{\odot}$ (Vogelsberger et al., 2014b; Schaye et al., 2015; Dubois et al., 2016; Weinberger et al., 2017b; Davé et al., 2019; Wellons et al., 2023; Dolag et al., 2025).

The hot, gaseous atmospheres that exist within massive haloes have long cooling times, naturally suppressing gas inflow towards central galaxies. But this suppression is insufficient without an additional source of energy such as AGN (Sijacki et al., 2007; van de Voort et al., 2011; Beckmann et al., 2017; Feldmann et al., 2023). Processes unrelated to AGN have been invoked as sources of feedback, including Type Ia supernovae, asymptotic giant branch (AGB) stellar winds, magnetic fields, (non-AGN) cosmic rays, thermal conduction, gravitational heating and morphological quenching (Martig et al., 2009; Dekel et al., 2009). Some key processes (e.g. thermal conduction) still do not feature in state-of-the-art cosmological boxes. While they all plausibly contribute to suppressing star formation in massive haloes, AGN feedback is thought to be dominant (Su et al., 2019).

### 6.1.1. Galaxy formation without AGN feedback

Models fail to capture key aspects of galaxy evolution when AGN feedback is omitted; see Fig. 32. Massive galaxies would convert too many baryons into stars. They would be gas-rich, actively star-forming, and would fail to 'quench' (Weinberger et al., 2018; Terrazas et al., 2020; Wellons et al., 2023). As

shown in the left-hand panel in Fig. 32, they would exhibit blue optical colours, while AGN feedback ensures that massive galaxies become gas-poor, passive and red (Di Matteo et al., 2005; Springel et al., 2005; Somerville et al., 2008; Gabor et al., 2011). By suppressing cooling in massive galaxies, AGN feedback sets the high end of the stellar mass- and luminosity functions (Scannapieco and Oh, 2004; Bower et al., 2006; Vogelsberger et al., 2014a). Feedback from AGN may also mediate the morphological transition from discs to ellipticals (Dubois et al., 2016; Sparre and Springel, 2017), likely by preventing discs from re-forming after major mergers. The morphology of the stellar component would also likely be affected, as the lack of AGN feedback facilitates the formation of compact stellar cores that do not as easily form in the presence of strong AGN feedback (Peirani et al., 2017; Choi et al., 2018; van der Vlugt and Costa, 2019; Cochrane et al., 2023). Differences arise also in the chemical enrichment properties of galaxy populations (e.g. the $\alpha$-enhancement of stellar populations), which can be better reconciled with observational constraints if AGN feedback is present (Taylor and Kobayashi, 2015; Segers et al., 2016). Discrepancies are evident also at the scale of the CGM. The absence of strong feedback in massive galaxies results in denser gaseous atmospheres (Davies et al., 2020a; Voit et al., 2024), producing more X-ray emission than is observed (Choi et al., 2015). At high redshift, the absence of AGN feedback would prevent bright quasars from being unobscured in the UV (Ni et al., 2020; Vito et al., 2022; Bennett et al., 2024), preventing the formation of extended $Ly\alpha$ nebulae (Costa et al., 2022).

### 6.1.2. Modelling AGN feedback

In cosmological volumes of $\gtrsim 50^3 h^{-3} \text{ Mpc}^3$, spatial resolution has improved over the past decade from ~ 1 kpc to
~ 100 pc (e.g., Nelson et al., 2019). Despite these advances, this resolution still falls short of resolving SMBH accretion and feedback processes, which occur on scales of ≲ 10⁻² pc (see Fig. 1). In state-of-the-art cosmological simulations, the various effects of AGN feedback described in §6.1 are the outcome of ‘sub-grid’ models.

The purpose of AGN feedback sub-grid models is to parametrize the effects of unresolved processes, such as winds, jets, and AGN radiation, on the scales that simulations can resolve. Since they operate at scales of ~ 100 pc, the sub-grid models used in large cosmological simulations cannot be straightforwardly connected to smaller-scale processes like jets and accretion disk winds, which are launched on sub-pc scales (see Fig. 10). As a result, these models rely on assumptions about how such interactions affect gas on much larger scales. Moreover, sub-grid models are often partly designed with the need to overcome problems of numerical origin in view (e.g. Booth and Schaye, 2009).

A common approach to simulating AGN feedback is to continuously inject thermal energy around sink particles representing SMBHs (Di Matteo et al., 2005; Springel et al., 2005) and is standard for modeling feedback in quasar phases, i.e., “quasar feedback”, in simulations like MAGNETICUM (Steinborn et al., 2015; Dolag et al., 2025), ILLUSTRISTNG (Weinberger et al., 2017b), and NEWHORIZON (Dubois et al., 2021). Yet, continuous thermal energy injection can suffer from numerical cooling losses, dispersing energy into dense gas where it dissipates before a significant effect occurs (Booth and Schaye, 2009; Weinberger et al., 2018). To address this, Booth and Schaye (2009) proposed injecting thermal energy in bursts, raising the gas temperature above a threshold ($T \gtrsim 10^8$ K) to ensure longer cooling times. This model, adopted in simulations such as COSMOOWLS (Le Brun et al., 2014) and EAGLE (Schaye et al., 2015), leads to feedback occurring in intense, periodic bursts.

Other models abandon thermal energy injection altogether and, instead, inject energy in kinetic form to represent quasar feedback. This approach can produce stronger feedback, provided the energy is not immediately thermalized (Costa et al., 2020). Injecting momentum can cause acceleration, while offsetting cooling losses (e.g., Choi et al., 2012). Notably, Choi et al. (2015) and Farcy et al. (2025) show that kinetic quasar feedback can yield a similarly realistic $M_{BH} - M_{\star}$ relation as thermal feedback, while exerting a much stronger influence on halo gas and the AGN host galaxy.

When continuous thermal injection is used to model quasar feedback, simulations indeed often rely on enhanced feedback during low accretion rate phases. A clear example is the ILLUSTRISTNG simulation, where AGN feedback at high accretion rates is modeled with continuous thermal feedback, and, at low accretion rates, is modeled through (i) the injection of kinetic energy and (ii) a pulsating rather than continuous injection. This approach is key to enabling the quenching of massive galaxies in the simulation (Weinberger et al., 2018; Terrazas et al., 2020), as shown in the right-most panel of Fig. 32. Weinberger et al. (2017b) emphasize this model remains agnostic about underlying microphysics, although wide-angle winds and jets are highlighted as possible physical motivations for their

approach.

Large-scale simulations are often performed with multiple such sub-grid models (e.g. Rennehan et al., 2024; Huško et al., 2024). NEWHORIZON (Dubois et al., 2021) and SIMBA (Davé et al., 2019), for instance, include both a ‘wind mode’, proceeding via bipolar injection of momentum (in SIMBA) or thermal energy (in NEWHORIZON), and a ‘jet mode’ producing higher-velocity, kinetically-dominated, collimated outflows. SIMBA further includes an ‘X-ray mode’ whereby a combination of thermal and kinetic energy is deposited into gas to mimic X-ray heating and radiation pressure, i.e. without explicit radiative transfer. Understanding to what extent the flows produced by such models at the current resolution scale of cosmological simulations ($\gtrsim 100$ pc) can result from the impact of jets, radiation and winds with the ISM on smaller scales is an important future milestone in the field.

## 6.2. Physical mechanisms

From a theoretical perspective, the last decade has witnessed growing consensus that AGN feedback has to proceed via ‘energy-driven’ large-scale outflows (King and Pounds, 2015). Importantly, this conclusion appears to hold regardless of whether these outflows are powered by jets or winds. This insight was brought about by advances in analytic theoretical modeling (King et al., 2011; Zubovas and King, 2012; Faucher-Giguère and Quataert, 2012; Thompson et al., 2015; Hartwig et al., 2018), higher-resolution simulations probing the impact of jets and winds on multiphase media (Wagner et al., 2012, 2013; Mukherjee et al., 2016; Ward et al., 2024; Sivasankaran et al., 2025), but also by new ‘cross-talk’ between analytic models and numerical simulations for winds (Nayakshin and Power, 2010; Costa et al., 2014a) and jets (Bourne and Sijacki, 2017; Talbot et al., 2021; Huško et al., 2022).

### 6.2.1. Energy- vs. momentum-driven outflows

A hot bubble forms when wind or jet fluid passes through a strong shock that converts the bulk of its kinetic energy into thermal energy (see Fig. 33). For winds and jets with velocities $\gtrsim 0.1c$, the temperature of the post-shock fluid can exceed ~ $10^8$ K. If cooling losses are negligible in the shocked wind/jet fluid, its energy can, in principle, be fully transferred to ambient gas. The resulting outflow is termed ‘energy-driven’ or ‘energy-conserving’. In an ‘energy-driven’ outflow, gas expands, doing ‘PdV’ work on intervening ambient gas, sweeping it up and driving it outwards at speeds $\gtrsim 1000 \text{ km s}^{-1}$ (King et al., 2011; Zubovas and King, 2012; Faucher-Giguère and Quataert, 2012; Costa et al., 2014a; Richings and Faucher-Giguère, 2018).

Much of the debate in the preceding decade has centred around whether these hot bubbles indeed retain their thermal energy or whether they cool catastrophically, to produce ‘momentum-driven’ or ‘momentum-conserving’ outflows. In a ‘momentum-driven’ outflow, the large-scale outflow is no longer accelerated via adiabatic expansion of an inner shocked bubble, relying instead on (conserved) residual momentum from the wind (King, 2003). Partially-radiative outflows also exist. Here, the shocked wind fluid cools more slowly, on a
<figure><img src="image_34.png" /><figcaption>Figure 33: Energy-driven outflow powered by a spherical wind with velocity $v = 10,000 \text{ km s}^{-1}$ propagating through an ambient medium of constant density. The top panel shows gas temperature and the bottom panel shows gas density (at a later time). The backward facing reverse shock (black circle) converts the bulk of the wind’s kinetic energy to thermal energy, producing a hot bubble that expands in an energy-conserving fashion. The denser swept-up ambient gas shell bounded at large radii by the forward shock (white circle) cools radiatively at late times (bottom panel). Source: produced from simulation results of</figcaption></figure>

Costa et al. (2020).

timescale comparable to the outflow time (Faucher-Giguère and Quataert, 2012).

Whether outflows are energy- or momentum-driven determines the strength of the outflow, its kinematics and phase structure. Energy conservation between a small-scale wind, with a kinetic wind power $L_w$, and the swept-up ISM, expelled at velocity $v_{\text{out}}$ and momentum flux $\dot{P}_{\text{w}}$, gives

$$ \dot{P}_{\text{out}} \sim L_{\text{w}}/v_{\text{out}}. \qquad (2) $$

Momentum conservation alone gives

$$ \dot{P}_{\text{out}} \sim L_{\text{w}}/v_{\text{w}}. \qquad (3) $$

If there is significant entrainment of ambient gas, the outflowing mass must slow down, giving $v_{\text{out}} \ll v_w$ and thus $L_w/v_w \ll L_w/v_{\text{out}}$. Therefore energy-driven outflows exert stronger forces on ambient gas than momentum-driven outflows (King, 2005). Indeed, observational efforts have centred on identifying outflows with enhanced momentum fluxes (see §6.4). The energy coupling efficiencies are also different. While energy-driven winds inject energy at a rate $\sim L_w$ into ambient gas, the energy coupling efficiency is lower by a factor $v_{\text{out}}/v_w = 10^{-1} (v_{\text{out}}/10^3 \text{ km s}^{-1})(v_w/10^4 \text{ km s}^{-1})^{-1}$.

It is based on momentum-driven solutions that King (2003) reproduces both the normalisation and the slope of the observed $M_{\text{BH}} - \sigma$ relation. The model predicts

$$ M_{\text{BH}} = 3.9 \times 10^8 \left( \frac{f_{\text{gas}}}{0.17} \right) \left( \frac{\sigma}{200 \text{ km s}^{-1}} \right)^4, \qquad (4) $$

in striking agreement with e.g., Kormendy and Ho (2013), who find the empirical relationship to be

$$ M_{\text{BH}} = 3 \times 10^8 \left( \frac{\sigma}{200 \text{ km s}^{-1}} \right)^{4.4}. \qquad (5) $$

The momentum-driven outflow picture has since been challenged. Faucher-Giguère and Quataert (2012) argue that two temperature effects expected in the shocked wind fluid prevent the catastrophic cooling required for a momentum-driven solution. Using hydrodynamic simulations, Costa et al. (2020) show with hydrodynamic simulations that fast winds often thermalise outside their cooling radius, which causes them not to ‘skip’ a momentum-driven phase.

Nayakshin (2014) argues that even purely energy-driven outflows may produce a $M_{\text{BH}} - \sigma$ relation as in Eqn. 4. In a medium comprising gas clouds and low density channels, energy-driven bubbles take paths of least resistance, pushing the cold gas via ram pressure (see Fig. 34). In a more complex setting such as a porous medium (Wagner et al., 2013; Ward et al., 2024) or the cosmological environment of a bright quasar (Costa et al., 2014a), energy-driven outflows couple their energy less efficiently than in homogeneous media, where their impact is maximal (Costa et al., 2014a; Ward et al., 2024). Crucially, a momentum flux $\sim L/c$ appears to be insufficient on scales $\gtrsim 100$ pc. This conclusion has been drawn repeatedly in the last $\sim 15$ years. For instance, Debuhr et al. (2011, 2012) require momentum fluxes $\sim 10L/c$ to effectively impact the galactic ISM and reproduce a realistic $M_{\text{BH}} - \sigma$ relation in simulations of major mergers. In their cosmological simulations, Costa et al. (2014a, 2018) also find that momentum fluxes $\sim 10L/c$ are required to power strong outflows in quasar-host galaxies.

When pieced together, the theoretical developments of the last decade point to AGN feedback via energy-driven bubbles. Whether mediated via (inclined) jets or winds, a low-density, hot $\gtrsim 10^9$ K, energy-conserving fluid is likely to be shaping the evolution of massive galaxies. What the last decade has also taught us is that this outflow phase is almost impossible to detect directly (Nims et al., 2015). Observations of AGN outflows typically probe much denser and colder gas phases. The recent development in our understanding of how this dense phase comes about is reviewed next.

### 6.2.2. Formation of multi-phase structure

Observations of AGN-driven outflows (see §6.4) indicate efficient entrainment of warm ($\sim 10^4$ K) and cold, molecular ($T \sim 100$ K) gas. Gas passing a strong shock with velocity $\gtrsim 200 \text{ km s}^{-1}$ shock heats the gas to a temperature $\gtrsim 10^6$ K. Collisions of ambient gas with outflows at velocities $> 1000 \text{ km s}^{-1}$ naturally produces hot gas. Central to research in the last decade has been the origin of the cool gas component. Pre-existing ambient cold gas clouds should be rapidly destroyed, before they are fully entrained in the outflow (Klein et al., 1994; Marinacci et al., 2010; Scannapieco and Brüggen, 2015; Zhang et al., 2017). The likely small ($\lesssim 10$ pc) cloud sizes places stringent resolution demands that cannot be met in cosmological boxes, requiring smaller-scale idealised simulations. In the
<figure><img src="image_35.png" /></figure>

Figure 34: In an energy-driven outflow, ambient gas is driven out via adiabatic expansion of a hot, energy-conserving bubble. These outflows impact galaxies on different spatial and temporal scales. Large-scale bubbles propagate through halo material, clearing out the low-density CGM. Within the ISM, these outflows trigger shock waves that compress the ISM. Shear flows between the energy-driven bubble and ISM clouds destroy these clouds, although, under certain conditions, cold material can reform within the outflow. Source: reproduced from Fig. 2 of Ward et al. (2024).

context of AGN-driven outflows, two main origin scenarios for cold gas outflows have emerged from these studies:

**Formation of cold gas in shock fronts:** The propagation of a strong forward shock compresses and accelerates ambient gas (see Fig. 33). Enhanced radiative cooling in swept-up ambient gas causes the condensation of a warm, dense gas phase with $T \sim 10^4$ K (Zubovas and King, 2014; Costa et al., 2014a). Using hydrodynamic simulations performed with an on-the-fly chemical network, Richings and Faucher-Giguère (2018) show that cooling proceeds down to temperatures $\sim 100$ K, at which point the outflow becomes dominated by molecular gas. In this scenario, newly-formed cold gas is comoving with the fast hot outflow and is thus less affected by destructive shear. But there are several difficulties: (i) dust, needed to catalyse the formation of molecules, may be destroyed by the strong shock and thermal sputtering (Ferrara and Scannapieco, 2016; Barnes et al., 2020) and (ii) this mechanism requires the energy-driven bubble to be effectively trapped by the ambient medium (Nims et al., 2015; Ward et al., 2024) and does not operate efficiently if hot gas vents out

(Costa et al., 2015).

**Cloud destruction, mixing and re-cooling:** Fig. 34 shows a recent example of a simulated AGN-driven outflow, produced by the collision between a fast ($v_w \sim 10^4$ km s$^{-1}$) wind and a clumpy disc (from Ward et al., 2024). There is coexistence over six orders of magnitude range in gas temperatures and densities. Here, the coldest, densest gas does not form in shocks, but within mixing layers comprising cloud debris and outflow material that cool radiatively on an outflow timescale (Marinacci et al., 2010; Armillotta et al., 2016; Gronke and Oh, 2018). While it produces fast-moving cold gas, the hot AGN-driven bubble retains most of its energy, remaining energy-conserving (Ward et al., 2024).

### 6.3. The impact of AGN on galaxy evolution

AH12 distinguish two AGN feedback modes: (i) a 'super-wind mode' (often termed 'quasar-' or 'wind' mode) which operates in luminous AGN and expels cold gas from galaxies, and (ii) a 'radio mode' which is mediated via jets and maintains halo gas hot. More physical 'sub-grid' models for AGN winds
<figure><img src="image_36.png" /></figure>

Figure 35: Simulation of relativistic low-power jet interacting with inhomogeneous disc galaxy. Such jets, like quasar winds, inflate energy-driven bubbles that carve out cavities in galactic nuclei and clear out halo material. Even in configurations where jets propagate along the disc plane, the outflow becomes bipolar on larger scales. Source: reproduced from Fig. 16 of Mukherjee et al. (2018).

and jets, and influx IFU observations conducted over the last ~10 years blur this neat two-mode division.

### 6.3.1. Feedback channels

Whether driven by jets or winds, AGN-driven outflows produce multi-faceted effects operating on different timescales:

**Central gas cavities:** At the smallest scales resolved in hydrodynamic simulations, outflows clear out material, carving low-density cavities in the galactic nucleus (Gabor and Bournaud, 2014; Costa et al., 2014a, 2020; Torrey et al., 2020; Raouf et al., 2023; Mercedes-Feliz et al., 2023). The size of these cavities, which is of order ~100 pc, is sensitive to AGN wind strength and the AGN outburst duration. The size of central cavities is also very sensitive on numerical resolution (Bourne et al., 2015; Curtis and Sijacki, 2016) and exactly how the AGN feedback is modelled.

**Ejection of the ISM:** Outflows take paths of least resistance (Costa et al., 2014a; Bieri et al., 2017; Nelson et al., 2019). Even if energy is communicated in the form of a collimated jet at small scales, outflows on larger scales are shaped by the properties of the ambient medium through which they propagate (Wagner et al., 2012; Tanner and Weaver, 2022). This means it is possible for jets propagating along the plane of a galactic disc to produce bipolar bubbles (Fig. 35) that escape along the disc poles (e.g. Mukherjee et al., 2016, 2018; Talbot et al., 2022). The typical outcome of simulations following the propagation of jets and winds in systems comprising galactic discs is a weak ejective feedback effect (Gabor and Bournaud, 2014; Costa et al., 2020; Torrey et al., 2020). In configurations where gas has a high covering factor around the AGN (e.g., the aftermath of a merger), outflows should clear out

larger ISM masses (Faucher-Giguère and Quataert, 2012), although it remains unclear whether feedback via ejection of the ISM will become dominant even in this case (e.g. Costa et al., 2018; Mercedes-Feliz et al., 2023).

**Clearing out the halo:** Once they break out of their host galaxy, AGN-driven outflows push into tenuous, galactic halo gas (Fig. 36). Here the effect of outflows is to reduce the density of the CGM and/or raise its temperature (Anglés-Alcázar et al., 2017; Costa et al., 2018). In simulations, the expulsion of the CGM coincides with quenching (Davies et al., 2020a; Costa et al., 2020; Oppenheimer et al., 2020; Terrazas et al., 2020; Voit et al., 2024). After the CGM is cleared out, accretion from the halo stalls and star formation in the central galaxy drops due to reduced availability of gas fuel on a timescale of ~100 Myr–1 Gyr (Oppenheimer et al., 2020),

### 6.3.2. 'Radio' vs. 'Quasar' and 'Ejective' vs. 'Preventive' feedback

Mechanisms such as quasar-powered winds produce both 'ejective' and 'preventive' feedback<fnref n="41" /> (Costa et al., 2020; Zinger et al., 2020). Winds powered at small scales produce energy-conserving bubbles that expand both into the ISM and the CGM (Davies et al., 2020a). As they propagate into the more isotropic and tenuous halo, these bubbles isotropise and push out significant amounts of halo gas. The density decrease caused by halo gas ejection leads to longer cooling timescales, which suppress long-term gas accretion onto the central galaxy, producing strong preventive feedback (e.g. Davies et al., 2020a). The tendency for energy-conserving bubbles to escape through low density channels means that it is even possible that quasar-mode feedback is predominantly preventive. This question remains open to debate. Answering it will require improvements both in the modelling of the ISM and AGN feedback processes in numerical simulations.

At low Eddington ratios, feedback was pictured to be communicated via radio jets. Their high degree of collimation was assumed to imply that jets produce no effect on the interstellar medium. Instead, the reasoning went, they interact primarily with halo gas, where they inflate over-pressurised cavities that expand, do work and heat halo material. The effect was thought to be mainly 'preventive'. While the impact of jets on suppressing cooling in hot haloes is established (e.g. Weinberger et al., 2017a), jets, when highly inclined with respect to the disc, can have a similar effect as quasar winds (Wagner et al., 2012; Mukherjee et al., 2018; Tanner and Weaver, 2022; Talbot et al., 2022), powering bubbles that clear out gas from the nucleus, but expand to take paths of least resistance. This is especially true for low velocity, low power jets, which are more easily confined, and deposit their energy within the host galaxy ISM, where-as higher power jets deposit most of their energy

<footnote n="41">⁴¹The term 'ejective' feedback refers to the removal of material from the ISM, whilst the term 'preventative' feedback refers to regulating the rate at which gas can accrete onto the host galaxy from the wider halo, to ultimately form part of future ISM.</footnote>
<figure><img src="image_37.png" /><figcaption>Figure 36: In cosmological simulations, massive galaxies quench when AGN-driven outflows clear out the CGM. Shown here is the fraction of total halo mass in the CGM as a function of halo mass in the state-of-the-art simulations EAGLE (left-hand panel) and ILLUSTRiSTNG (right-hand panel). The colour gives specific SFRs. This drops when the CGM mass fraction drops as well as when the SMBH mass increases (not shown here). Source: adapted from Fig. 2 of Davies et al. (2020a).</figcaption></figure>

in the surrounding halo. As Wagner et al. (2016) put it, the interaction of young jets in a radio galaxy with the ISM leads to quasar-mode, energy-driven, mechanical, negative or positive feedback. In summary, quasars winds and jets can provide both preventative ('maintenance') and ejective feedback, and it is not trivial to separate observed AGN populations into contributing to distinct feedback modes (see discussion in Harrison and Ramos Almeida, 2024).

## 6.4. Observational work on AGN feedback

As recently reviewed by Harrison and Ramos Almeida (2024), there has been considerable effort over the last three decades to constrain the mechanisms of AGN feedback, measure AGN-driven outflow properties, and search for evidence of feedback 'in action'. This requires a multi-faceted set of observations, covering detailed multi-wavelength observations of individual targets (e.g., to understand the physics of how AGN drive multi-phase outflows) through to population studies (e.g., to look for global evidence of AGN feedback). As detailed below, this work has lead to considerable progress over the last decade, in particular, in understanding the challenges involved in interpreting observations of AGN-driven outflows, and finding 'smoking gun' evidence for global AGN feedback.

### 6.4.1. Searching for the impact of AGN feedback

Despite many earlier attempts to directly search for the impact of feedback on galaxies hosting luminous AGN (see review in Harrison, 2017), the last decade of work has resulted in the consensus, that there is no strong reason to assume the impact of AGN feedback on a star-forming galaxy's global properties is observable, while the SMBH is still active. This is, at least partly, because AGN feedback also suppresses BH growth and

the AGN will turn itself off, potentially before impacting upon the host galaxy. In the CGM expulsion scenario, it also takes > 100 Myr for the effects of AGN feedback to become evident. Consequently, there will be no correlation between the state of the AGN at that time with the episode that brought about quenching. The gas reservoir will indeed become depleted, long after a bright AGN has switched off. Indeed, state-of-the-art cosmological simulations like ILLUSTRiSTNG or EAGLE typically predict that the brightest AGN reside in gas-rich star-forming galaxies, while gas-poor, quenched galaxies have faint AGN (Scholtz et al., 2018; Piotrowska et al., 2022; Ward et al., 2022). Therefore, the fact that gaseous reservoirs are not typically suppressed around the brightest AGN is not in contradiction with theoretical predictions. According to these models, it is the SMBH mass, which is proportional to integrated injected feedback energy, that predicts whether a galaxy is quenched (Terrazas et al., 2016). Indeed, the signatures of AGN feedback on a population level are likely better found in the galaxy population as a whole, irrespective of a currently active SMBH. The latest results from JWST of very early (z ~3–4), massive quenched galaxies, are an interesting population to test against different models for AGN feedback in relation to the timescale of growing the SMBHs and shutting down star formation (e.g., D'Eugenio et al., 2024; Nanayakkara et al., 2024; Scholtz et al., 2024).

On the other hand, localised and transient feedback effects on host galaxy properties can be seen in specific galaxies hosting in-situ luminous AGN. For example, a local region of depleted molecular gas, or a change in the distributions of young stellar populations in the vicinity of an outflow or jet is sometimes observed (Harrison and Ramos Almeida, 2024). Nonetheless, it is likely that multiple AGN episodes would ultimately be re-
sponsible for *global* star formation quenching. Furthermore, a global role of regulative AGN feedback in the most massive galaxies, is observationally evidenced by combining radio observations (as a tracer of AGN) and X-ray observations of the intracluster, or intragroup medium. Observations over the last decade have confirmed earlier work (and expanded to higher redshift), that the required amount of heating is available from AGN to offset cooling for massive early-type galaxies (McDonald et al., 2019; Hardcastle et al., 2019).

### 6.4.2. *Observed multi-phase outflows*

Theoretical models naturally predict a multi-phase gas outflow structure (see Table 3 for the different gas phases). However, it is difficult to know if the different *observed* outflows are all related to the same underlying cause. Larger statistical studies by necessity often focus on a single survey / waveband, but there have been attempts to build a more holistic picture of multi-phase outflows in AGN. For example, Fiore et al. (2017) investigated outflow properties for molecular, ionised, broad absorption line (BAL; see Footnote 30), and X-ray absorption line outflows (Ultra-Fast Outflows; UFOs; see Tombesi et al., 2010) for 94 AGN. However, the same objects are not covered by the different outflow phases. Fig. 37 shows the maximum velocity of the different types of outflows in this sample (also including new BALQSO measurements), as a function of the dust sublimation radius with point size scaled by bolometric luminosity (K. Leighly, priv. comm.). There are clear correlations between the outflow and AGN parameters, but the relationships may differ between the various types of outflows.

Until recently, nearly all outflow studies focused solely on a single tracer, providing an incomplete picture of the multi-phase content of the outflows (Cicone et al., 2018). Even among single tracer studies, there may be more than one causal factor driving the outflows. For example, using SimBAL, Choi et al. (2022a) demonstrated that BALQSOs split into two subpopulations with distinct outflow properties when considering emission line properties (i.e., using Eigenvector 1 as a proxy for Eddington ratio; Boroson and Green, 1992). These single-tracer studies are still the only viable way to effectively undertake large statistical studies of outflows (e.g., Mullaney et al., 2013; Förster Schreiber et al., 2014; Avery et al., 2021). This is because the cross section of observations of different phases of outflows drastically reduces the number of sources with information available at the required wavebands.

It is critical to have spatially-matched information across multiple different phases of outflows, to understand how they are related, and gain a more complete understanding of their energy and mass content. This is expensive in terms of telescope time and person power, since the analyses require expertise across a wide range of techniques. However, the last few years has seen some progress in this area, typically by studying individual sources with detailed multi-wavelength information. These have helped pinpoint many important aspects of the properties and drivers of multi-phase outflows. For example, ISM modifications that depend on radio jets (e.g. Girdhar et al., 2022), or where there is spatial correlation between tracers of different phases of outflows (Lansbury et al., 2018; Ven-

turi et al., 2023; Speranza et al., 2024). Advances have recently been made towards building up larger samples of carefully selected samples with multiple tracers of outflow phases (e.g., Fluetsch et al., 2021; Yamada et al., 2021; Riffel et al., 2023; Speranza et al., 2024), but it is yet unclear how much these studies can be extrapolated to the general population.

Understanding how efficiently an AGN can couple their energy output to the host galaxy ISM is crucial for assessing the potential impact on the host galaxy, and testing different feedback models. However, a lot of work over the last few years has highlighted that with the limitations of typical observations, it is extremely difficult to measure outflow properties, such as energetics, and momentum rates reliably (e.g., Rose et al., 2018; Davies et al., 2020b; Lutz et al., 2020; Harrison and Ramos Almeida, 2024). When simulations are to mimic observational experiments, many further challenges become apparent, including the potential to miss most of the energy content (without a tracer of the hottest phases) and to mis-interpret energy-driven outflows as momentum-driven outflows (Ward et al., 2024). Nonetheless, the future is promising, as cutting-edge simulations, can now be used more directly with observations to help determine what can, and can not, be inferred from specific observations that will have incomplete information (e.g., limited gas phases, limited spatial information etc.).

Understanding the physical drivers of multi-phase outflows (e.g., jets; radiation pressure; accretion-disk winds) is also challenging observationally, but improvements in observations have resulted in significant progress. Recently, the advent of more sensitive radio surveys like the LOFAR Two-metre Sky Survey (LoTSS; Shimwell et al., 2019, 2022) has allowed studies to expand at least in terms of including radio information (Rankine et al., 2021; Petley et al., 2022; Kukreti and Morganti, 2024). For a typical synchrotron source, LoTSS is almost an order of magnitude deeper than the Faint Images of the Radio Sky at Twenty cm (FIRST; Becker et al., 1995), and delves into the canonically ‘radio-quiet’ population (see § 4.1.4. While this does not yet provide spatially matched information, LoTSS is currently undergoing post processing to provide sub-arcsecond resolution imaging of sources with $S_{144MHz} > 10$ mJy (Morbito et al., 2022a). Images where radio emission is spatially resolved can help provide information on the interaction between jets and multi-phase outflows (e.g., Tadhunter et al., 2014; Girdhar et al., 2022; Venturi et al., 2023). Advances in simulations are also helping provide direct observational predictions for different outflow driving mechanisms (e.g., Meenakshi et al., 2024; Ward et al., 2024).

The last decade has seen a shift towards more multi-wavelength studies which look at different outflow gas tracers in the bigger picture context, but still for relatively small samples. Large statistical studies are limited to one or two outflow tracers at a time. To really make progress on how outflows in different phases are related to each other and their origin, we need larger samples with spatially resolved information on multi-phase outflows, and more direct comparisons with simulations to try to unpick their underlying cause(s).
<table><caption>Table 3: Gas phases in which outflows are observed, and their associated tracers and properties.</caption><thead><tr><th>Gas phase</th><th>Temp [K]</th><th>Density [particles cm<sup>-3</sup>]</th><th>Typical tracers</th><th>Scales</th></tr></thead><tbody><tr><td>Hot ionised</td><td>10<sup>5</sup> - 10<sup>8</sup></td><td>10<sup>6</sup> - 10<sup>8</sup></td><td>X-ray absorption lines</td><td>&lt;1 pc</td></tr><tr><td>Warm ionised</td><td>10<sup>3</sup> - 10<sup>5</sup></td><td>10<sup>2</sup> - 10<sup>5</sup></td><td>emission lines in UV, optical, NIR</td><td>&lt;1 pc to ~1 Mpc</td></tr><tr><td>Neutral atomic</td><td>10<sup>2</sup> - 10<sup>3</sup></td><td>10<sup>1</sup> - 10<sup>2</sup></td><td>H<sub>I</sub> absorption, Na I D absorption, [C II] emission lines</td><td>~1 pc to ~1 Mpc</td></tr><tr><td>Molecular</td><td>10<sup>1</sup> - 10<sup>3</sup></td><td>&gt;10<sup>3</sup></td><td>OH, CO, H<sub>2</sub> lines</td><td>~1 pc to ~1 Mpc</td></tr></tbody></table>

<figure><img src="image_38.png" /><figcaption>Figure 37: Top: Outflow velocity as a function of normalized radius (dust sublimation radius $R_{sub} = 0.2L_{46}^{1/2}$ pc; Laor and Draine, 1993). Grey, purple, blue, and yellow markers show objects from Fiore et al. (2017), specifically, ultrafast outflows observed in the X-ray band, high-ionization BALQSO outflows from the literature, ionized outflows detected in emission lines, and outflows detected via molecular transitions, respectively. Magenta circles show results from analysis of low redshift FeLoBAL quasars Choi et al. (2022a) using SimBAL Leighly et al. (2018). The size of the point scales with the bolometric luminosity, which spans $43.1 < \log(L_{bol}/\text{erg s}^{-1}) < 48.0$. The black circles show objects with inferred $L_{bol}/L_{KE} > 0.5\%$ that have the potential for significant host-galaxy impact through AGN feedback (Hopkins and Elvis, 2010). Bottom: Number of outflows per decade of $R/R_{sub}$. Low-redshift FeLoBAL outflows span the largest range, and appear to be concentrated near the torus and at larger distances from the central engine. Kindly produced by K. Leighly.</figcaption></figure>

## 6.5. Summary of the key drivers of progress

From a theoretical point of view, the last decade has seen a significant transition from a simplified approach and understanding of AGN feedback, to more sophisticated models focused on understanding the underlying physical mechanisms, advances driven by improvements in facilities and techniques. Early models and simulations implemented simplistic AGN feedback prescriptions as a required means to present for the formation of too many massive galaxies. However, advances in both understanding and simulation development has seen more focus on testing various physically-motivated AGN feedback prescriptions (e.g., radiative driving on dust, accretion disk winds, low power jets). Furthermore, the field is moving on from assuming two simplistic, and isolated, modes of feedback

(i.e., an “ejective” mode associated quasars and a “preventative” mode associated with radio lets from low accretion systems). It is now clear that multiple feedback channels are possible, from a variety of different mechanisms, across different accretion rates. Furthermore, sophisticated modelling of AGN-driven outflows now suggest that energy-driven outflows are more relevant than momentum-driven outflows, with both accretion-disk winds and (low power, inclined) jets able to produce similar effects.

The predicted dominant hot phase of outflows remains a challenge to test directly with observations (it is too tenuous); yet the last few years are seeing more attempts to use simulations to directly predict observational tracers of the sub-dominant colder phases (e.g., by incorporating additional models to enable predictions of emission lines from warm ionised and cold molecular gas phases). The general “closing of the gap” between observations and simulations is proving to be crucial for both understanding the limitations of observations, and testing specific AGN feedback models.

From an observational point of view, the explosion of multi-wavelength data, through a greater collection of data particularly IFU spectroscopy with new facilities (Gemini; VLT-MUSE; VLT-VIMOS), has driven progress. A key facet of these advances, is moving to a more complete census of typical (more common) AGN populations, such as those which are ‘radio quiet’, facilitated by a new generation of radio surveys such as those undertaken by LOFAR. Another key facet, is the ability to obtain a detailed (sensitive and spatially-resolved) multi-phase picture of the ISM. Whilst still in the early days, we are already having a much better census of how AGN impact on the multi-phase ISM across a diversity of AGN samples (including extreme quasars, those hosting low power jets, and those with low power winds). Nonetheless, work is still required to understand how these measurements can definitely inform the long-term impact on the ISM in these host galaxies, as we overcome challenges with systematic uncertainties for moving from measured to derived properties.

The last decade has also seen the community move on from assuming there is a simple “smoking gun” piece of evidence for suppressed star formation directly associated with luminous AGN, a conceptual advance. A greater understanding of the likely timescales involved has seen more work searching for indirect evidence of AGN feedback, such as studying the population of non-active quenched galaxies. Using such population-level observations of galaxy populations, to test specific predictions (which are not ‘tuned’) from different AGN feedback models will be an important area of progress over the coming
years.

# 7. Discussion of the scientific progress

## 7.1. Reflections on the last decade of progress

The last decade has witnessed huge advances in our understanding of “what drives the growth of black holes?”, as detailed in the last four sections (see §3–§6). Here we briefly summarise some of the key areas of significant progress, making direct reference to the progress drivers defined in Table 1.

New *facilities* (category 1) have been a strong driver of progress in each of the scientific areas, allowing many key conceptual advances. For example, new facilities have allowed for significant progress in shaping our understanding of the nature of the obscuring medium in AGN (ALMA; VLTI-GRAVITY; VLTI-MATISSE; see §3.3), the time-variable nature of accretion from close photometric and spectroscopic monitoring (ZTF; ASAS-SN; spectroscopic follow up; see §3.2 & §3.5), the identification of faint high-redshift candidate AGN (JWST; see §5.4), the nature and impact of AGN-driven outflows and feedback (LOFAR and many new IFU facilities; see §6.4), and the development of higher spatial and temporal resolution simulations with improved sub-grid physics (see §3.7 and §6.1–6.3), amongst other discoveries. It is unsurprising that significant progress has been made through the exploitation of new facilities since astronomy is, largely, a facility driven discipline and facilities are defined and approved on the very basis that they will address the key scientific issues in the research field.

Increases in the quantity and quality of *data* (category 2) from existing facilities and the development of new *techniques* (category 3) have also been significant drivers of scientific progress. These two drivers often go in parallel since the increase in data will often necessitate the employment of new techniques to handle the larger datasets. For example, the increase in the availability of optical spectroscopic data from the SDSS survey, along with developments in the modelling of the physics of accretion disks, has enabled greater insight on how the observed properties of quasars connect to their underlying accretion properties (see §5.2). Our understanding that AGN are typical “events” in the broader lifetimes of galaxies has also largely come about from the development of techniques to take account of limitations in the Eddington-ratio parameter space probed by existing data, along with the realisation that the different variability timescales of AGNs and galaxies (itself assisted by improved hydro-dynamical simulations of accretion disks), which showed that different *apparent* correlations will be identified (or not) depending on what observed quantities are compared (see §4.2 & §4.5). The development of simulation tools to directly predict the observed properties of AGN-driven outflows, another example of techniques, also greatly helped the interpretation of outflow and “AGN feedback” signatures from spatially resolved spectroscopy, and ultimately showed that the impact of outflows on the host-galaxy environment is more localised and less “immediate” than the community initially assumed (see §6.4).

The last two scientific drivers, *community* (category 4) and *conceptual* (category 5), are often more difficult to assess. The

definition of *community* is broad in this context and essentially refers to approaches that aim to engage and receive feedback from a broader community rather than from an individual research team. A “pure” example of community driven progress is the Galaxy-Zoo project (Lintott et al., 2008), which engaged the general public (“Citizen Science”) to assist in the classification of the morphology of galaxies, which led to the identification of unusual systems such as the “Voorwerpje”, now understood to be the light echoes of past AGN episodes (see §4.2.6). A more typical example of community driven progress is the organisation of scientific workshops (occasionally multi or inter disciplinary; Burns et al., 2025) that bring together groups of scientists with complimentary skills and scientific backgrounds to facilitate broad discussion along with the development of new approaches to analyse and interpret data. As mentioned in §2.3, the Lorentz Center workshops are good examples of this type of community engagement, as are the Aspen workshops which require participants to commit to attending for an extended duration to allow for more in-depth discussion; community engagement was also a major driver behind the development of our series of black-hole workshops (see Appendix A). Many examples of key conceptual advances have been highlighted above. Perhaps the key factor with this progress driver is that these advances are typically achieved in association with one of the other four drivers rather than being a “Eureka” moment of pure inspiration without any other input. Indeed, the examples highlighted above represent conceptual advances but were driven by new *facilities*, increases in existing *data*, the development of new *techniques*, or facilitated through broader *community* engagement.

## 7.2. A view from citation network analyses

Citation networks are powerful tools to track scientific progress through the identification of research clusters, such as those defined in Fig. 3. These clusters are essentially research “communities” (i.e., the community with which an individual researcher would likely identify). However, these clusters are not static but change over time as new papers are published, citation patterns shift, and new research connections are made. We can employ the citation network analysis tools introduced in §2.1 to track how the clusters shown in Fig. 3 have progressed and evolved with time, a measure of community progress.

In Fig. 38a we show the changes in the number of publications per year for the 13 research clusters over the 2003–2023 period. A number of features are evident from this figure. First, over decadal-long timescales, the number of publications within each research cluster have increased: the typical increase from the 2003–2012 period to the 2013–2023 period is a factor of ≈ 1.5. Three clusters (C0, C2, and C9) have clearly undergone strong growth over the last decade with particularly rapid increases in the number of papers since 2015. Each of these clusters are tightly connected to fundamental BH research and the identification of the first gravitational waves, and the consequential need to interpret and understand these data are strong drivers of the impressive growth in these clusters. However, it is notable that the majority of the papers in these clusters are
<figure><img src="image_39.png" /><figcaption>Figure 38: (a) Annual number of publications per year for the 13 research clusters identified in §2.1 over 2003–2023 using the same colour coding as that adopted in Fig. 3; see Table B.5 for the cluster labels. Several clusters with a strong rise in publication rate since ≈ 2015 are highlighted along with C1, the cluster most closely associated with the research in this review. (b) Total number of publications for each research cluster (“Full”: blue), along with the number of publications in the six selected “astrophysics” journals (“6-Journal”: orange; see Footnote 42).</figcaption></figure>

not published in typical “astrophysical” journals,<fnref n="42" /> as shown in Fig. 38b: the researchers in these clusters probably feel more associated to theoretical physics and particle physics than astronomy and astrophysics. Although less clear from Fig. 38a, C12, which focuses on transient and TDE-related research, has also grown rapidly over the last decade almost certainly due to the increased availability of transient photometric data from the PTF, ZTF, and ASAS-SN facilities. The C1 cluster, that is most closely associated with this review, has grown at the median rate over the last decade. However, it underwent much more dramatic growth in the previous decade likely driven by the discovery of black-hole–galaxy scaling relations in the mid-late 1990’s and the demonstration that “AGN feedback” is required in cosmological simulations to explain the observed properties of galaxies in the early-mid 2000’s.

In Fig. 39 we have sub divided the CND from Fig. 3 into two
decadal-long periods: 2003–2012 and 2013–2024. The larger
clusters (C0–C6) were already well defined in the first decade
and have kept a similar morphology (a measure of clustering
strength) through the last decade. However, it is notable that C0
was originally much more disparate and offset from the other
clusters but has become better connected over the last decade.
This is largely due to the detection of gravitational waves and
the EHT constraints on the physical properties of the SMBHs
in Sgr A* and M 87, which bridge interests in both theoreti-
cal physics and astronomy; indeed, see the relative “domains”
of the EHT papers on Sgr A* and M 87 in Fig. 3. C9 has be-
come much better defined over the last decade, reflecting it’s
growth from originally a small cluster to now a medium size
cluster (see Fig. 38a). Most of the other smaller clusters are
more diffuse and filamentary but reside in approximately the
same positions across both CNDs. The exception is C8, which
has a distinct new extension towards C0 due to the SMBH con-
straints on Sgr A* from the EHT.

We finish this brief analysis with a couple of related notes.

First, as shown by analyses of the publication output across Physics over the last century, increases in the number of published papers are driven by increases in the number of researchers working within a given cluster, rather than individual researchers becoming more productive (Sinatra et al., 2015; Fanelli and Larivière, 2016). Of course, there may be some variations between sub fields and clusters, although overall any significant rises in publication output over the median trend will be driven by researchers “flocking” to these research fields, possibly due to greater funding opportunities. Second, we should be cautious about forecasting future growth based on any dramatic increases in publication within a cluster as they are unlikely to be long lived and may even decline with time as other areas of research become more appealing to both researchers and funding agencies.

7.3. *Anticipated progress over the next decade*

Identifying the drivers of scientific progress over the last decade not only helps us retrospectively understand where we have been able to make progress and why, it also helps us to provide guidance for generating future progress more efficiently. Facilities are a major driver of the scientific progress in astronomy and the largest new facilities are decided and defined a decade prior to starting operations so we already have a good idea of what the next decade will bring in terms of facilities.

A non-exhaustive list of expected major new facilities in-
cludes several extremely large ground-based optical–infrared
telescopes with > 30 m diameter primary mirrors (e.g., the Eu-
ropean ELT; the Thirty Meter Telescope), the square-kilometer
radio array (SKA), the *Euclid* and *Roman* near-IR space tele-
scopes, the *Athena* and *AXIS* X-ray telescopes, the Cherenkov
Telescope Array (CTA) Gamma-ray facility, the *LISA* and *Ein-
stein* gravitational-wave observatories, ground-based large ded-
icated multi-object spectroscopic telescopes (DESI; MOONS;
PSF; 4MOST; WEAVE), and the LSST survey with the Vera
Rubin telescope. These are just some of the major new facili-
ties and doesn't include the current facilities which are already

<footnote n="42">⁴²In this analysis we used the following six journals: *Mon Not R Astron Soc*; *ApJ*; *ApJL*; *Astron. J*; *Astron. Astrophys*; and *New Astron. Rev*.</footnote>
<figure><img src="image_40.png" /><figcaption>Figure 39: The temporal evolution of the citation network diagram presented in Fig. 3, split into the last decade (2013–2024: bottom) and the prior decade (2003–2012: top), with the 13 identified clusters colour coded and highlighted; see Table B.5 for the cluster labels.</figcaption></figure>

making great advances and yet still maturing such as the LOFAR radio telescope, the eROSITA observatory, and the JWST space telescope. In terms of computing facilities, graphics processing units (GPUs) are likely to become more prevalent over the next decade offering much greater performance than central processing units (CPUs) particularly for parallel processing. Artificial intelligence (AI)-based processors will also likely become available, along with improvements in the capacity and performance of data storage and memory.

We can expect great discoveries and potential conceptual advances from all of these facilities, which open up new parameter space across the broad frequency–sensitivity–area–source statistic plane. Given the prominent growth rates of clusters C0, C2, C9, and C12 we highlight LISA which will provide the first sensitive constraints on gravitational waves from SMBHs and the Vera Rubin telescope which will fully open up the astronomical community to unprecedented large-scale sensitive multi-band optical photometric monitoring and the identification of transient phenomena. On the simulation side, developments in computing architecture over the next decade should allow for the first simulations on exascale supercomputers (> $10^{18}$ double-precision operations/s) and provide at least an order of magnitude improvement in the box-size–

simulation resolution ratio over current state-of-the-art simulations. Improvements in sub-grid modelling will be possible thanks to enhanced numerical resolution and to closer, cross-scale connections to ever more comprehensive and physically-accurate small-scale simulations. Such work is already underway (Costa et al., 2020; Rennehan et al., 2024; Koudmani et al., 2024; Huško et al., 2024) and future large-scale models will certainly tap into predictions from small-scale radiative transfer and general-relativistic-magnetohydrodynamic simulations (Ricarte et al., 2023; Kaaz et al., 2024), simulations performed with more accurate gravity solvers (Genina et al., 2024), and more *ab initio* simulations of BH seed formation (Regan and Volonteri, 2024). Innovative machine learning approaches to sub-grid modelling, currently in their infancy, may come to play an important role in porting small-scale simulation results and predictions to large-scale, coarse-resolution simulations (e.g. Stachenfeld et al., 2021).

It is also undeniable that astronomy has become more data rich over the last decade with the current generation of facilities already providing such vast quantities of data that it is becoming impractical for any individual to store all of their data (e.g., LOFAR generates tens–hundreds of terrabytes of data per day; Manzano et al., 2024). In some cases it has even become impractical for the observatory to retain all of the raw data as it is more effective to simply re-take the observations. This increase in data will drive further development of data-analysis techniques capable of handling such large datasets. Over the last decade, we have seen a large increase in the adoption of machine-learning approaches to analyse large datasets (see Footnote 3), which will surely become even more prevalent over the next decade. Given the massive rise in the ability of generative artificial intelligence (gAI) over just the last few years, it also seems inevitable that gAI will play a more dominant role in data analyses over the next decade.

A great, largely untapped, potential for significant scientific progress over the next decade is improved interactions between different research communities, potentially through multi or inter disciplinary workshops (Burns et al., 2025). The CNDs shown in Figs. 3 & 39 show connections between different “black hole” research communities, some of which are comparatively isolated and could benefit from having broader scientific discussions with other research communities. This can be easily achieved through the facilitation of focused workshops, aimed to address key questions by bringing together experts from different communities (even beyond the astronomical community) who can bring new knowledge, skills, and approaches to address key questions. Compared to the large budgets and long timescales required to develop new facilities, focused workshops provide extremely cost effective and immediate solutions to making significant scientific progress.

## 8. Conclusions

In this review we have provided an overview of our understanding of the growth of SMBHs to address the broad question: what drives the growth of black holes? This question is somewhat rhetorical as, has been established for several decades, the
growth of black holes is the result of the accretion of gas onto the black hole. However, the broader implication of this question is how does the gas reach the vicinity of the black hole, what was the earlier origin of that gas, are there are specific environments more conducive to black-hole growth than others, and what influence does the black hole have on the evolution of the host galaxy? We have shown that huge progress has been made in addressing each of these questions over the last decade.

We now have a much greater physical understanding of the accretion process, driven by both theoretical and observational advances, and have opened up the time-domain axis in AGN astronomy, identifying instances of changing accretion states. Our understanding of the obscuring medium has undergone a revolution and led to the development of more physical models of the AGN “torus”, largely driven by high spatial resolution observations as well as a better understanding of the role of galaxy scale structures in AGN obscuration. We have a better understanding of the connection between the AGN and its hosting galaxy and larger-scale environment and have demonstrated that AGN are “events” within the broader life cycle of galaxies which can significantly impact the long-term evolution of the galaxies through energetic “outflows” (accretion driven winds and jets). We also have a significantly greater census of AGN activity across the luminosity-obscuration-redshift plane, with sample sizes now numbering in the millions, and we are starting to witness the very earliest phases of SMBH growth.

We have also evaluated the drivers of this scientific progress, using the AH12 review as a “snapshot” of the field from a decade ago. Inevitably much of the scientific progress is driven by new facilities and more abundant data across the whole of the broad population in luminosity, redshift, and wavelength space, which has required the development of new techniques to analyse the ever increasing data volume. This has led to several conceptual advances such as the examples given in the previous paragraph. We have also shown how scientific progress can be driven through community developments, both in citizen science and the professional research community, and utilised citation network analyses to demonstrate the full breadth of the black-hole research network. A key component of community engagement is scientific workshops such as our Durham–Dartmouth black-hole workshops, which form the backbone of this (and the original AH12) review. With that final concluding thought we thank you for reading this review and we hope to see you at the next black-hole workshop focused on understanding “What drives the growth of black holes?”.

## Acknowledgements

We dedicate this review to the memory of Richard Bower who sadly passed away on 6th January 2023. Richard was not only a brilliant scientist and polymath, making significant contributions to theory, observation, and instrumentation in astronomy and beyond, he was also a warm and entertaining colleague and friend. Richard was also one of the greatest patrons of our black-hole workshop series, missing just two meetings, lighting up a session with his knowledge, passion, and personality. His influence was particularly great in our 2012 workshop on “What

is the role of AGN in the evolution of galaxies?” at Dartmouth College where he became the de facto arbiter on the outcome of scientific debates, giving his thumbs up or thumbs down in agreement or not (pollicie verso)!

We thank the 187 participants of our Iceland meeting in 2022 for the great contributions to the discussions and debates, which reinvigorated our passion for scientific conferences following 2 years of online meetings during the Covid-19 pandemic. We also want to thank everyone who has attended at least one of our Durham–Dartmouth black hole meetings over the last decade, particularly those who contributed to the great discussions which have helped to guide the field and lead to scientific progress. A special thanks to Roger Blandford FRS, our surprise attendee to the Iceland meeting, for presenting the awards to our long-serving patrons at the Iceland meeting and for the great after-dinner speech on the drivers and progress over the last 50 years of black-hole research.

We thank the referee for constructive and positive feedback and the following people for commenting on earlier drafts of this review and/or providing data and figures: Alastair Basden, Xiaohui Fan, Vicky Fawcett, Claire Greenwell, Stephanie Juneau, Dale Kocevski, Rohit Kondapally, Karen Leighly, Rusen Lu, James Matthews, John Paice, Krisztina Perger, Cristina Ramos-Almeida, Amy Rankine, Jan-Torge Schindler, Jan Scholtz, Matthew Temple, Carolin Villforth, and Sarah White. We extend a special thanks to Matthew Temple for producing Figure 27. DMA thanks the Science Technology Facilities Council (STFC) for support through the Durham consolidated grant (ST/T000244/1). JA, CMH, and LKM acknowledge funding from United Kingdom Research and Innovation (UKRI) grants MR/T020989/1, MR/V022830/1, and MR/Y020405/1, respectively. MH acknowledges funding from the Swiss SNSF Starting Grant (218032). RCH acknowledges support from NASA via grant numbers 80NSSC22K0862 and 80NSSC23K0485, and is grateful to the Durham Department of Physics and Center for Extragalactic Astronomy for hosting a very enjoyable and productive sabbatical visit. We thank the Matariki Network of Universities, of which both Dartmouth College and Durham University are partners, for financial support of our institutional meetings.

## Appendix A. Durham–Dartmouth black-hole workshops

Over the last decade we have held six black-hole growth workshops, our Durham–Dartmouth black-hole workshop series. Brief details on each workshop, including the number of participants, are given in Table A.4. Overall, 512 participants have attended at least one of our workshops. A small group of elite participants have attended four or more workshops: James Aird, Richard Bower, Ric Davies, Chris Harrison, Dale Kocevski, James Mullaney, David Rosario, Carolin Villforth, and Belinda Wilkes. Five of these, our greatest workshop patrons, also attended the Iceland conference and received a “long-service” award, presented by our special guest, Roger Blandford.

An illustration of these workshops, including the challenging series of natural hazards related to any of our workshops with
<figure><img src="image_41.png" /><figcaption>Figure A.40: The long and interesting road to the first decade of our black-hole growth workshops, illustrated by the workshop posters (top) and punctuated with natural disasters (bottom). The 2014, 2016, and 2022 workshop posters were designed by George Lansbury and inspired by a 1950's boxing promotion poster, a horror B movie, and a disaster movie, respectively. The 2018 workshop poster was designed by Allegra Santis-Rosario, inspired by science fiction and the "Restaurant at the edge of the Universe" book by Douglas Adams. The natural-disaster images show from left to right: the Iceland volcanic-eruption on Eyjafjallajökull, taken on 17 April 2010 (two days before the original start of the 2010 workshop) by Henrik Thorburn, the ash-cloud distribution over the Northern hemisphere constructed by the Metoffice over 14-25 April 2010 with Eyjafjallajökull (red) and Durham (blue) highlighted and labelled, and a Covid-19 coronavirus schematic by CDC on Unsplash to illustrate the global pandemic which lasted from early 2020 to early 2022.</figcaption></figure>

<table><caption>Table A.4: Basic details of our series of black-hole workshops along with the workshop dates, locations, and number of participants.</caption><thead><tr><td>Workshop title</td><td>Dates</td><td>Location</td><td>N</td></tr></thead><tbody><tr><td>What drives the growth of black holes?</td><td>26–29 July 2010</td><td>Durham</td><td>119</td></tr><tr><td>What is the role of AGN in the evolution of galaxies?</td><td>30 July–3 August 2012</td><td>Dartmouth</td><td>64</td></tr><tr><td>The fate of the gas in galaxies: AGN versus star formation</td><td>28 July–1 August 2014</td><td>Durham</td><td>97</td></tr><tr><td>Hidden monsters: obscured AGN and connections to galaxy evolution</td><td>8–12 August 2016</td><td>Dartmouth</td><td>98</td></tr><tr><td>Are AGN special? The environmental dependence and global impact of AGN activity</td><td>30 July–3 August 2018</td><td>Durham</td><td>112</td></tr><tr><td>What drives the growth of black holes? A decade of reflection</td><td>26–30 September 2022</td><td>Reykjavík</td><td>187</td></tr></tbody></table>
the words “What drives” in the title, is shown in Fig. A.40! Our first “What drives the growth of black holes?” workshop was originally scheduled for 19-22 April 2010 at Durham University. However, the Icelandic volcano Eyjafjallajökull, which had first started erupting one month earlier, drove large quantities of fine ash into the atmosphere which got blown across a large fraction of the Northern hemisphere and, critically, all of the UK; see Fig. A.40. This caused European airspace to shut down, forcing us to reschedule the workshop to July, thankfully an eruption-free month. Two determined participants still managed to travel to Durham for the original week in April!

A series of workshops followed on a two-year cycle alternating between Dartmouth College and Durham University. In these meetings we focused on specific subject areas covered in the original “What drives the growth of black holes?” workshop and AH12 review. Early on in the planning we decided that our sixth workshop, 10 years on from the first workshop, would be a celebration of the scientific progress made over the last decade. There could only be one location for this meeting: Iceland, to face Eyjafjallajökull and challenge nature to further disrupt our workshop series! Nature, instead, gave us a 2-year global pandemic with the Covid-19 coronavirus! But, eventually, after failed attempts of holding the meeting in 2020 and then 2021, with some angst and trepidation we successfully scheduled and held the sixth workshop in September 2022 at the Harpa Concert and Conference venue in Reykjavík, Iceland. This meeting was eponymously named after the original workshop and scheduled to have the same scientific sessions but with the additional focus of investigating the scientific progress made over the last decade: “What drives the growth of black holes? A decade of reflection”. The conference photo from this meeting is shown in A.41 and a photo taken following the long-service awards at the conference-dinner venue, the Gamla Bíó theatre, is shown in A.42.

## Appendix B. Citation Network Analysis

Citation network analysis is the application of graph theory and network science to bibliographic data. Citation networks consist of nodes that represent individual papers, with directed edges between nodes representing citations *from* a citing *to* a cited paper. These networks are used to analyse how the scientific literature is structured into specific disciplines, fields, and topics. These structures are self-organising: papers on the same or similar topics cluster together into densely interacting literatures due to the tendency of authors to cite other papers relevant to their own (de Solla Price, 1965). To understand how a literature clusters into different areas, modularity maximisation (Newman and Girvan, 2004) via the Leiden Algorithm(Traag et al., 2019) is a common and well validated approach (Klavans and Boyack, 2017). By these methods, citation networks have been used to produce ‘maps’ of literatures in a number of different areas, including sustainability sciences (Kajikawa et al., 2014), oxytocin research (Leng and Leng, 2021; Leng et al., 2022), and energy storage (Mejía and Kajikawa, 2020)). Beyond its use in literature survey, it is also used to test for cita-

<figure><img src="image_42.png" /><figcaption>Figure A.41: Iceland conference photo taken inside the Harpa concert and conference venue, Reykjavík, in September 2022.</figcaption></figure>

<figure><img src="image_43.png" /><figcaption>Figure A.42: Iceland conference dinner entertainment team on stage at the Gamla Bíó theatre following the long-service awards ceremony. From left to right: Ryan Hickox, Roger Blandford FRS, David Alexander, and Leah Morabito.</figcaption></figure>