[1,2,3]Ming-Yang Zhuang

1]Kavli Institute for Astronomy and Astrophysics, Peking University, Beijing,
100871, China

2]Department of Astronomy, School of Physics, Peking University, Beijing,
100871, China

3]Department of Astronomy, University of Illinois Urbana-Champaign, Urbana,
61801, IL, USA

# Evolutionary Paths of Active Galactic Nuclei and Their Host Galaxies

<EMAIL_ADDRESS>Luis C. Ho  [  [  [

###### Abstract

The tight correlations between the masses of supermassive black holes (BHs)
and the properties of their host galaxies suggest that BHs coevolve with
galaxies. However, what is the link between BH mass ($M_{\rm BH}$) and the
properties of the host galaxies of active galactic nuclei (AGNs) in the nearby
Universe? We measure stellar masses ($M_{*}$), colors, and structural
properties for $\sim 11,500$ $z\leq 0.35$ broad-line AGNs, nearly 40 times
larger than that in any previous work. We find that early-type and late-type
AGNs follow a similar $M_{\rm BH}-M_{*}$ relation. The position of AGNs on the
$M_{\rm BH}-M_{*}$ plane is connected with the properties of star formation
and BH accretion. Our results unveil the evolutionary paths of galaxies on the
$M_{\rm BH}-M_{*}$ plane: objects above the relation tend to evolve more
horizontally with substantial $M_{*}$ growth; objects on the relation move
along the local relation; and objects below the relation migrate more
vertically with substantial $M_{\rm BH}$ growth. These trajectories suggest
that radiative-mode feedback cannot quench the growth of BHs and their host
galaxies for AGNs that lie below the relation, while kinetic-mode feedback
hardly suppress long-term star formation for AGNs situated above the relation.
This work provides important constraints for numerical simulations and offers
a framework for studying the cosmic coevolution of supermassive BHs and their
host galaxies.

##  1 Main

The discovery of tight correlations between the masses of supermassive BHs and
the properties of their host galaxies suggests that BHs and galaxies
coevolve[1]. No consensus has been reached, however, as to how the BH-galaxy
relations arose or how they evolved over time. Popular scenarios ascribe a
critical role for energy feedback from the AGN in regulating the growth of the
BH and its host galaxy[2, 3], but yet some contend that the scaling relations
merely reflect the natural consequences of the hierarchical assembly of BHs
and their host galaxies via merging events[4].

Considerations of when and where the BH-galaxy scaling relations were first
established and how they evolved over time face the immediate practical
difficulty that direct $M_{\rm BH}$ measurements currently are limited only to
the nearest galaxies whose BH sphere-of-influence can be resolved[1, 5]. To
avoid the contamination and disturbances from the active nucleus, direct
$M_{\rm BH}$ measurements also deliberately and nearly universally have
avoided galaxies with substantial nuclear activity. Consequently, the BH-
galaxy scaling relations have been established securely only in the most local
galaxies, ones whose central BH is essentially quiescent and whose overall
star formation activity is also largely muted. This severe limitation can be
bypassed by exploiting indirect methods of estimating $M_{\rm BH}$ using the
properties of the broad-line region, either through reverberation mapping
experiments or empirical estimators of $M_{\rm BH}$ based on single-epoch
spectroscopy calibrated against reverberation-mapped AGNs[6, 7]. Securing
reliable estimates of $M_{\rm BH}$ is only the first step. Obtaining robust
host galaxy parameters with which to compare $M_{\rm BH}$ presents an equally,
if not more, daunting challenge.

Definitive conclusions on the nature of the $M_{\rm BH}$-host scaling
relations of AGNs, including whether and how they evolve with redshift, remain
elusive[8, 9, 10, 11]. The study of cosmic evolution of $M_{\rm BH}$-host
scaling relations using high-$z$ AGNs relies critically on our understanding
of the relation of their low-redshift counterparts. Do nearby AGNs differ with
respect to local inactive galaxies? Do AGN hosts with different morphologies
or bulge types behave the same? Progress has been stymied by the small samples
(typically a few tens, at most a few hundred) and the systematics of host
measurements. To anchor our understanding of the local BH-host galaxy scaling
relations, we urgently need a large sample of nearby AGNs that covers a wide
parameter space and has reliable host measurements, one that can serve as the
definitive benchmark for evolutionary studies.

We derive the structural and photometric properties of the host galaxies of
11,474 $z\leq 0.35$ type 1 AGNs, by performing simultaneous, multi-band
decomposition of Pan-STARRS1 $grizy_{\rm P1}$ images that considers the
wavelength-dependent variation of galaxy structure (see Methods). This sample
allows us to revisit this topic with a more comprehensive, less biased
population of AGNs that spans more than 4 orders of magnitude in AGN
bolometric luminosity ($L_{\rm bol}=10^{41.9}-10^{46.4}\,{\rm erg\,s^{-1}}$),
almost 5 orders of magnitudes in BH mass ($M_{\rm
BH}=10^{5.1}-10^{10.0}\,M_{\odot}$), and a much broader distribution of
Eddington ratios ($\lambda_{\rm E}=10^{-3.5}-10^{0.7}$).

Figure 1: The relation between $M_{\rm BH}$ and $M_{*}$ for $z\leq 0.35$ type
1 AGNs. Objects with late-type, disk-dominated and early-type, bulge-dominated
host morphologies are shown in green and magenta, respectively. Best-fit
relations from G20[5] for local early-type galaxies (black dashed line), late-
type galaxies (black dotted line), and both types combined (black solid line;
intrinsic scatter of 0.81 dex indicated as shaded gray region). We plot the
best-fit relation for the early-type AGNs as a magenta dashed line and for the
late-type AGNs as green dotted line. Data are presented as median values and
typical uncertainty is shown in the lower-right corner. Contours indicate the
distribution of 10%, 30%, 60%, and 90% of the entire sample, respectively.
Dots are individual objects located outside the 90% contour.

The $M_{\rm BH}-M_{*}$ Relation of AGNs in the Nearby Universe

We consider the $M_{\rm BH}-M_{*}$ relations of Greene et al. (2020;
hereinafter G20[5]) in the local Universe as the reference framework. G20 use
dynamical $M_{\rm BH}$ measurements of local, inactive galaxies spanning more
than 6 orders of magnitude in $M_{\rm BH}$ and a wide range of morphological
types. They provide separate relations for early-type and late-type galaxies,
and for both galaxy types combined.

AGNs generally follow the local relations (Fig. 1). To better facilitate
detailed quantitative comparison with G20, we present linear regression fits
to our AGNs. A fit to the early-type, bulge-dominated AGNs using linmix[12]
(the independent and dependent variables are switched when performing the fits
in Equations (1–4), to better trace the orientation of the contours) gives

${\log\,(M_{\rm BH}}/M_{\odot})=(1.88\pm
0.05)\log\,(M_{*}/10^{11}M_{\odot})+(7.89\pm 0.40),$ (1)

with a total scatter of 0.64 dex. For the late-type, disk-dominated AGNs,

${\log\,(M_{\rm BH}}/M_{\odot})=(1.81\pm
0.06)\log\,(M_{*}/10^{11}M_{\odot})+(7.61\pm 0.50),$ (2)

with a total scatter of 0.59 dex. AGNs with late-type hosts follow a relation
close to that of late-type inactive galaxies, with a slightly steeper slope.
Although the BHs in early-type AGNs are on average $\sim 0.6$ dex
systematically more massive than late-type AGNs ($10^{7.8}M_{\odot}$ versus
$10^{7.2}M_{\odot}$), they are still significantly less massive ($-0.7$ dex at
$M_{*}=10^{11}M_{\odot}$) than those in early-type, bulge-dominated inactive
galaxies, which are predominantly quiescent galaxies with very weak star
formation and nuclear activity. Interestingly, the relations of early-type and
late-type AGNs are quite similar. At $M_{*}=10^{11}M_{\odot}$, the difference
in $M_{\rm BH}$ is merely 0.28 dex, smaller than the typical uncertainty of
the single-epoch virial $M_{\rm BH}$ estimates (0.4 dex).

The classical bulges of inactive galaxies, whose formation physics likely
parallels that of ellipticals, follow the same relation as ellipticals, while
pseudo bulges, which share properties akin to those of galactic disks, do not
correlate strongly with $M_{\rm BH}$ and tend to be offset below the relation
of classical bulges[1]. In terms of total stellar mass, late-type galaxies
exhibit $\sim 1.2$ dex lower $M_{\rm BH}$ than early-type galaxies (G20). For
active galaxies, however, definitive conclusions on whether they follow the
relations of inactive galaxies have remained elusive. While some[13] contend
that local AGNs follow a similar correlation between $M_{\rm BH}$ and bulge
stellar mass as inactive galaxies, with no difference between classical and
pseudo bulges, others[14, 15] maintain that nearby AGNs in general lie offset
below the scaling relation of normal classical bulges and ellipticals,
regardless of whether the bulge of the AGN host belongs to the pseudo or
classical variety.

Figure 2: The relation between $M_{\rm BH}$ and $M_{*}$ for $z\leq 0.35$ AGNs.
Colors represent (a) rest-frame optical color $(g^{\prime}-r^{\prime})_{0}$,
(b) Eddington ratio ($\lambda_{\rm E}$), and (c) relative Fe ii strength.
Best-fit relations from G20 for local early-type galaxies (black dashed line),
late-type galaxies (black dotted line), and both types combined (black solid
line; intrinsic scatter of 0.81 dex indicated as shaded gray region). Squares
indicate median values of objects within a box of size 0.2 dex in $M_{*}$ and
0.3 dex in $M_{\rm BH}$ and source number $\geq 5$, after performing locally
weighted regression smoothing using the Python package LOESS[148]. The above
trends persist regardless of whether the data are smoothed (Extended Fig. 2)
or the choice of $M_{\rm BH}$ estimator (Extended Fig. 3). Figure 3: The
relation between $M_{\rm BH}$ and $M_{*}$ for $z\leq 0.35$ AGNs. Colors
represent (a) the median lookback time, (b) growth fraction of $M_{*}$, and
(c) growth fraction of $M_{\rm BH}$. The growth fraction is defined as $\Delta
M/M$ at $z=0$. The host galaxies are assumed to be able to sustain star
formation at the current rate. $M_{\rm BH}$-dependent duty cycle of the BH
accretion is adopted to estimate the mass gain of BHs (see Methods).

Connections with Star Formation and BH Accretion Properties

The large sample size and reliable image decomposition of this study allow us
to investigate, systematically for the first time, the star formation and BH
accretion properties of AGNs distributed across the $M_{\rm BH}-M_{*}$ plane.
In light of the similarity between early-type and late-type AGNs on this plane
(Fig. 1), we do not separate them in the following analysis. First, we confirm
that active galaxies, as in inactive galaxies in the color$-M_{*}$ diagram
(see Extended Fig. 1), exhibit a clear positive correlation between rest-frame
$(g^{\prime}-r^{\prime})_{0}$ color and $M_{*}$, with redder colors found in
more massive hosts (Fig. 2a). We also see a clear correlation between
$\lambda_{\rm E}$ and $M_{\rm BH}$, such that more massive BHs have lower
$\lambda_{\rm E}$ (Fig. 2b), in line with trends established in nearby
galaxies[16] that arise mostly from sample selection bias in magnitude-limited
surveys. Intriguingly, besides the aforementioned expected trends, we find
that $(g^{\prime}-r^{\prime})_{0}$ decreases systematically (color becomes
bluer) toward higher $M_{\rm BH}$ at fixed $M_{*}$, while $\lambda_{\rm E}$
increases significantly toward higher $M_{*}$ at fixed $M_{\rm BH}$. As
$\lambda_{\rm E}$ is derived from $M_{\rm BH}$, it suffers from the large
uncertainty ($\sim 0.4$ dex) associated with single-epoch estimates of $M_{\rm
BH}$ and thus could be driven potentially by unknown systematics. To mitigate
against this concern, we use an independent, empirical indicator of Eddington
ratio, $R_{\rm Fe\,\textsc{ii}}$, defined as the ratio between the equivalent
width of the broad optical Fe ii emission (integrated over 4434–4684 Å) to
that of broad H$\beta$[17]. The surrogate indicator $R_{\rm Fe\,\textsc{ii}}$
also increases toward higher $M_{*}$ at fixed $M_{\rm BH}$ (Fig. 2c),
confirming the robustness of the trends observed based on $\lambda_{\rm E}$
directly. The above trends also persist regardless of whether the data are
smoothed (Extended Fig. 2) or the choice of $M_{\rm BH}$ estimator (Extended
Fig. 3).

Taking into account their relative positions with respect to the local $M_{\rm
BH}-M_{*}$ relation, our results show that objects lying below the local
relation tend to have higher $\lambda_{\rm E}$, while objects lying above the
local relation have characteristically higher specific star formation rate
(sSFR). This supports previous works based on much smaller samples[18, 19, 15,
20]. Within this context, it is noteworthy that the class of AGNs known as
narrow-line Seyfert 1s (NLS1s), which are suspected of experiencing high
accretion rates, also appears to lie offset below the BH-galaxy scaling
relations of inactive galaxies[21], although these claims are
controversial[22].

To understand the cumulative effect of star formation and BH accretion after
the AGN population has evolved to $z=0$, we estimate the amount of mass gain
of the BH and of the stars in the host galaxy (see Methods). Fig. 3 redisplays
the $M_{\rm BH}-M_{*}$ relation, now color-coded by (a) the lookback time and
the relative growth fraction ($\Delta M/M$) till $z=0$ of (b) $M_{*}$ and (c)
$M_{\rm BH}$. The period available for growth ranges from 0.5 to 3.5 Gyr. We
find that the stellar masses of AGNs generally will experience limited growth
(fraction $\lesssim 10\%$), except for objects located on or above the G20
relation of early-type galaxies, which may grow up to $35\%$. Objects lying on
or below the G20 relation of late-type galaxies can grow their BHs
substantially, by as much as $200\%$ for the most extreme outliers.

Combining the growth trajectories of $M_{*}$ and $M_{\rm BH}$, we
schematically outline in Fig. 4 the predicted evolutionary paths for sources
located in different positions across the $M_{\rm BH}-M_{*}$ plane. Three
clear trends emerge: (1) objects above the local relation experience
significant $M_{*}$ growth and migrate mostly horizontally; (2) objects
already on the local relation move along it by growing both $M_{\rm BH}$ and
$M_{*}$ in lockstep for less massive objects, or remain on the relation with
little change in either $M_{\rm BH}$ or $M_{*}$ for massive objects; (3)
objects below the local relation enhance $M_{\rm BH}$ substantially and mainly
ascend vertically toward the relation. Furthermore, objects that depart the
most from the local relation have a greater tendency to evolve toward it. Our
results imply that the location of an AGN on the $M_{\rm BH}-M_{*}$ plane is
closely linked to the star formation history of its host galaxy on the one
hand, and the accretion history of the BH on the other hand.

Figure 4: Evolutionary paths toward $z=0$ of objects on the $M_{\rm BH}-M_{*}$
plane. Evolutionary paths are indicated by the direction and length of each
arrow. The red, green, and blue ellipses represent typical locations of
$z\gtrsim 6$ quasars[44, 11], $z\approx 2$ submillimeter galaxies, SMGs and
nearby luminous infrared galaxies, LIRGs[29, 9, 28], and nearby narrow-line
Seyfert 1s, NLS1s[21, 20]. Best-fit relations from G20 are shown for local
early-type galaxies (black dashed line), late-type galaxies (black dotted
line), and both types combined (black solid line; intrinsic scatter of 0.81
dex indicated as shaded gray region). Figure 5: Comparison of the current
distribution of low-$z$ AGNs (blue) and the expected distribution (red) when
they have evolved to $z=0$. Contours indicate the distribution of 10%, 30%,
60%, 90%, and 95% of the entire sample, respectively. The blue solid line
gives the fit to AGNs of all morphological types, and the red solid line is
the fit to the expected counterparts of the AGNs at $z=0$. Best-fit relations
from G20 for local early-type galaxies (black dashed line), late-type galaxies
(black dotted line), and both types combined (black solid line; intrinsic
scatter of 0.81 dex indicated as shaded gray region).

To visualize the relative positional shift of AGNs toward $z=0$ on the $M_{\rm
BH}-M_{*}$ plane, we compare their distributions in Fig. 5, together with the
fits to the entire sample of AGNs (regardless of morphology) before and after
applying the evolutionary correction. The original distribution (total scatter
0.72 dex)

$\log\,(M_{\rm BH}/M_{\odot})=(2.20\pm
0.04)\log\,(M_{*}/10^{11}M_{\odot})+(7.86\pm 0.30)$ (3)

evolves to $z=0$ to become (total scatter 0.51 dex)

$\log\,(M_{\rm BH}/M_{\odot})=(1.45\pm
0.03)\log\,(M_{*}/10^{11}M_{\odot})+(7.81\pm 0.21).$ (4)

The fit to the original distribution yields a rather steep slope of 2.20,
likely a consequence of mixing hosts of early and late types, which have
slightly different normalizations. After applying the evolutionary correction,
the slope of 1.45 lies between the values of inactive early-type ($1.33\pm
0.12$) and late-type ($1.61\pm 0.24$) galaxies[5].

Intriguingly, the zero point for the AGNs remains offset below that of
inactive galaxies, by 0.5 dex at $M_{*}=10^{11}M_{\odot}$. This difference is
likely contributed by many factors, such as sample size ($\sim 11,500$ in this
work versus $\sim 120$ in G20), potential selection biases, and normalization
and relative fractions of different morphological types. The upper envelope of
the distribution of objects after applying the evolution correction overlaps
with early-type inactive galaxies at $z=0$. The relatively small fraction of
these objects suggests that the BHs in the majority of the local inactive
ellipticals and classical bulges that define the local relation have not been
active since $z=0.35$. As both cosmic star formation history and BH accretion
history peak at $z\approx 2$[1, 23], these massive objects likely experienced
most of their mass growth at high redshifts.

Implications for the Coevolution of Supermassive BHs and their Host Galaxies

Our results reveal clear indications that the level of BH accretion and star
formation of a galaxy is related to its position on the $M_{\rm BH}-M_{*}$
plane. We illustrate these trends graphically by overlaying in Fig. 4 the
typical locations of objects that are found to exhibit various degrees of
departures from the canonical locus of inactive galaxies. We do not consider
the effects of mergers[24], which play a minor role in triggering AGN activity
in the nearby Universe[25].

Why are there undermassive BHs, and how do they grow? Objects with
undermassive BHs on the $M_{\rm BH}-M_{*}$ plane invariably have late-type,
disk-dominated morphologies (Fig. 1). In light of their relatively blue colors
(Extended Fig. 4), it is likely that they contain substantial amounts of cold
gas. Low-redshift AGNs similar to those studied here typically have $\sim
10^{9}\,M_{\odot}$ of molecular gas[26], which maintains healthy levels of
star formation[27]. In this context, a parallel may be drawn with another
class of actively star-forming, massive galaxies notable for having BHs
undermassive with respect to their stellar masses, namely the nearby luminous
infrared galaxies (LIRGs)[28] and their distant counterparts, the
submillimeter galaxies (SMGs) at $z\approx 2-4$[29, 9]. These, too, possess a
significant stellar disk component[30, 31], and they are naturally rich in
dust and gas[32, 33]. The comparison can be further extended to the NLS1s,
which also share some similarities insofar as their BHs appear to be
undermassive[21, 20] and their host galaxy morphologies are conspicuously
late-type[34]. Statistical studies of the cold interstellar medium of low-
redshift quasars, among them members that qualify as NLS1s, detect abundant
cold gas that sustains vigorous ongoing star formation[35, 36]. In the gas-
rich environments that typify these classes of sources, gas inflow leads to
star formation and BH accretion. Although both processes draw from a common
fuel reservoir, the two are not tightly synchronized because they occur on
different spatial and temporal scales. Star formation can take place over a
wide range of spatial scales, so long as the gas is gravitationally unstable,
while BH accretion occurs only in the innermost region of a galaxy after gas
has lost sufficient angular momentum. Star formation proceeds on a local
dynamical timescale, which depending on the location within a galaxy, on the
order of $10^{8}$ yr[37], while AGN lifetimes are $\sim 1-2$ orders of
magnitude shorter[38]. Hence, star formation outpaces BH growth, resulting in
more rapid increase of $M_{*}$ compared to $M_{\rm BH}$. Of course, other
factors can also stunt the growth of the BH[39, 40], such as strong supernova
feedback or insufficient early BH merger events.

Recent studies of the host galaxies of luminous AGNs conclude that, far from
being gas-deficient, they have at least as much molecular gas as normal star-
forming galaxies of the same stellar mass[41]. Closer inspection of the subset
with adequate follow-up observations shows the gas distribution to be
compact[42], and the star formation efficiency to be elevated[35]. The
relatively constant molecular gas fraction and especially the correlation
between molecular gas mass and gas fraction with high $\lambda_{\rm E}$[41]
support the idea that the gas supply to the BH, at fixed $M_{\rm BH}$,
increases with increasing $M_{*}$. This would trigger a higher BH accretion
rate and thus larger mass growth along the horizontal direction toward larger
$M_{*}$, as observed (Figures 2b and 3c). Super-Eddington accretion associated
with the late-stage mergers of LIRGs can grow undermassive BHs by up to an
order of magnitude[28], as found in simulations[3], and restore them to the
local $M_{\rm BH}-M_{*}$ relation of early-type galaxies. The predicted
evolutionary path of systems with undermassive BHs is qualitatively consistent
with recent simulations of BH growth[39, 40], wherein at early times, when
$M_{*}$ is low, the growth of the BH lags that of the stars because of
supernova feedback, while later on the BH catches up its host after the gas
reservoir has become stabilized at higher $M_{*}$. This resonates with the
$M_{*}$-dependent ratio of $\dot{M}_{\rm BH}$ to SFR in X-ray–selected AGNs,
for which $\dot{M}_{\rm BH}$/SFR increases toward higher $M_{*}$[43].

At the opposite extreme, we draw attention to the subset of AGNs that lies
clearly above the $M_{\rm BH}-M_{*}$ relation of inactive early-type galaxies.
Although such objects account for $\lesssim 10\%$ of our sample, they may be
just the tip of an iceberg representing the most extreme of the galaxy
population. These overmassive BHs are reminiscent of $z\gtrsim 6$ quasars[44,
11], which, in view of the limited time ($<1$ Gyr) available for their rapid
formation and growth, pose a challenge for possible formation scenarios[45].
Although the requirement is less stringent for our low-$z$ AGNs, they may
share similar growth channels. Moreover, at least some of these objects may be
the remnants of less extreme high-$z$ quasars.

Studies of local massive elliptical galaxies—the direct descendants of the
most massive, distant quasars—have established that they experienced
substantial size growth compared to their precursor red nuggets[46] at
$z\approx 2$, through the accumulation of an extended outer envelope[47]. In
the popular two-phase formation scenario of early-type galaxies[48], the
transformation of red nuggets to today’s massive ellipticals proceeded largely
through a series of minor (but see [49]), dry mergers. The relatively blue
colors and hence high sSFR of the AGNs with overmassive BHs (Fig. 2a), which
necessarily must gain stellar mass in order to converge to the $M_{\rm
BH}-M_{*}$ relation of inactive early-type galaxies (Fig. 3b) suggest that
their growth evidently was not completely dissipationless and involved some
level of star formation.

The growth trajectories for the majority of the AGNs in our sample follow the
local $M_{\rm BH}-M_{*}$ relation. Applying an evolutionary correction affects
the bulk of the objects rather modestly, with significant changes occurring
only at the outskirts of the distribution (Fig. 5). The mutual dependence of
$L_{\rm bol}$ and SFR on the molecular gas reservoir[50] may account for the
lockstep growth of $M_{\rm BH}$ and $M_{*}$.

Our results are based on nearby AGNs. Abundant gas and frequent mergers at
higher redshifts should lead to rapid evolution with more pronounced growth on
the $M_{\rm BH}-M_{*}$ plane. For the obscured type 2 AGNs, it is difficult to
directly measure BH mass for a large sample to investigate their $M_{\rm
BH}-M_{*}$ relation and whether they behave similarly to the type 1 AGNs.
However, based on 14 Seyfert 2 AGNs with near-infrared broad-line detections
(see Methods), we find that they are located on the local relation and follow
the same trends as their unobscured counterparts (Extended Fig. 5). This
similarity suggests a common evolutionary scenario for both type 1 and type 2
AGNs.

Implications for AGN Feedback

AGNs with blue host galaxies have relative slowly accreting BHs, whereas AGNs
with rapidly accreting BHs tend to be hosted by redder galaxies (Fig. 2). Are
we witnessing the effects of negative AGN feedback? Radiative-mode feedback
has been an appealing mechanism proposed for quenching star formation in
galaxies[2], although linking it to observable consequences has proven not
straight-forward[51]. If negative AGN feedback operates effectively for the
general population of active galaxies as viewed through their position on the
$M_{\rm BH}-M_{*}$ plane, we expect $(g^{\prime}-r^{\prime})_{0}$ to correlate
with $\lambda_{\rm E}$: more highly accreting BHs, having more powerful AGNs,
inhibit star formation and turn the hosts red. This basic expectation is not
borne out by our sample. Indeed, the exact opposite is seen (Extended Fig. 6).
This disfavors the proposition that strong, instantaneous, negative feedback
from accreting BHs quenches star formation, at least in nearby galaxies, in
concert with the conclusions of previous studies of AGN hosts based on
measurements of gas content, SFR, and star formation efficiency[35, 26, 41,
36, 27].

What about the AGNs with undermassive BHs with high $\lambda_{\rm E}$ and red
$(g^{\prime}-r^{\prime})_{0}$ colors (Fig. 2)? At first glance, these traits
appear to implicate negative AGN feedback. We argue, however, that AGN
feedback only minimally affects the interstellar medium in the host galaxies
of this population of objects. The consequences of AGN winds on the host
galaxy depend on the properties of the shocks that they impart on the
interstellar medium[52]. If the wind cools efficiently on a timescale shorter
than its flow time, most of its kinetic energy is lost. As momentum is
conserved, a momentum-driven outflow injects only $\sim 10\%$ of the binding
energy of the bulge to the interstellar medium, which is hardly capable of
clearing the surrounding gas and preventing accretion. It has been
demonstrated that an outflow is momentum-driven when $M_{\rm BH}$ is smaller
than a critical mass[52], which scales with the fourth power of the velocity
dispersion (or, equivalently, mass) of the bulge, very close to the observed
local relation[53]. If we assume a bulge-to-total ratio of 0.2, as typical of
late-type AGN host galaxies[20] and nearby inactive disk galaxies[54], most
objects with a predicted increase of $M_{\rm BH}$ larger than $20\%$ still lie
below the $M_{\rm BH}-M_{\rm bulge}$ relation of early-type galaxies. This
implies that $M_{\rm BH}$ is smaller than the critical mass and that the
outflows are momentum-driven, rendering them incapable of shutting down BH
growth or curtailing star formation (Extended Fig. 7). This is consistent with
the results from a recent detailed analysis of the properties of the X-ray
outflows in I Zwicky 1, a well-studied NLS1 that is located far below the
local relation[55].

When $M_{\rm BH}$ reaches the critical mass, as is the case for AGNs with
overmassive BHs above the early-type relation, the outflow becomes energy-
driven, and the adiabatically expanding bubble from kinetic feedback
efficiently expels the inflowing gas and suppresses further BH growth.
However, our analysis suggests that $M_{*}$ increases significantly even after
substantial gain in $M_{\rm BH}$ (Fig. 3b). This implies that although AGN
feedback can suppress the growth of the BH, it fails to prevent long-term star
formation on galactic scales, seemingly at odds with the conclusions of recent
cosmological simulations[3, 56].

##  2 Methods

###  2.1 Low-$z$ Type 1 AGN Sample

We analyze 14,574 type 1 AGNs with $z\leq 0.35$ selected by Liu et al.
(2019)[57] from the seventh data release (DR7)[58] of the Sloan Digital Sky
Survey (SDSS)[59] to have broad H$\alpha$ emission. The redshift cut ensures
that H$\alpha$ lies within the SDSS spectral coverage. This constitutes the
largest low-$z$ type 1 AGN sample studied in terms of their host galaxy
properties. It surpasses other similar samples by a factor of $\sim 3$ in
size, and it covers a wider range of AGN luminosities owing to the selection
based on the detection of broad H$\alpha$ instead of H$\beta$, as well as the
lower threshold on line width used to define broad emission lines[60, 61, 62,
63]. The emission lines are measured after subtracting a model for the
continuum, which includes a broken power law for the AGN accretion disk, a
Balmer continuum, if present, an optical Fe ii pseudo-continuum, and the host
galaxy starlight represented by stellar population templates[57]. Liu et al.’s
catalog provides basic properties such as the AGN continuum luminosity at 5100
Å ($L_{5100}$), the fluxes, line widths, and equivalent widths of strong
permitted and forbidden lines, as well the strength of the broad Fe ii
emission blend at 4434–4684 Å. We estimate the AGN bolometric luminosity
($L_{\rm bol}$) from $L_{5100}$ with a bolometric correction of $L_{\rm
bol}=9.8L_{5100}$[64].

We obtain BH masses ($M_{\rm BH}$) from the full width at half-maximum (FWHM)
of broad H$\beta$ and $L_{5100}$, using the virial mass estimator based on
single-epoch spectra[65]. This estimator takes into account different
calibration coefficients for pseudo and classical bulges, which have virial
factors of $f=3.2$ and 6.3, respectively[66]. We adopt the calibration for
pseudo bulges for hosts with reliable late-type, disk-dominated morphology,
and we use the zero point for classical bulges for hosts that are robustly
classified as early-type, bulge-dominated systems. Hosts without reliable
morphologies get assigned the zero point for both bulge types combined. When
broad H$\beta$ is absent, we convert the FWHM of broad H$\alpha$ to that of
broad H$\beta$ using the empirical relation between the two line widths[67].
The derived BH masses are on average $0.17\pm 0.26$ dex lower for late-type
hosts and $0.23\pm 0.28$ dex higher for early-type hosts compared to those
derived from broad H$\alpha$ alone with $f=4.47$[68]. Switching to the
H$\alpha$-based $M_{\rm BH}$ estimator does not qualitatively affect our main
conclusions. Extended Fig. 3 shows a different version of Fig. 2 without
performing smoothing and with $M_{\rm BH}$ estimated from broad H$\alpha$.

Recent reverberation mapping observations of nearby type 1 AGNs suggest that
the scatter in the relation between broad-line region radius and luminosity
depends systematically on Eddington ratio[69, 70, 71] (but see counter
argument[72]), $\lambda_{\rm E}\equiv L_{\rm bol}/L_{\rm E}$, where $L_{\rm
E}=1.26\times 10^{38}\,(M_{\rm BH}/M_{\odot})$ is the Eddington luminosity. If
this dependence on Eddington ratio is not taken into account, the virial BH
masses estimated from single-epoch spectra, which rely on the radius-
luminosity relation, will be underestimated for AGNs with low $\lambda_{\rm
E}$ and overestimated for those with high $\lambda_{\rm E}$. This effect might
impact our results, in view of the variation of Eddington ratio with
morphological type observed in our sample. To mitigate this potential
systematic effect, we conservatively assign an uncertainty of 0.4 dex to the
$M_{\rm BH}$ estimates.

The final sample after excluding objects with unreliable measurements on host
properties ($\sim 78\%$ of the original) consists of 11,474 moderately
luminous to powerful AGNs, which are characterized by $L_{\rm
bol}=10^{41.9}-10^{46.4}$ erg s-1 (median $10^{44.6}$ erg s-1), $M_{\rm
BH}=10^{5.1}-10^{10.0}\,M_{\odot}$ (median $10^{7.6}M_{\odot}$), and
$\lambda_{\rm E}=10^{-3.5}-10^{0.7}$ (median $10^{-1.2}$). It still well
represents the parent sample (Supplementary Fig. 1), with a major difference
being the significant reduction (by $\sim 50\%$) of the number of objects with
$L_{\rm bol}\gtrsim 10^{45.4}$ erg s-1. This sample is nearly 40 times larger
than that in any previous work in the nearby ($z\leq 0.35$) Universe[73, 74,
75, 76, 14, 19, 77, 15, 78, 13, 79, 80, 81].

###  2.2 Low-$z$ Inactive Sample

To investigate systematically the AGN host galaxies and the potential biases
that can be introduced from the presence of a bright nucleus, it is vital to
select a sample of inactive galaxies in the same redshift range that covers a
similar parameter space in host galaxy properties as the AGN sample. We start
from the SDSS data release 16 spectroscopic catalog (specObj table;
https://data.sdss.org/sas/dr16/sdss/spectro/redux/specObj-dr16.fits) and limit
to galaxies within a radius of 2.5′ of each AGN. This is the minimum size of
the image cutout of the AGNs, one that yields similar image properties, in
terms of depth and point-spread function (PSF), for both the galaxy and AGN
sample. We exclude duplicate observations of the same objects and
spectroscopically identified AGNs by omitting objects with CLASS=QSO,
SUBCLASS=AGN or AGN BROADLINE. The above SUBCLASS criteria remove objects
classified as (broad-line) Seyferts and low-ionization nuclear emission-line
regions (LINERs)[82, 16] based on the [O iii]/H$\beta$ versus [N ii]/H$\alpha$
line-intensity ratio diagnostic diagram[83], to leave a clean sample of
inactive galaxies free of potential effects of the central AGN on the galaxy
properties. After selecting objects with $0.01\leq z\leq 0.35$ and ZWARNING=0,
we obtain 12,236 inactive galaxies, among which 11,576 are retained after
excluding 538 objects without image coverage (located in a different image
tile) and 122 objects affected by image masks. As we are mainly concerned
about the reliability of photometric and structural properties of the AGN host
galaxies, the final inactive galaxy sample consists of 11,017 objects
($95.2\%$) with robust (measurement/error $>3$) magnitudes in all five bands,
half-light radii, and Sérsic indices, as well as successful spectral analysis
from spectral energy distribution (SED) fitting. The inactive sample, apart
from the absence of an active nucleus, ideally should have exactly the same
distribution of galaxy properties as the sample of active galaxies. This is
difficult to realize in practice, but the above-constructed inactive galaxy
sample suffices for our primary goal of probing the typical photometric
parameter space spanned by the AGN host galaxy sample.

###  2.3 Reduction of the Pan-STARRS1 Data

We use the coadded five-band ($grizy_{\rm P1}$) images (pixel size 0.25′′)
from Pan-STARRS1 (PS1)[84]. We retrieve the stack images, weight images, and
mask images, as well as the source catalog that includes coordinates, PSF
magnitudes ($m_{\rm PSF}$), and Kron[85] magnitudes ($m_{\rm Kron}$) from the
second data release (DR2)[86] of PS1. Each image cutout is centered on the
position of the AGN with an initial size of $5^{\prime}\times 5^{\prime}$, and
the size increases in steps of 1′ according to the number of stars available
for PSF construction. We exclude 92 objects that lack an image in at least one
band, that are affected by close bright stars and foreground galaxies, or that
have at least 50% of the pixels masked in the central 1 arcsec2. For the
remaining 14,482 AGNs, we subtract a two-dimensional background model and
generate segmentation maps and source catalogs using SExtractor[87]. It has
been demonstrated that PSF models generated by stacking nearby stars perform
best in describing the light profiles and recovering the magnitudes of point-
like sources by comparing the behavior of different codes for PSF
construction, finding no significant spatial variation within the $\sim$30
arcmin2 field-of-view of the stacked PS1 images[27]. Following their
procedure, we use SWarp[88] to construct PSF models by stacking isolated
point-like sources in the image cutout, which are selected by requiring
$m_{\rm PSF}<21$ mag and $m_{\rm PSF}-m_{\rm Kron}<0.05$ mag in the $i_{\rm
P1}$ band[89]. The number of applicable stars is determined by gradually
increasing the size of the image cutout so that the magnitudes of the PSF
models are brighter than $m_{\rm PSF}$ of the AGN from the PS1 catalog
(conservative overestimation) in all five bands, to ensure enough signal-to-
noise ratio to model the light profile of the nucleus.

We summarize the image properties of our sample in Supplementary Fig. 2 and
Supplementary Table 1, which gives the statistics of the FWHM of the PSF
obtained from fitting a circular, two-dimensional Moffat[90] function, the
$3\sigma$ surface brightness depth, and the $5\sigma$ point source depth,
where $\sigma$ is the standard deviation of 100 flux measurements made using a
circular aperture of diameter 15′′ for extended sources and 2 times the PSF
FWHM for point sources, randomly placed in the source-masked background of the
image.

###  2.4 Two-dimensional Simultaneous Multiwavelength Fits

The morphology and structure of a galaxy can change significantly with
wavelength as a result of variations in stellar population, metallicity
gradient, or dust attenuation. Using GALFITM[91, 92], it has been demonstrated
that both the Sérsic[93] index ($n$) and half-light radius ($R_{e}$) of a
galaxy vary smoothly from the ultraviolet to the near-infrared, with
increasing $n$ and decreasing $R_{e}$ toward longer wavelengths[94]. GALFITM
is a multiwavelength version of GALFIT[95, 96], which includes all the
functions of the original code and models the wavelength dependence of
parameters using a series of Chebyshev polynomials. By simultaneously fitting
multiple images at different wavelengths, GALFITM can improve the accuracy,
reliability, and signal-to-noise ratio of the derived parameters[91, 97].

We use GALFITM to model the unresolved nucleus with the PSF and the underlying
host galaxy with a single Sérsic component. The choice of wavelength
dependence of model parameters follows that of previous work[27]. We keep the
ellipticity and position angle constant across all bands, but allow the
magnitudes of the PSF and Sérsic components to vary with wavelength as a
fourth-order polynomial, while $n$ and $R_{e}$ change linearly with
wavelength. The positions of the PSF and Sérsic components are tied together,
under the assumption that the supermassive BH is located at the center of the
galaxy.

We also perform a set of decompositions with the positions of the PSF and
Sérsic components untied, allowing their central coordinates to vary by up to
10 pixels. Consistent results are obtained in the $i_{\rm P1}$ band for
$m_{\rm S\acute{e}rsic}$, $R_{e}$, and $n$, with no systematic differences in
the median and a scatter of $\lesssim 0.06$ mag, 0.03 dex, and 0.07 dex,
respectively. We caution that leaving too much freedom to the central
positions of the two components may lead to artificial offsets in hosts with
fainter surface brightness[98]. The segmentation map and source properties
returned by SExtractor are used to generate the source mask and to set up the
initial guesses for the model parameters and fitting region. We use the weight
maps from PS1 DR2, which properly account for the uncertainty of the image
using information from the image stacking and data reduction. All the sources
inside the fitting region are masked, except for the AGN and possible blended
companions identified by SExtractor with DEBLEND_NTHRESH=32 and
DEBLEND_MINCONT=0.005, which, if present, are fitted simultaneously. The
fitting region is set to be the larger value between 50 pixels and $10$ times
A_IMAGE, which is the length of the semi-major axis of the isophotal
elliptical shape of the source. If a companion is present, the source extent
and the separation between the AGN and the companion are summed up to
determine the final fitting region. We model the companion with a PSF function
if its CLASS_STAR > 0.9, and as a Sérsic function otherwise. We use $m_{\rm
PSF}$ from the PS1 catalog and $m_{\rm S\acute{e}rsic}$
$=-2.5\,\log\,(10^{-m_{\rm Kron}/2.5}-10^{-m_{\rm PSF}/2.5})$ as initial
guesses of the magnitudes for the PSF and Sérsic components, respectively. For
the structural parameters of the Sérsic model, we take values from the output
of SExtractor for $R_{e}$, axis ratio, and position angle, but we restrict the
initial guess of $n$ to 2. We constrain $R_{e}$ to vary within $0.5-500$
pixels and $n$ within values of $0.3-8$.

The simultaneous, multi-object, multiwavelength decomposition is illustrated
in Supplementary Fig. 3 for a late-type host galaxy, and in Supplementary Fig.
4 for an early-type system. In Supplementary Fig. 3, the prominent emission
from the nucleus is well extracted and clearly visible in the surface
brightness radial profiles, despite its small contribution (as low as 3%) to
the total flux. The hole in the central region, evident in the data $-$ model
and residuals images, indicates that the nucleus is oversubtracted. As
discussed in the following subsection, this is more severe in objects with
lower AGN-to-host contrast but does not significantly affect the measured
galaxy properties. Although substructures from the bulge, bar, and spiral arms
can leave noticeable residuals, a single-component Sérsic model can recover
the total magnitude of the galaxy as robustly as more detailed (e.g.,
bulge+disk) models[99, 97].

Successful fits were obtained for a total of 12,903 objects. We cannot derive
reliable magnitudes (magnitude error $>0.36$ mag, equivalent to flux/flux
error $<3$) for 1,579 sources, which are mainly powerful ($L_{\rm bol}\approx
10^{45.0}$ erg s-1) AGNs at higher redshift ($z\approx 0.3$).

###  2.5 Verification of Results

Although a single-component Sérsic model can deliver a reliable total
magnitude for the galaxy, the situation is more complicated with the presence
of a prominent active nucleus. The structure and color of a galaxy can
significantly affect the magnitude of the nucleus, and the bright nucleus can
also bias the measurements of the galaxy parameters. Moreover, the formal
errors returned by fitting codes always underestimate the true uncertainties,
regardless of the goodness of the fit[100, 101]. We perform the most
comprehensive input-output tests to date to derive realistic uncertainties and
study the potential systematics introduced by the presence of a bright
nucleus. We produce two sets of mock AGN images: one from model galaxies
generated from the best-fit parameters, and the other from actual observations
of low-$z$ inactive galaxies with structural and photometric properties
carefully matched to each AGN. The second approach is often employed in the
literature[102, 79], but only using control samples matched in redshift and
stellar mass, neglecting the effect of color and morphology of the AGN host
galaxies.

####  2.5.1 Mock AGNs from Idealized Galaxies

We quantify the uncertainty of each model parameter by constructing idealized
mock AGN images using the best-fit parameters of each object while mimicking
the conditions of its imaging properties. Three sources of noise are applied
independently to the data: (1) Poisson noise associated with the source, (2)
Gaussian noise of the individual background pixels, and (3) Gaussian noise of
the constant background after background subtraction prior to the
decomposition. The last source of noise is related to the extent of the source
and is measured as the standard deviation of 100 ellipses randomly distributed
in the source-masked, background-dominated regions (the background constitutes
more than $95\%$ of all pixels), with semi-major axis 3 times $R_{e}$ and axis
ratio and position angle fixed to the best-fit values of the AGN host; these
apertures enclose $96\%$ and $79\%$ of the total light for Sérsic profiles
with $n=1$ and $n=4$, respectively. We generate 100 realizations of mock AGNs
for each object in each band, and we repeat the decomposition using GALFITM
with the same configuration file and constraints as for the real observations.
The final measurement and its uncertainty are taken as the median value and
standard deviation of the 100 results. We exclude 580/12,903 ($4.5\%$) objects
for the following reasons: (1) they have unreliable magnitudes for both the
PSF and Sérsic components; (2) the fitting results for the magnitudes,
$R_{e}$, and $n$ of the observed and mock AGNs differ by more than 3 times the
uncertainty of the mock AGNs; and (3) the magnitude of the Sérsic component is
unreliable in at least two bands for either the real observations or the mock
AGNs, which would render the SED analysis highly uncertain. The uncertainties
derived from the idealized model AGNs are added in quadrature to the final
error budget of each parameter. The same procedure is adopted also for the
low-$z$ inactive galaxy sample, which resulted in 11,479 objects.

####  2.5.2 Mock AGNs from Real Galaxies

Model mismatch due to the presence of galactic substructures, such as the
bulge, bar, ring, lens, or spiral arms, even if they are not well resolved at
relatively higher redshifts, can lead to significant systematic parameter bias
from our adoption of a simplified two-component (PSF + Sérsic) model. It is of
fundamental importance to understand how the final parameters of each AGN host
galaxy may be affected by this possible source of systematic bias. To
accomplish this goal, we generate mock AGNs by manually adding an artificial
point source to the center of each real galaxy selected from the low-$z$
inactive sample. We then study the derived parameters before and after adding
the point source. Although we do not consider explicitly galaxy substructures,
we note that $n$ does correlate with galaxy concentration and bulge-to-total
ratio[103].

Galaxies that are similar to each AGN are selected based on the error-weighted
difference (${\rm\mid AGN-gal\mid}/\sqrt{\rm error(AGN)^{2}+error(gal)^{2}}$)
of four parameters: observed $(g-i)_{\rm P1}$ color, $R_{e}$, $n$, and axis
ratio. The error of axis ratio is not used but instead the largest difference
between the galaxy and the AGN cannot exceed 0.1. If the error is too small,
we adopt a fiducial value of 0.1 mag for color and 10% for $R_{e}$ and $n$. An
integer factor of $2-5$ is multiplied to the error to increase the number of
inactive galaxies if it is fewer than 20. We exclude 33 AGNs that do not have
enough inactive galaxies even after applying an error boost factor of 5. If
more than 20 objects satisfy our criteria, we choose the top 20 with the
smallest magnitude difference between the AGN and the inactive galaxy in the
$i_{\rm P1}$ band. Supplementary Fig. 5 illustrates our procedure for the
target AGN SDSS J011448.68$-$002946.1 (already shown in Supplementary Fig. 3),
in which we compare its $i_{\rm P1}$-band image to those of a collection of 20
inactive galaxies that have similar structural properties and color. In total,
we use 245,800 inactive galaxies (each galaxy can be used multiple times) to
generate mock images for 12,290 AGNs. The error boost factor is $\leq 3$ for
$\sim 95\%$ of the sources. The distributions of differences between AGNs and
inactive galaxies in observed-frame $(g-i)_{\rm P1}$ color, $R_{e}$, and $n$
confirm the similarity between the real AGNs and their mocks (Supplementary
Fig. 6). The median differences cluster near 0 with a small dispersion; the
tails of the distribution arise from a tiny fraction of sources with
relatively large parameter errors.

To generate the corresponding mock active galaxies for a given AGN target, we
add a point source to the center of each of the 20 inactive galaxies. The
magnitude of the point source in each band is calculated from $L_{\rm bol}$
and the composite quasar SED[104], independent of the robustness and
reliability of $m_{\rm PSF}$ from the image decomposition. To ensure a
consistent AGN-to-host contrast in each band, these magnitudes are further
subject to corrections according to the magnitude difference between the host
galaxy of the AGN target and each inactive galaxy.

We perform simultaneous multiwavelength decomposition of the mock AGNs
following the same procedure as employed for the real AGNs. Supplementary Fig.
7 shows the general trend of the difference ($\rm\Delta=AGN-galaxy$) between
the results measured in the mock AGNs and those measured in the inactive
galaxies used to generate them, as a function of $m_{\rm S\acute{e}rsic}$,
$R_{e}$, $n$, and $m_{\rm AGN}-$$m_{\rm S\acute{e}rsic}$, which indicates the
AGN-to-host contrast. We find that the $R_{e}$ of AGN host galaxies can be
measured reliably with no systematic bias ($\Delta\log R_{e}\approx 0.01$ dex;
scatter $\sim 0.03$ dex in all five bands) and with little dependence on the
properties of the galaxy or relative AGN brightness. By contrast, $m_{\rm
S\acute{e}rsic}$ and $n$ exhibit small systematic biases ($\Delta$$m_{\rm
S\acute{e}rsic}$$=0.07$ mag; $\Delta\log\,n=-0.11$ dex) and slightly larger
scatter compared to $R_{e}$ (90% of the sources are located within a
dispersion of $\sim 0.25$ mag for $\Delta$$m_{\rm S\acute{e}rsic}$ and $\sim
0.2$ dex for $\Delta\log\,n$). The scatter of $\Delta$$m_{\rm S\acute{e}rsic}$
and $\Delta\log\,n$ are affected by $m_{\rm S\acute{e}rsic}$ and the AGN-to-
host contrast, with larger scatter toward fainter host magnitude and
relatively smaller scatter in redder, more star-dominated bands. We also find
a tendency for $n$ to be underestimated in galaxies with larger $R_{e}$. Among
the parameters investigated here, $m_{\rm AGN}$ performs least favorably. The
AGN brightness is systematically overestimated when the nucleus is much
fainter than its host galaxy, in the most extreme cases by as much as $\sim 3$
mag when the AGN is only 1% as bright as the host. We apply these object-
tailored corrections to every AGN and combine them with the results from
previous subsection to derive the final value and uncertainty of each
parameter.

###  2.6 Comparison to Decomposition Results for Two-component Sérsic Models

To investigate in detail the effect on host magnitude when neglecting
substructures, we select all late-type ($n<2$) objects with sufficient
resolution (PSF FWHM $<2$ kpc in all five bands; resolution better than 1 kpc
in radius) to detect potential central substructures, such as a bulge. This
gives us 581 AGNs, $\sim 5\%$ of the parent sample. We perform a new set of
image decomposition for these objects, following exactly the same procedure
used for the parent sample but using, instead, a two-component Sérsic model,
where for the first component the Sérsic index set to $n=1$ to mimic emission
from a disk. The position of both components are tied together, and $R_{e}$
for the second component is constrained to be smaller than that of the first.

We find that in the $i_{\rm P1}$ band the total host magnitude (sum of the
fluxes of both Sérsic components) from the double-component fit is consistent
with that from the single-component fit, with a median difference of $-0.05$
mag ($5\%$; fluxes from the double-Sérsic fit are brighter) and a scatter of
0.16 mag ($15\%$). This indicates that the single-component can recover the
total magnitude of the galaxy reliably. Moreover, the median difference is in
excellent agreement with the median correction of $-0.06$ mag applied to
$m_{\rm S\acute{e}rsic}$ for the single-component decomposition, as estimated
from our mock AGNs based on real inactive galaxies. Nevertheless, systematic
biases may still be present even for the two-Sérsic models. It is much harder
to follow our procedure of generating mock AGNs with similar morphologies to
study systematic biases, because of the lack of a sufficient amount of real
inactive galaxies with similar bulge and disk parameters. Assuming negligible
systematic bias in host magnitude from the double-Sérsic component
decomposition, we estimate the stellar masses of host galaxies following the
procedures and configurations to be introduced below. We find consistent
results: the median difference in stellar mass is $\sim 0.05\pm 0.09$ dex,
which is smaller than the uncertainty of individual objects ($\sim 0.10$ dex).

We conclude that adding an extra host component has only a minor effect on the
magnitudes and stellar masses of a small fraction of AGN host galaxies in our
sample. For simplicity, throughout this work we adopt the results ($m_{\rm
S\acute{e}rsic}$, $R_{e}$, and $n$) from the single-component Sérsic models.

###  2.7 Stellar Masses of AGN Host Galaxies

We derive the stellar masses of the host galaxies from the AGN-removed
photometry, using version 2022.0 of the SED-fitting code CIGALE[105, 106]. The
flux densities of the five PS1 bands are first corrected for Galactic
foreground extinction using the Python package dustmaps[107], using a Milky
Way dust map[108] and converting $E(B-V)$ to extinction in the PS1 bands using
$R_{V}=3.1$[109]. An additional $10\%$ uncertainty is added in quadrature to
the flux densities to account for the uncertainty in absolute flux calibration
and other potential unknown factors. In light of the limited available
photometric coverage of our SEDs, we adopt a parametric “delayed” star
formation history combined with single stellar population models[110], a
nebular emission-line model[111], and a dust attenuation law (dustatt modified
starburst) based on the starburst attenuation curve[112]. Albeit simplistic,
parametric star formation histories have been demonstrated to yield robust
stellar masses with $\lesssim 10\%$ systematic bias and $<20\%$ scatter, even
when only optical photometry is available[113, 114]. The models and fitting
parameters for CIGALE are given in Supplementary Table 2. We scale the stellar
masses by a factor of 1.08[23] to convert the initial mass function from
Chabrier[115] assumed by CIGALE to that of Kroupa[116] adopted in this paper.
We also test the reliability of the stellar masses directly derived from the
observed integrated (AGN$+$host) SED, using SED-fitting as described in the
following section. In general, the stellar masses derived from integrating the
SEDs can be overestimated by as much as $\sim 0.5$ dex for AGNs of the
strength typically found in our sample.

To select the most reliable SED measurements, we use the difference between
the estimate of $M_{*}$ from the best-fitting template ($M_{*,\rm{best}}$) and
the mean of the probability density function distribution of the likelihood of
all templates ($M_{*,\rm{bayes}}$) and $\chi^{2}_{\nu}$. We apply a strict cut
to keep only objects with the most secure SEDs, by excluding 816 ($6.6\%$)
sources with $\mid\log\frac{M_{*,\rm{bayes}}}{M_{*,\rm{best}}}\mid>0.3$ and
$\chi^{2}_{\nu}>1$. For the remaining 11,474 AGNs, we obtain the rest-frame
SDSS $g^{\prime}-r^{\prime}$ color, $(g^{\prime}-r^{\prime})_{0}$, and adopt
$M_{*,\rm{bayes}}$ and its associated uncertainty as the final estimate of
stellar mass. The median uncertainty of $M_{*}$ is $\sim 0.10$ dex ($23\%$),
which is larger than the uncertainties of the input flux densities ($\sim
11\%$). This suggests that the uncertainty of $M_{*}$ is dominated by the
degeneracy of the templates.

G20[5] assume a “diet Salpeter” initial mass function by using mass-to-light
ratios from Ref.[117] to derive $M_{*}$. To match the initial mass function
adopted in this work, we scale down their stellar masses by a factor of 1.06.

###  2.8 Stellar masses of AGNs derived from AGN$+$host SEDs

It is common to derive the stellar masses of AGN host galaxies by fitting
their total (AGN$+$host), integrated SEDs. To test the reliability of this
method, we use the stellar masses and SEDs of the inactive galaxy sample to
generate mock SEDs for the AGNs by adding a nuclear flux to each band,
following the same procedure used to generate the mock AGN images. We consider
three typical AGN strengths ($L_{\rm bol}=10^{43},10^{44}$, and $10^{45}\,\rm
erg~{}s^{-1}$) and limit the inactive sample to objects with
$M_{*}=10^{9.5}-10^{12}\,M_{\odot}$ to cover a parameter space similar to that
of our actual AGN sample. We then fit their SEDs using exactly the same
configuration of CIGALE models, but with an additional AGN component generated
using SKIRTOR[118, 119]. We set the inclination angle to four discrete values
($i=0^{\circ}$, $10^{\circ}$, $20^{\circ}$, $30^{\circ}$) and allow the AGN
fraction (defined in the wavelength range 0.3–1.0 $\mu$m) to freely vary
($f_{\rm AGN}=0.001-0.9$), while fixing all other parameters to their default
values, since they are related to the properties of the AGN torus and are not
relevant here.

Supplementary Fig. 8 compares the stellar masses derived from the mock AGNs
and those measured in the inactive galaxies used to generate them (the ground
truth). We find that the differences between the two sets of measurements and
their scatter increase with increasing AGN $L_{\rm bol}$ (first column),
likely a consequence of the increasing degeneracy between the AGN and host
galaxy toward larger AGN-to-host contrast and bluer host galaxy color. As a
result, CIGALE tends to overestimate $f_{\rm AGN}$ toward bluer host galaxies
(vertical trend). We also find a turnover (“V”-shape) of the trend in
relatively red host galaxies: $\Delta M_{*}$ first decreases and then
increases with increasing AGN-to-host contrast. This can be understood as
underestimation of $f_{\rm AGN}$ at low AGN-to-host contrast and
overestimation of $f_{\rm AGN}$ toward high AGN-to-host contrast. As mass-to-
light ratio rises for galaxies with redder colors[117], stellar mass would be
overestimated if $f_{\rm AGN}$ is overestimated, and vice versa.

We also fit the AGN$+$host SEDs of the AGN sample using measurements of
$m_{\rm Kron}$ from the PS1 DR2 catalog. Supplementary Fig. 9 compares stellar
masses from fitting the AGN+host SEDs (SED) and those from SED-fitting to AGN
contribution-subtracted host SEDs from image decomposition (IMG). Consistent
with results of mock AGNs, stellar mass difference ($\Delta=$ SED$-$IMG) tends
to be more positive toward smaller $M_{*}$, higher AGN-to-host contrast and
bluer host color (Supplementary Figures 9a and 9b), but with larger scatter.
The larger scatter is partly due to the difference in fluxes between model
fitting and aperture photometry, as evident from a negative trend of $\Delta
M_{*}$ with $n$. At larger $n$, a higher fraction of light is extended, with
fluxes enclosed within an aperture of radius of 3 times $R_{e}$ contribute to
96% and 79% of the total fluxes for Sérsic profiles with $n=1$ and $n=4$,
respectively. This potential missing flux issue, together with the different
aperture sizes for different bands, renders fluxes measured from aperture
photometry less reliable and less self-consistent. Supplementary Fig. 9d
indicates that after controlling for the effect of $n$, the scatter of $\Delta
M_{*}$ is relatively small at low AGN-to-host contrast (consistent to the
middle panel in Supplementary Fig. 8a), while the scatter increases as a
combined effect of higher AGN-to-host contrast and large host color range
(consistent to the middle panels in Supplementary Figures 8b and 8c). These
results demonstrate the biases and uncertainties in using SED-fitting to
derive stellar mass of type 1 AGNs using AGN$+$host SEDs.

###  2.9 Color$-$Stellar Mass Diagram and Estimation of Specific Star
Formation Rate

Broad-band optical colors encode valuable information about galaxies. For
example, the SDSS $g^{\prime}-r^{\prime}$ color is an indicator of star
formation properties and is widely used to classify galaxies into red
(quiescent) and blue (star-forming) sequences[120, 121]. Using the rest-frame
SDSS $g^{\prime}-r^{\prime}$ color from SED fitting, we compare the
distributions of active and inactive galaxies on the color-stellar mass
diagram. Extended Fig. 1 splits the samples into two redshift bins, one local
($z<0.15$) and the other slightly more distant ($0.15\leq z\leq 0.35$). Note
that the contours of the inactive galaxies are only used to establish the
typical range of parameter space instead of giving a quantitative comparison
with their active counterparts, since other properties (e.g., morphology or
Sérsic index) correlate with color[122], and inactive galaxies are biased
toward luminous red galaxies at higher redshift. Nevertheless, we find that
the AGNs in the two redshift bins span a wide range of color, residing in both
blue and red galaxies. AGNs are hosted mostly in galaxies with stellar masses
$M_{*}\gtrsim 10^{10}\,M_{\odot}$, which mainly reflects the selection effect
that luminous AGNs are biased toward more massive BHs and hence galaxies of
larger stellar masses[123, 1]. Lower mass BHs, which live in less massive
galaxies, are less luminous and thus would be missed from a flux-limited
sample such as ours. Notably, $\sim 40\%$ of the AGNs in the lower redshift
bin and $\sim 53\%$ of the ones in the higher redshift bin have colors bluer
than $(g^{\prime}-r^{\prime})_{0}\approx 0.6$ mag, which indicates the
presence of significant star formation activity alongside BH accretion[124].

To further obtain a rough estimate of the SFR of the AGN host galaxies, we
investigate the correlation between $(g^{\prime}-r^{\prime})_{0}$ color and
specific SFR (${\rm sSFR\equiv SFR}/M_{*}$) in inactive galaxies. We retrieve
the SFRs of our inactive galaxy sample from the second version of the GALEX-
SDSS-WISE Legacy Catalog[125, 126], which were derived from SED fitting of
broad-band photometry from the ultraviolet, optical, and mid-infrared.
Supplementary Fig. 10 shows that the local galaxy population delineates a
sharp “elbow” feature[127, 128], such that $(g^{\prime}-r^{\prime})_{0}$ color
tightly traces sSFR for colors $(g^{\prime}-r^{\prime})_{0}\lesssim 0.65$ mag,
while $(g^{\prime}-r^{\prime})_{0}$ becomes almost insensitive to sSFR at
$(g^{\prime}-r^{\prime})_{0}\gtrsim 0.65$ mag, likely due to the degeneracy
between age, metallicity, and dust extinction[129]. We take advantage of this
strong empirical trend and use the $(g^{\prime}-r^{\prime})_{0}$ color as a
proxy for sSFR, as this color is consistently covered by the PS1 filters
across the entire redshift range of interest for our study. A fit to the star-
forming inactive galaxies, which we define by
$(g^{\prime}-r^{\prime})_{0}<0.65$ mag and ${\rm sSFR}>10^{-11}$ yr-1, with
the linear regression code linmix[12] gives

$\log\,({\rm sSFR}/{\rm yr}^{-1})=-(2.98\pm
0.06)(g^{\prime}-r^{\prime})_{0}-(8.61\pm 0.03),$ (5)

which has a total scatter of 0.26 dex. For $(g^{\prime}-r^{\prime})_{0}>0.65$
mag, we adopt a conservative estimate of ${\rm sSFR}=10^{-12}~{}{\rm
yr}^{-1}$.

###  2.10 Morphologies of the Host Galaxies

The Sérsic index $n$, a flexible, effective representation of a galaxy’s light
distribution, is often used to classify the morphologies of both active and
inactive galaxies into early-type, bulge-dominated (large $n$) and late-type,
disk-dominated (small $n$) systems[130, 131]. However, as the image resolution
can affect the measurement of $n$ even at low redshifts[132], it is necessary
to calibrate the Sérsic index with independent, more accurate morphology
classifications. Comparison of detailed morphological classification of
galaxies imaged with the Hubble Space Telescope[27] reveals that for PS1
images of galaxies at $z\approx 0.3$ the criterion $n=2$ suffices to separate
bulge-dominated from disk-dominated objects. For more nearby objects, we
crossmatch the inactive galaxy sample with the Galaxy Zoo 2 catalog[133, 134]
to obtain more detailed morphology. Galaxy Zoo 2 improves upon the original
Galaxy Zoo classification by measuring more detailed morphological features
(e.g., bulge, bar, spiral arms) in a subset of the brightest and largest
galaxies. In total we obtain reliable morphological classification for 2,259
massive ($M_{*}\geq 10^{10}\,M_{\odot}$) galaxies with $0.01\lesssim z\lesssim
0.25$, including 933 ellipticals, 478 Sa–Sb types (bulge-dominated disks), 686
Sc–Sd types (disk-dominated disks), and 162 edge-on disks. S0 galaxies are
more commonly classified as ellipticals, but with a nonnegligible possibility
as disk galaxies[133]. Supplementary Fig. 11 shows the distributions of
$i_{\rm P1}$-band Sérsic indices for galaxies of different morphological
types. Ellipticals and Sa–Sb galaxies tend to have a high Sérsic index broadly
clustered at $n\approx 4-5$, while Sc–Sd galaxies peak at $n\approx 1-2$,
although all types have significant tails toward lower and higher values of
$n$. To quantify the effectiveness with which galaxy morphology can be
classified crudely using the Sérsic index, we calculate the number ratio (the
higher the better) of Sc–Sd galaxies over galaxies of all other types (i.e.,
ellipticals, Sa–Sb, and edge-on disks). We find that $n\leq 2$ performs better
compared to $n\leq 2.5$, with a number ratio of 2.17 compared to 1.76. This is
consistent with studies of galaxies at slightly higher average redshift
($z\approx 0.3$), which conclude that $n=2$ provides a good separation between
disk-dominated (conservatively Sc–Sd) and bulge-dominated galaxies[27].
Moreover, the high frequency of pseudo bulges in Sc–Sd galaxies—at least $\sim
60\%$ and as high as $100\%$, depending on galaxy sample[135, 136]—justifies
our decision to use the bulge type-dependent calibration[65] to estimate
$M_{\rm BH}$.

Applying the $n=2$ criterion to the AGN sample yields 6,672 ($58\%$) early-
type, bulge-dominated galaxies, 4,116 ($36\%$) late-type, disk-dominated
galaxies, and 686 ($6\%$) galaxies without reliable morphological
classification owing to the large uncertainty in $n$ (measurement/uncertainty
$<3$). In general, early-type hosts are $\sim 0.2$ dex more massive and
dominate over the late-type hosts across the entire redshift range
(Supplementary Fig. 12). The number ratio between early-type and late-type
hosts slightly decreases toward higher redshift, from 2.1 at $z<0.1$ to 1.5 at
$z\geq 0.1$. Closer inspection of the colors of the two types reveals that
they span a similar color range from star-forming to quiescent, with
preference for bluer $(g^{\prime}-r^{\prime})_{0}$ color for late types
compared to early types across the entire range of stellar masses (Extended
Fig. 4), consistent with the results of nearby inactive galaxies[121]. A
number of recent studies have also reported that low-redshift, luminous AGNs
have a significant fraction of late-type host galaxies[137, 138, 139, 79],
highlighting the importance of internal secular processes in regulating AGN
activity[135, 140].

###  2.11 Mass Gain of the BHs and Their Host Galaxies to $z=0$

We obtain the lookback time using the median redshift, estimate sSFR using
Equation (5), and convert bolometric luminosity to BH accretion rate following

$\dot{M}_{\rm BH}=0.15\left(\frac{0.1}{\epsilon}\right)\left(\frac{L_{\rm
bol}}{10^{45}\ {\rm erg\,s^{-1}}}\right)\ M_{\odot}\,{\rm yr}^{-1},$ (6)

with a radiative efficiency $\epsilon=0.1$[141]. Since the sSFR and
$\dot{M}_{\rm BH}$ are median quantities computed over a large number of
objects, we assume that star formation can sustain the current sSFR during the
entire lookback time. We need to consider the duty cycle of the BH accretion,
which can be defined as the fraction of BHs that are actively accreting above
a certain Eddington ratio with respect to all BHs, including inactive ones.
The duty cycle depends on redshift and BH mass[142, 143, 144], with generally
higher values toward higher redshift and smaller BH mass. High Eddington ratio
leads to lower duty cycle, but BHs of lower mass experience a higher duty
cycle. We adopt the relation between duty cycle of AGNs with $\lambda_{\rm
E}>0.01$ and BH mass[144]. The BH mass coverage of the adopted relation spans
the range $M_{\rm BH}\approx 10^{7}-10^{9.6}\,M_{\odot}$, over which the duty
cycle decreases from $\sim 0.1$ to $6\times 10^{-4}$. For $M_{\rm
BH}<10^{7}\,M_{\odot}$, we simply fix the duty cycle to 0.1.

###  2.12 $M_{\rm BH}-M_{*}$ Relation in Local Type 2 AGNs

It is difficult to obtain BH mass measurements in type 2 AGNs, which lack
broad emission lines in optical spectra. Recent observations[145, 146] have
reported detection of near-infrared broad lines in Seyfert 2s, which can be
used to derive BH masses in these systems and offer us the opportunity to
compare their location and evolution on the $M_{\rm BH}-M_{*}$ plane. We
crossmatch the type 2 AGNs presented in Ref.[146] with the PS1 DR2 catalog and
Dark Energy Survey (DES) Data Release 2 catalog[147] to obtain
optical+infrared SEDs. In the end, 14 $z<0.1$ type 2 AGNs have suitable
photometries in five bands $grizy_{\rm P1}$ ($grizY_{\rm DES}$). Assuming that
AGN contamination to the optical and near-infrared continuum is negligible in
type 2 AGNs, we obtain their stellar masses, $(g^{\prime}-r^{\prime})_{0}$
color, and sSFR following the same procedure as used for the host galaxies of
type 1 AGNs. Their BH masses are derived with a virial factor of 6.3, same as
our choice for early-type galaxies. We find that local type 2 AGNs generally
follow the same relation as inactive galaxies and type 1 AGNs, with the caveat
that the current sample is severely limited by small-number statistics
(Extended Fig. 5). Their evolution toward $z=0$, as indicated with red arrows,
is more limited because of the low redshift ($z<0.1$), but it is still similar
to that of type 1 AGNs (Fig. 4).

###  2.13 Cosmology

Throughout the paper, we adopt a cosmology with $H_{0}=70$ km s-1 Mpc-1,
$\Omega_{\rm m}=0.3$, and $\Omega_{\Lambda}=0.7$.

Data Availability

Pan-STARRS1 images are available from the Pan-STARRS1 data archive
https://outerspace.stsci.edu/display/PANSTARRS/Pan-
STARRS1+data+archive+home+page. The basic properties of the final AGN sample,
PS1 flux densities, and the derived physical properties of the host galaxies
are available from https://doi.org/10.12149/101278[149].

Code Availability Astropy, dustmaps, LOESS, Matplotlib, Numpy, and Scipy are
available from python package index (PyPI) https://pypi.org/. CIGALE is
available from https://cigale.lam.fr/. GALFITM is available from
https://www.nottingham.ac.uk/astronomy/megamorph/. linmix is available from
https://linmix.readthedocs.io/en/latest/. SExtractor and SWarp are available
from https://www.astromatic.net/software/.

Acknowledgments This work was supported by the National Key R&D Program of
China (2022YFF0503401), the National Science Foundation of China (11721303,
11991052, 12011540375, 12233001) and the China Manned Space Project (CMS-
CSST-2021-A04, CMS-CSST-2021-A06).

Author Contributions Statement M.-Y. Z. and L. C. H. developed and discussed
the idea. M.-Y. Z. reduced the data, performed the analyses, and wrote the
manuscript. L. C. H. revised the manuscript.

Competing Interests Statement The authors declare that they have no competing
financial interests.

Extended Data

Extended Fig. 1: Rest-frame SDSS $g^{\prime}-r^{\prime}$ color
[$(g^{\prime}-r^{\prime})_{0}$] versus stellar mass ($M_{*}$) for the inactive
galaxy sample (red) and AGN sample (blue). Objects with $z<0.15$ are shown in
panel (a), and objects with $0.15\leq z\leq 0.35$ are shown in panel (b).
Contours indicate the distribution of 10%, 30%, 60%, and 90% of the entire
sample, respectively. Dots are individual objects located outside the 90%
contour. The lower-right corner of each panel shows the typical uncertainties.
Extended Fig. 2: The relation between $M_{\rm BH}$ and $M_{*}$ for $z\leq
0.35$ type 1 AGNs without smoothing. Colors represent (a) rest-frame optical
color $(g^{\prime}-r^{\prime})_{0}$, (b) Eddington ratio ($\lambda_{\rm E}$),
and (c) relative Fe ii strength. Best-fit relations from G20 for local early-
type galaxies (black dashed line), late-type galaxies (black dotted line), and
both types combined (black solid line; intrinsic scatter of 0.81 dex indicated
as shaded gray region). This figure shows the version of Figure 2 without
smoothing to demonstrate that the trends are still present in the original
data. Extended Fig. 3: The relation between $M_{\rm BH}$ and $M_{*}$ for
$z\leq 0.35$ type 1 AGNs without smoothing and with $M_{\rm BH}$ derived from
broad H$\alpha$ emission with a virial factor $f=4.47$. Colors represent (a)
rest-frame optical color $(g^{\prime}-r^{\prime})_{0}$, (b) Eddington ratio
($\lambda_{\rm E}$), and (c) relative Fe ii strength. Best-fit relations from
G20 for local early-type galaxies (black dashed line), late-type galaxies
(black dotted line), and both types combined (black solid line; intrinsic
scatter of 0.81 dex indicated as shaded gray region). This figure demonstrates
that switching to the H$\alpha$-based $M_{\rm BH}$ estimator does not
qualitatively affect our main conclusions. Extended Fig. 4: Rest-frame SDSS
$g^{\prime}-r^{\prime}$ color [$(g^{\prime}-r^{\prime})_{0}$] versus stellar
mass ($M_{*}$) for AGNs. AGNs with early-type and late-type morphologies are
shown in magenta and green, respectively. Objects with $z<0.15$ are shown in
panel (a), and objects with $0.15\leq z\leq 0.35$ are shown in panel (b).
Contours indicate the distribution of 10%, 30%, 60%, and 90% of the entire
sample, respectively. Dots are individual objects located outside the 90%
contour. The lower-right corner of each panel shows the typical uncertainties.
Extended Fig. 5: The relation between $M_{\rm BH}$ and $M_{*}$ for local
($z<0.1$) type 2 AGNs. BH masses are derived from broad emission lines in the
near-infrared (see Methods). Their evolution to $z=0$ is indicated by the
direction and length of each arrow. Best-fit relations from G20 for local
early-type galaxies (black dashed line), late-type galaxies (black dotted
line), and both types combined (black solid line; intrinsic scatter of 0.81
dex indicated as shaded gray region). Typical uncertainties are shown in the
lower-right corner. Extended Fig. 6: Rest-frame optical color
$(g^{\prime}-r^{\prime})_{0}$ versus $\lambda_{\rm E}$ for AGNs. AGNs with
early-type, late-type, and both types combined morphologies are shown in
magenta, green, and black, respectively. Contours indicate the distribution of
10%, 30%, 60%, and 90% of the entire sample, respectively. Dots are individual
objects located outside the 90% contour. Thick solid and dashed lines
represent 50th, 16th, and 84th percentiles of the objects. The Spearman
correlation coefficient $\rho$ and $p$-value are shown in the lower-left
corner. Typical uncertainties are shown in the upper-right corner. Extended
Fig. 7: The relation between $M_{\rm BH}$ and $M_{\rm bulge}$ for $z\leq 0.35$
AGNs with $M_{\rm BH}$ growth fraction larger than $20\%$. Colors represent
growth fraction of $M_{\rm BH}$. $M_{\rm bulge}$ is predicted from $M_{*}$
with a bulge-to-total ratio of 0.2. We show the best-fit relations from G20
for local early-type galaxies (black dashed line), late-type galaxies (black
dotted line), and both types combined (black solid line; intrinsic scatter of
0.81 dex indicated as shaded gray region).

## References

  *   * [1] Kormendy, J. & Ho, L. C.  Coevolution (Or Not) of Supermassive Black Holes and Host Galaxies.  _ARA &A_ 51, 511–653 (2013).
  * [2] Di Matteo, T., Springel, V. & Hernquist, L.  Energy input from quasars regulates the growth and activity of black holes and their host galaxies.  _Nature_ 433, 604–607 (2005).
  * [3] Weinberger, R. _et al._ Supermassive black holes and their feedback effects in the IllustrisTNG simulation.  _MNRAS_ 479, 4056–4072 (2018).
  * [4] Peng, C. Y.  How Mergers May Affect the Mass Scaling Relation between Gravitationally Bound Systems.  _ApJ_ 671, 1098–1107 (2007).
  * [5] Greene, J. E., Strader, J. & Ho, L. C.  Intermediate-Mass Black Holes.  _ARA &A_ 58, 257–312 (2020).
  * [6] Vestergaard, M. & Peterson, B. M.  Determining Central Black Hole Masses in Distant Active Galaxies and Quasars. II. Improved Optical and UV Scaling Relationships.  _ApJ_ 641, 689–709 (2006).
  * [7] Wang, J.-G. _et al._ Estimating Black Hole Masses in Active Galactic Nuclei Using the Mg II $\lambda$2800 Emission Line.  _ApJ_ 707, 1334–1346 (2009).
  * [8] Peng, C. Y. _et al._ Probing the Coevolution of Supermassive Black Holes and Galaxies Using Gravitationally Lensed Quasar Hosts.  _ApJ_ 649, 616–634 (2006).
  * [9] Alexander, D. M. _et al._ Weighing the Black Holes in z $\approx$ 2 Submillimeter-Emitting Galaxies Hosting Active Galactic Nuclei.  _AJ_ 135, 1968–1981 (2008).
  * [10] Schulze, A. & Wisotzki, L.  Accounting for selection effects in the BH-bulge relations: no evidence for cosmological evolution.  _MNRAS_ 438, 3422–3433 (2014).
  * [11] Neeleman, M. _et al._ The Kinematics of z $\gtrsim$ 6 Quasar Host Galaxies.  _ApJ_ 911, 141 (2021).
  * [12] Kelly, B. C.  Some Aspects of Measurement Error in Linear Regression of Astronomical Data.  _ApJ_ 665, 1489–1506 (2007).
  * [13] Bennert, V. N. _et al._ A Local Baseline of the Black Hole Mass Scaling Relations for Active Galaxies. IV. Correlations Between MBH and Host Galaxy $\sigma$, Stellar Mass, and Luminosity.  _ApJ_ 921, 36 (2021).
  * [14] Reines, A. E. & Volonteri, M.  Relations between Central Black Hole Mass and Total Galaxy Stellar Mass in the Local Universe.  _ApJ_ 813, 82 (2015).
  * [15] Kim, M. & Ho, L. C.  Evidence for a Young Stellar Population in Nearby Type 1 Active Galaxies.  _ApJ_ 876, 35 (2019).
  * [16] Ho, L. C.  Nuclear Activity in Nearby Galaxies.  _ARA &A_ 46, 475–539 (2008).
  * [17] Boroson, T. A. & Green, R. F.  The Emission-Line Properties of Low-Redshift Quasi-stellar Objects.  _ApJS_ 80, 109 (1992).
  * [18] Urrutia, T. _et al._ Spitzer Observations of Young Red Quasars.  _ApJ_ 757, 125 (2012).
  * [19] Sun, M. _et al._ Evolution in the Black Hole—Galaxy Scaling Relations and the Duty Cycle of Nuclear Activity in Star-forming Galaxies.  _ApJ_ 802, 14 (2015).
  * [20] Zhao, Y. _et al._ The Diverse Morphology, Stellar Population, and Black Hole Scaling Relations of the Host Galaxies of Nearby Quasars.  _ApJ_ 911, 94 (2021).
  * [21] Mathur, S., Fields, D., Peterson, B. M. & Grupe, D.  Supermassive Black Holes, Pseudobulges, and the Narrow-line Seyfert 1 Galaxies.  _ApJ_ 754, 146 (2012).
  * [22] Greene, J. E. & Ho, L. C.  The MBH-$\sigma$∗ Relation in Local Active Galaxies.  _ApJ_ 641, L21–L24 (2006).
  * [23] Madau, P. & Dickinson, M.  Cosmic Star-Formation History.  _ARA &A_ 52, 415–486 (2014).
  * [24] Graham, A. W. & Sahu, N.  Appreciating mergers for understanding the non-linear Mbh-M∗,spheroid and Mbh-M∗,galaxy relations, updated herein, and the implications for the (reduced) role of AGN feedback.  _MNRAS_ 518, 2177–2200 (2023).
  * [25] Steinborn, L. K. _et al._ Cosmological simulations of black hole growth II: how (in)significant are merger events for fuelling nuclear activity?  _MNRAS_ 481, 341–360 (2018).
  * [26] Zhuang, M.-Y. & Ho, L. C.  The Interplay between Star Formation and Black Hole Accretion in Nearby Active Galaxies.  _ApJ_ 896, 108 (2020).
  * [27] Zhuang, M.-Y. & Ho, L. C.  The Star-forming Main Sequence of the Host Galaxies of Low-redshift Quasars.  _ApJ_ 934, 130 (2022).
  * [28] Farrah, D. _et al._ Stellar and black hole assembly in $z<0.3$ infrared-luminous mergers: intermittent starbursts versus super-Eddington accretion.  _MNRAS_ 513, 4770–4786 (2022).
  * [29] Borys, C. _et al._ The Relationship between Stellar and Black Hole Mass in Submillimeter Galaxies.  _ApJ_ 635, 853–863 (2005).
  * [30] Veilleux, S. _et al._ A Deep Hubble Space Telescope H-Band Imaging Survey of Massive Gas-rich Mergers.  _ApJ_ 643, 707–723 (2006).
  * [31] Chen, C.-C. _et al._ JWST Sneaks a Peek at the Stellar Morphology of $z\sim 2$ Submillimeter Galaxies: Bulge Formation at Cosmic Noon.  _ApJ_ 939, L7 (2022).
  * [32] Greve, T. R. _et al._ An interferometric CO survey of luminous submillimetre galaxies.  _MNRAS_ 359, 1165–1183 (2005).
  * [33] Shangguan, J. _et al._ Interstellar Medium and Star Formation of Starburst Galaxies on the Merger Sequence.  _ApJ_ 870, 104 (2019).
  * [34] Orban de Xivry, G. _et al._ The role of secular evolution in the black hole growth of narrow-line Seyfert 1 galaxies.  _MNRAS_ 417, 2721–2736 (2011).
  * [35] Shangguan, J., Ho, L. C., Bauer, F. E., Wang, R. & Treister, E.  AGN Feedback and Star Formation of Quasar Host Galaxies: Insights from the Molecular Gas.  _ApJ_ 899, 112 (2020).
  * [36] Xie, Y., Ho, L. C., Zhuang, M.-Y. & Shangguan, J.  The Infrared Emission and Vigorous Star Formation of Low-redshift Quasars.  _ApJ_ 910, 124 (2021).
  * [37] Kennicutt, J., Robert C.  The Global Schmidt Law in Star-forming Galaxies.  _ApJ_ 498, 541–552 (1998).
  * [38] Martini, P.  Ho, L. C. (ed.) _QSO Lifetimes_.  (ed.Ho, L. C.) _Coevolution of Black Holes and Galaxies_ , 169 (2004).  astro-ph/0304009.
  * [39] Habouzit, M., Volonteri, M. & Dubois, Y.  Blossoms from black hole seeds: properties and early growth regulated by supernova feedback.  _MNRAS_ 468, 3935–3948 (2017).
  * [40] Çatmabacak, O. _et al._ Black hole-galaxy scaling relations in FIRE: the importance of black hole location and mergers.  _MNRAS_ 511, 506–535 (2022).
  * [41] Koss, M. J. _et al._ BAT AGN Spectroscopic Survey. XX. Molecular Gas in Nearby Hard-X-Ray-selected AGN Galaxies.  _ApJS_ 252, 29 (2021).
  * [42] Molina, J. _et al._ Compact Molecular Gas Distribution in Quasar Host Galaxies.  _ApJ_ 908, 231 (2021).
  * [43] Aird, J., Coil, A. L. & Georgakakis, A.  X-rays across the galaxy population - III. The incidence of AGN as a function of star formation rate.  _MNRAS_ 484, 4360–4378 (2019).
  * [44] Decarli, R. _et al._ An ALMA [C II] Survey of 27 Quasars at $z>5.94$.  _ApJ_ 854, 97 (2018).
  * [45] Inayoshi, K., Visbal, E. & Haiman, Z.  The Assembly of the First Massive Black Holes.  _ARA &A_ 58, 27–97 (2020).
  * [46] Daddi, E. _et al._ Passively Evolving Early-Type Galaxies at $1.4\lesssim z\lesssim 2.5$ in the Hubble Ultra Deep Field.  _ApJ_ 626, 680–697 (2005).
  * [47] van Dokkum, P. G. _et al._ The Growth of Massive Galaxies Since z = 2.  _ApJ_ 709, 1018–1041 (2010).
  * [48] Oser, L., Ostriker, J. P., Naab, T., Johansson, P. H. & Burkert, A.  The Two Phases of Galaxy Formation.  _ApJ_ 725, 2312–2323 (2010).
  * [49] Davari, R. H., Ho, L. C., Mobasher, B. & Canalizo, G.  Detection of Prominent Stellar Disks in the Progenitors of Present-day Massive Elliptical Galaxies.  _ApJ_ 836, 75 (2017).
  * [50] Zhuang, M.-Y., Ho, L. C. & Shangguan, J.  Black Hole Accretion Correlates with Star Formation Rate and Star Formation Efficiency in Nearby Luminous Type 1 Active Galaxies.  _ApJ_ 906, 38 (2021).
  * [51] Ward, S. R., Harrison, C. M., Costa, T. & Mainieri, V.  Cosmological simulations predict that AGN preferentially live in gas-rich, star-forming galaxies despite effective feedback.  _MNRAS_ 514, 2936–2957 (2022).
  * [52] King, A. & Pounds, K.  Powerful Outflows and Feedback from Active Galactic Nuclei.  _ARA &A_ 53, 115–154 (2015).
  * [53] Tremaine, S. _et al._ The Slope of the Black Hole Mass versus Velocity Dispersion Correlation.  _ApJ_ 574, 740–753 (2002).
  * [54] Gao, H., Ho, L. C., Barth, A. J. & Li, Z.-Y.  The Carnegie-Irvine Galaxy Survey. VIII. Demographics of Bulges along the Hubble Sequence.  _ApJS_ 244, 34 (2019).
  * [55] Ding, Y., Li, R., Ho, L. C. & Ricci, C.  Accretion Disk Outflow during the X-Ray Flare of the Super-Eddington Active Nucleus of I Zwicky 1.  _ApJ_ 931, 77 (2022).
  * [56] Davé, R. _et al._ SIMBA: Cosmological simulations with black hole growth and feedback.  _MNRAS_ 486, 2827–2849 (2019).
  * [57] Liu, H.-Y. _et al._ A Comprehensive and Uniform Sample of Broad-line Active Galactic Nuclei from the SDSS DR7.  _ApJS_ 243, 21 (2019).
  * [58] Abazajian, K. N. _et al._ The Seventh Data Release of the Sloan Digital Sky Survey.  _ApJS_ 182, 543–558 (2009).
  * [59] York, D. G. _et al._ The Sloan Digital Sky Survey: Technical Summary.  _AJ_ 120, 1579–1587 (2000).
  * [60] Shen, Y. _et al._ A Catalog of Quasar Properties from Sloan Digital Sky Survey Data Release 7.  _ApJS_ 194, 45 (2011).
  * [61] Oh, K. _et al._ A New Catalog of Type 1 AGNs and its Implications on the AGN Unified Model.  _ApJS_ 219, 1 (2015).
  * [62] Lyke, B. W. _et al._ The Sloan Digital Sky Survey Quasar Catalog: Sixteenth Data Release.  _ApJS_ 250, 8 (2020).
  * [63] Rakshit, S., Stalin, C. S. & Kotilainen, J.  Spectral Properties of Quasars from Sloan Digital Sky Survey Data Release 14: The Catalog.  _ApJS_ 249, 17 (2020).
  * [64] McLure, R. J. & Dunlop, J. S.  The cosmological evolution of quasar black hole masses.  _MNRAS_ 352, 1390–1404 (2004).
  * [65] Ho, L. C. & Kim, M.  A Revised Calibration of the Virial Mass Estimator for Black Holes in Active Galaxies Based on Single-epoch H$\beta$ Spectra.  _ApJ_ 809, 123 (2015).
  * [66] Ho, L. C. & Kim, M.  The Black Hole Mass Scale of Classical and Pseudo Bulges in Active Galaxies.  _ApJ_ 789, 17 (2014).
  * [67] Greene, J. E. & Ho, L. C.  Estimating Black Hole Masses in Active Galaxies Using the H$\alpha$ Emission Line.  _ApJ_ 630, 122–129 (2005).
  * [68] Woo, J.-H., Yoon, Y., Park, S., Park, D. & Kim, S. C.  The Black Hole Mass-Stellar Velocity Dispersion Relation of Narrow-line Seyfert 1 Galaxies.  _ApJ_ 801, 38 (2015).
  * [69] Du, P. _et al._ Supermassive Black Holes with High Accretion Rates in Active Galactic Nuclei. V. A New Size-Luminosity Scaling Relation for the Broad-line Region.  _ApJ_ 825, 126 (2016).
  * [70] Du, P. & Wang, J.-M.  The Radius-Luminosity Relationship Depends on Optical Spectra in Active Galactic Nuclei.  _ApJ_ 886, 42 (2019).
  * [71] Dalla Bontà, E. _et al._ The Sloan Digital Sky Survey Reverberation Mapping Project: Estimating Masses of Black Holes in Quasars with Single-epoch Spectroscopy.  _ApJ_ 903, 112 (2020).
  * [72] Fonseca Alvarez, G. _et al._ The Sloan Digital Sky Survey Reverberation Mapping Project: The H$\beta$ Radius-Luminosity Relation.  _ApJ_ 899, 73 (2020).
  * [73] Treu, T., Woo, J.-H., Malkan, M. A. & Blandford, R. D.  Cosmic Evolution of Black Holes and Spheroids. II. Scaling Relations at z=0.36.  _ApJ_ 667, 117–130 (2007).
  * [74] Bennert, V. N. _et al._ Cosmic Evolution of Black Holes and Spheroids. IV. The M BH-L sph Relation.  _ApJ_ 708, 1507–1527 (2010).
  * [75] Falomo, R., Bettoni, D., Karhunen, K., Kotilainen, J. K. & Uslenghi, M.  Low-redshift quasars in the Sloan Digital Sky Survey Stripe 82. The host galaxies.  _MNRAS_ 440, 476–493 (2014).
  * [76] Matsuoka, Y., Strauss, M. A., Price, I., Ted N. & DiDonato, M. S.  Massive Star-forming Host Galaxies of Quasars on Sloan Digital Sky Survey Stripe 82.  _ApJ_ 780, 162 (2014).
  * [77] Bentz, M. C. & Manne-Nicholas, E.  Black Hole-Galaxy Scaling Relationships for Active Galactic Nuclei with Reverberation Masses.  _ApJ_ 864, 146 (2018).
  * [78] Ishino, T. _et al._ Subaru Hyper Suprime-Cam view of quasar host galaxies at $z<1$.  _PASJ_ 72, 83 (2020).
  * [79] Li, J. _et al._ The Sizes of Quasar Host Galaxies in the Hyper Suprime-Cam Subaru Strategic Program.  _ApJ_ 918, 22 (2021).
  * [80] Li, J. I. H. _et al._ The Sloan Digital Sky Survey Reverberation Mapping Project: The MBH-Host Relations at $0.2\lesssim z\lesssim 0.6$ from Reverberation Mapping and Hubble Space Telescope Imaging.  _ApJ_ 906, 103 (2021).
  * [81] Son, S., Kim, M., Barth, A. J. & Ho, L. C.  Relation between Black Hole Mass and Bulge Luminosity in Hard X-Ray Selected Type 1 AGNs.  _Journal of Korean Astronomical Society_ 55, 37–57 (2022).
  * [82] Heckman, T. M.  An Optical and Radio Survey of the Nuclei of Bright Galaxies - Activity in the Normal Galactic Nuclei.  _A &A_ 87, 152 (1980).
  * [83] Baldwin, J. A., Phillips, M. M. & Terlevich, R.  Classification parameters for the emission-line spectra of extragalactic objects.  _PASP_ 93, 5–19 (1981).
  * [84] Chambers, K. C. _et al._ The Pan-STARRS1 Surveys.  _arXiv e-prints_ arXiv:1612.05560 (2016).
  * [85] Kron, R. G.  Photometry of a complete sample of faint galaxies.  _ApJS_ 43, 305–325 (1980).
  * [86] Flewelling, H. A. _et al._ The Pan-STARRS1 Database and Data Products.  _ApJS_ 251, 7 (2020).
  * [87] Bertin, E. & Arnouts, S.  SExtractor: Software for source extraction.  _A &AS_ 117, 393–404 (1996).
  * [88] Bertin, E. _et al._ Bohlender, D. A., Durand, D. & Handley, T. H. (eds) _The TERAPIX Pipeline_.  (eds Bohlender, D. A., Durand, D. & Handley, T. H.) _Astronomical Data Analysis Software and Systems XI_ , Vol. 281 of _ASP Conf. Ser._ , 228 (San Francisco, CA: ASP, 2002).
  * [89] Farrow, D. J. _et al._ Pan-STARRS1: Galaxy clustering in the Small Area Survey 2.  _MNRAS_ 437, 748–770 (2014).
  * [90] Moffat, A. F. J.  A Theoretical Investigation of Focal Stellar Images in the Photographic Emulsion and Application to Photographic Photometry.  _A &A_ 3, 455 (1969).
  * [91] Häußler, B. _et al._ MegaMorph - multiwavelength measurement of galaxy structure: complete Sérsic profile information from modern surveys.  _MNRAS_ 430, 330–369 (2013).
  * [92] Vika, M. _et al._ MegaMorph - multiwavelength measurement of galaxy structure. Sérsic profile fits to galaxies near and far.  _MNRAS_ 435, 623–649 (2013).
  * [93] Sersic, J. L.  _Atlas de Galaxias Australes_ (1968).
  * [94] Kelvin, L. S. _et al._ Galaxy And Mass Assembly (GAMA): Structural Investigation of Galaxies via Model Analysis.  _MNRAS_ 421, 1007–1039 (2012).
  * [95] Peng, C. Y., Ho, L. C., Impey, C. D. & Rix, H.-W.  Detailed Structural Decomposition of Galaxy Images.  _AJ_ 124, 266–293 (2002).
  * [96] Peng, C. Y., Ho, L. C., Impey, C. D. & Rix, H.-W.  Detailed Decomposition of Galaxy Images. II. Beyond Axisymmetric Models.  _AJ_ 139, 2097–2129 (2010).
  * [97] Häußler, B. _et al._ GALAPAGOS-2/GALFITM/GAMA - Multi-wavelength measurement of galaxy structure: Separating the properties of spheroid and disk components in modern surveys.  _A &A_ 664, A92 (2022).
  * [98] Zhuang, M.-Y. & Shen, Y.  Characterization of JWST NIRCam PSFs and Implications for AGN+Host Image Decomposition.  _arXiv e-prints_ arXiv:2304.13776 (2023).
  * [99] Casura, S. _et al._ Galaxy And Mass Assembly (GAMA): bulge-disc decomposition of KiDS data in the nearby Universe.  _MNRAS_ 516, 942–974 (2022).
  * [100] Häussler, B. _et al._ GEMS: Galaxy Fitting Catalogs and Testing Parametric Galaxy Fitting Codes: GALFIT and GIM2D.  _ApJS_ 172, 615–633 (2007).
  * [101] Gao, H. & Ho, L. C.  An Optimal Strategy for Accurate Bulge-to-disk Decomposition of Disk Galaxies.  _ApJ_ 845, 114 (2017).
  * [102] Fan, L. _et al._ Structure and Morphology of X-Ray-selected Active Galactic Nucleus Hosts at $1<z<3$ in the CANDELS-COSMOS Field.  _ApJ_ 784, L9 (2014).
  * [103] Simard, L., Mendel, J. T., Patton, D. R., Ellison, S. L. & McConnachie, A. W.  A Catalog of Bulge+disk Decompositions and Updated Photometry for 1.12 Million Galaxies in the Sloan Digital Sky Survey.  _ApJS_ 196, 11 (2011).
  * [104] Richards, G. T. _et al._ Spectral Energy Distributions and Multiwavelength Selection of Type 1 Quasars.  _ApJS_ 166, 470–497 (2006).
  * [105] Boquien, M. _et al._ CIGALE: a python Code Investigating GALaxy Emission.  _A &A_ 622, A103 (2019).
  * [106] Yang, G. _et al._ Fitting AGN/Galaxy X-Ray-to-radio SEDs with CIGALE and Improvement of the Code.  _ApJ_ 927, 192 (2022).
  * [107] Green, G. M.  dustmaps: A Python interface for maps of interstellar dust.  _The Journal of Open Source Software_ 3, 695 (2018).
  * [108] Schlegel, D. J., Finkbeiner, D. P. & Davis, M.  Maps of Dust Infrared Emission for Use in Estimation of Reddening and Cosmic Microwave Background Radiation Foregrounds.  _ApJ_ 500, 525–553 (1998).
  * [109] Schlafly, E. F. & Finkbeiner, D. P.  Measuring Reddening with Sloan Digital Sky Survey Stellar Spectra and Recalibrating SFD.  _ApJ_ 737, 103 (2011).
  * [110] Bruzual, G. & Charlot, S.  Stellar population synthesis at the resolution of 2003\.  _MNRAS_ 344, 1000–1028 (2003).
  * [111] Inoue, A. K.  Rest-frame ultraviolet-to-optical spectral characteristics of extremely metal-poor and metal-free galaxies.  _MNRAS_ 415, 2920–2931 (2011).
  * [112] Calzetti, D. _et al._ The Dust Content and Opacity of Actively Star-forming Galaxies.  _ApJ_ 533, 682–695 (2000).
  * [113] Mitchell, P. D., Lacey, C. G., Baugh, C. M. & Cole, S.  How well can we really estimate the stellar masses of galaxies from broad-band photometry?  _MNRAS_ 435, 87–114 (2013).
  * [114] Ciesla, L. _et al._ Constraining the properties of AGN host galaxies with spectral energy distribution modelling.  _A &A_ 576, A10 (2015).
  * [115] Chabrier, G.  Galactic Stellar and Substellar Initial Mass Function.  _PASP_ 115, 763–795 (2003).
  * [116] Kroupa, P.  On the variation of the initial mass function.  _MNRAS_ 322, 231–246 (2001).
  * [117] Bell, E. F., McIntosh, D. H., Katz, N. & Weinberg, M. D.  The Optical and Near-Infrared Properties of Galaxies. I. Luminosity and Stellar Mass Functions.  _ApJS_ 149, 289–312 (2003).
  * [118] Stalevski, M., Fritz, J., Baes, M., Nakos, T. & Popović, L. Č.  3D radiative transfer modelling of the dusty tori around active galactic nuclei as a clumpy two-phase medium.  _MNRAS_ 420, 2756–2772 (2012).
  * [119] Stalevski, M. _et al._ The dust covering factor in active galactic nuclei.  _MNRAS_ 458, 2288–2302 (2016).
  * [120] Baldry, I. K. _et al._ Quantifying the Bimodal Color-Magnitude Distribution of Galaxies.  _ApJ_ 600, 681–694 (2004).
  * [121] Blanton, M. R. & Moustakas, J.  Physical Properties and Environments of Nearby Galaxies.  _ARA &A_ 47, 159–210 (2009).
  * [122] Blanton, M. R. _et al._ The Broadband Optical Properties of Galaxies with Redshifts $0.02<z<0.22$.  _ApJ_ 594, 186–207 (2003).
  * [123] Magorrian, J. _et al._ The Demography of Massive Dark Objects in Galaxy Centers.  _AJ_ 115, 2285–2305 (1998).
  * [124] Trump, J. R. _et al._ A Census of Broad-line Active Galactic Nuclei in Nearby Galaxies: Coeval Star Formation and Rapid Black Hole Growth.  _ApJ_ 763, 133 (2013).
  * [125] Salim, S. _et al._ GALEX-SDSS-WISE Legacy Catalog (GSWLC): Star Formation Rates, Stellar Masses, and Dust Attenuations of 700,000 Low-redshift Galaxies.  _ApJS_ 227, 2 (2016).
  * [126] Salim, S., Boquien, M. & Lee, J. C.  Dust Attenuation Curves in the Local Universe: Demographics and New Laws for Star-forming Galaxies and High-redshift Analogs.  _ApJ_ 859, 11 (2018).
  * [127] Weinmann, S. M., van den Bosch, F. C., Yang, X. & Mo, H. J.  Properties of galaxy groups in the Sloan Digital Sky Survey - I. The dependence of colour, star formation and morphology on halo mass.  _MNRAS_ 366, 2–28 (2006).
  * [128] Huang, S., Haynes, M. P., Giovanelli, R. & Brinchmann, J.  The Arecibo Legacy Fast ALFA Survey: The Galaxy Population Detected by ALFALFA.  _ApJ_ 756, 113 (2012).
  * [129] Brinchmann, J. _et al._ The physical properties of star-forming galaxies in the low-redshift Universe.  _MNRAS_ 351, 1151–1179 (2004).
  * [130] Ravindranath, S. _et al._ The Evolution of Disk Galaxies in the GOODS-South Field: Number Densities and Size Distribution.  _ApJ_ 604, L9–L12 (2004).
  * [131] Gabor, J. M. _et al._ Active Galactic Nucleus Host Galaxy Morphologies in COSMOS.  _ApJ_ 691, 705–722 (2009).
  * [132] Dewsnap, C. _et al._ GALFIT-ing AGN Host Galaxies in COSMOS: HST vs. Subaru.  _ApJ_ 944, 137 (2023).
  * [133] Willett, K. W. _et al._ Galaxy Zoo 2: detailed morphological classifications for 304 122 galaxies from the Sloan Digital Sky Survey.  _MNRAS_ 435, 2835–2860 (2013).
  * [134] Hart, R. E. _et al._ Galaxy Zoo: comparing the demographics of spiral arm number and a new method for correcting redshift bias.  _MNRAS_ 461, 3663–3682 (2016).
  * [135] Kormendy, J. & Kennicutt, J., Robert C.  Secular Evolution and the Formation of Pseudobulges in Disk Galaxies.  _ARA &A_ 42, 603–683 (2004).
  * [136] Gao, H., Ho, L. C., Barth, A. J. & Li, Z.-Y.  The Carnegie-Irvine Galaxy Survey. IX. Classification of Bulge Types and Statistical Properties of Pseudo Bulges.  _ApJS_ 247, 20 (2020).
  * [137] Kim, M., Ho, L. C., Peng, C. Y., Barth, A. J. & Im, M.  Stellar Photometric Structures of the Host Galaxies of Nearby Type 1 Active Galactic Nuclei.  _ApJS_ 232, 21 (2017).
  * [138] Kim, M., Barth, A. J., Ho, L. C. & Son, S.  A Hubble Space Telescope Imaging Survey of Low-redshift Swift-BAT Active Galaxies.  _ApJS_ 256, 40 (2021).
  * [139] Zhao, D., Ho, L. C., Zhao, Y., Shangguan, J. & Kim, M.  The Role of Major Mergers and Nuclear Star Formation in Nearby Obscured Quasars.  _ApJ_ 877, 52 (2019).
  * [140] Martín-Navarro, I., Shankar, F. & Mezcua, M.  Rejuvenation triggers nuclear activity in nearby galaxies.  _MNRAS_ 513, L10–L14 (2022).
  * [141] Frank, J., King, A. & Raine, D.  _Accretion power in astrophysics._ (Cambridge: Cambridge Univ. Press, 1992).
  * [142] Aversa, R., Lapi, A., de Zotti, G., Shankar, F. & Danese, L.  Black Hole and Galaxy Coevolution from Continuity Equation and Abundance Matching.  _ApJ_ 810, 74 (2015).
  * [143] Delvecchio, I. _et al._ The Evolving AGN Duty Cycle in Galaxies Since z $\sim$ 3 as Encoded in the X-Ray Luminosity Function.  _ApJ_ 892, 17 (2020).
  * [144] Ananna, T. T. _et al._ BASS. XXX. Distribution Functions of DR2 Eddington Ratios, Black Hole Masses, and X-Ray Luminosities.  _ApJS_ 261, 9 (2022).
  * [145] Lamperti, I. _et al._ BAT AGN Spectroscopic Survey - IV: Near-Infrared Coronal Lines, Hidden Broad Lines, and Correlation with Hard X-ray Emission.  _MNRAS_ 467, 540–572 (2017).
  * [146] Ricci, F. _et al._ BASS. XXIX. The Near-infrared View of the Broad-line Region (BLR): The Effects of Obscuration in BLR Characterization.  _ApJS_ 261, 8 (2022).
  * [147] Abbott, T. M. C. _et al._ The Dark Energy Survey Data Release 2.  _ApJS_ 255, 20 (2021).
  * [148] Cappellari, M. _et al._ The ATLAS3D project - XX. Mass-size and mass-$\sigma$ distributions of early-type galaxies: bulge fraction drives kinematics, mass-to-light ratio, molecular gas fraction and stellar initial mass function.  _MNRAS_ 432, 1862–1893 (2013).
  * [149] Zhuang, M.-Y. & Ho, L. C.  Basic Information and Derived Properties of the Final AGN Sample.  _National Astronomical Data Center of China_ https://doi.org/10.12149/101278 (2023).

See pages 1-9 of Supplementary_Information