#  Little Red Dots: Rapidly Growing Black Holes Reddened by Extended Dusty
Flows

Zhengrong Li Kavli Institute for Astronomy and Astrophysics, Peking
University, Beijing 100871, China  Department of Astronomy, School of Physics,
Peking University, Beijing 100871, China  Kohei Inayoshi Kavli Institute for
Astronomy and Astrophysics, Peking University, Beijing 100871, China  Kejian
Chen Kavli Institute for Astronomy and Astrophysics, Peking University,
Beijing 100871, China  Department of Astronomy, School of Physics, Peking
University, Beijing 100871, China  Kohei Ichikawa Global Center for Science
and Engineering, Faculty of Science and Engineering, Waseda University, 3-4-1,
Okubo, Shinjuku, Tokyo 169-8555, Japan  Department of Physics, School of
Advanced Science and Engineering, Faculty of Science and Engineering, Waseda
University, 3-4-1, Okubo, Shinjuku, Tokyo 169-8555, Japan  Luis C. Ho Kavli
Institute for Astronomy and Astrophysics, Peking University, Beijing 100871,
China  Department of Astronomy, School of Physics, Peking University, Beijing
100871, China

###### Abstract

The James Webb Space Telescope (JWST) observations have revolutionized
extragalactic research, particularly with the discovery of little red dots
(LRD), which we propose are dust-reddened broad-line active galactic nuclei
(AGNs). Their unique v-shape spectral feature observed through JWST/NIRCam
challenges us to discern the relative contributions of the galaxy and AGN. We
study a spectral energy distribution (SED) model for LRDs from rest-frame UV
to infrared bands. We hypothesize that the incident radiation from an AGN,
characterized by a typical SED, is embedded in an extended dusty medium with
an extinction law similar to those seen in dense regions such as Orion Nebula
or certain AGN environments. The UV-optical spectrum is described by dust-
attenuated AGN emission, featuring a red optical continuum at
$\lambda>4000~{}{\rm\AA}$ and a flat UV spectral shape established through a
gray extinction curve at $\lambda<3000~{}{\rm\AA}$, due to the absence of
small-size grains. There is no need for additional stellar emission or AGN
scattered light. In the infrared, the SED is shaped by an extended dust and
gas distribution ($\gamma<1$; $\rho\propto r^{-\gamma}$) with a characteristic
gas density of $\simeq 10-10^{3}~{}{\rm cm}^{-3}$, which allows relatively
cool dust temperatures to dominate the radiation, thereby shifting the energy
peak from near- to mid-infrared bands. This model, unlike the typical AGN hot
torus models, can produce an infrared SED flattening that is consistent with
LRD observations through JWST MIRI. Such a density structure can arise from
the coexistence of inflows and outflows during the early assembly of galactic
nuclei. This might be the reason why LRDs emerge preferentially in the high-
redshift universe younger than one billion years.

Galaxy formation (595); High-redshift galaxies (734); Quasars (1319);
Supermassive black holes (1663)

##  1 introduction

The James Webb Space Telescope (JWST) has dramatically advanced the study of
high-redshift ($z>4$) active galactic nuclei (AGNs) by unveiling low-
luminosity AGN populations more representative than luminous quasars (e.g.,
Onoue et al., 2023). Among these discoveries, very red and compact objects,
so-called little red dots (LRDs), have attracted considerable attention (e.g.,
Barro et al., 2023; Kocevski et al., 2023; Labbe et al., 2023; Matthee et al.,
2023). Spectroscopic observations on photometrically identified LRDs have
revealed that more than 70% exhibit broad hydrogen Balmer lines (e.g., Greene
et al., 2024), indicating the presence of massive black holes (BHs). The
abundance of these AGNs in cosmic volume (Kokorev et al., 2023; Kocevski et
al., 2024; Akins et al., 2024) is significantly higher than what is expected
from ground-based quasar surveys (e.g., Matsuoka et al., 2018, 2023; Niida et
al., 2020; He et al., 2023), and provides insights on the cosmological
evolution of quasar luminosity functions (Li et al., 2023), the impact on
cosmic reionization (e.g., Dayal et al., 2024; Madau et al., 2024), and the
implications to rapid spins of these high-redshift BHs through the measurement
of the radiation-to-mass conversion efficiency (Inayoshi & Ichikawa, 2024).

Despite their significance, the properties of these newly identified LRDs
remain puzzling, particularly concerning the origin of their characteristic
v-shaped spectral energy distribution (SED) (e.g., Barro et al., 2023;
Kocevski et al., 2023; Labbe et al., 2023). These SEDs feature a red continuum
in the rest-frame optical spectrum at $\gtrsim 4000~{}{\rm\AA}$ and a blue
excess in the rest-frame UV range. Observational studies using SED template
fitting suggested that the red optical part could either be due to dust-
reddened AGNs or starburst galaxies, while the blue UV component might
represent light from unobscured galaxies or scattered light from embedded
AGNs.

Furthermore, rest-frame near- to far-infrared observations offer additional
insights. Data from the JADES, COSMOS-Web and PRIMER-COSMOS surveys show that
most LRDs are not detected in the JWST/MIRI F1800W bands, and even for the
detected sources, the MIRI SEDs tend to flatten from the rest-frame optical to
the near-infrared (NIR) at 1–3 $\mu$m (e.g., Williams et al. 2024; Akins et
al. 2024; see also Wang et al. 2024). Those observations appear to contradict
expectations that the SED should be very red throughout the optical to NIR
spectrum due to re-emission from hot dust components heated by the AGNs.
Consequently, these results have led to arguments disfavoring the scenario
accounting for dust-obscured AGNs. While a certain fraction of LRDs show
extremely red SEDs from NIR to mid-infrared (MIR) (Lyu et al., 2024; Pérez-
González et al., 2024), the reasons for such diversity at $\lambda\gtrsim
1~{}\mu\mathrm{m}$ are not well understood. Alternatively, if considering
dust-obscured galaxies instead of AGNs, cold dust associated with starbursts
would dominate the far-infrared (FIR) light. However, current ALMA
observations provide no strong evidence to support this scenario (Labbe et
al., 2023; Williams et al., 2024; Akins et al., 2024; Casey et al., 2024).

It is worth noting that LRDs are extremely compact (effective radius of
$R_{\mathrm{eff}}<100$ pc at $z\sim 5$; Labbe et al., 2023; Furtak et al.,
2023) and almost point-like objects even with supreme spatial resolution by
JWST. A small fraction of red sources with a v-shaped SED show extended
morphology in images of shorter-wavelength filters, which may indicate the
presence of an underlying host galaxy (Killi et al., 2023) or diffuse ionized
gas. However, most LRDs do not have detectable extended components, allowing
us to place upper limits on the mass and size of the galaxies. Comparing the
BH mass estimated from broad-line emission spectroscopy to the upper limits on
galaxy mass, the BH-to-galaxy mass ratio in LRDs (e.g., Pacucci et al., 2023;
Kokorev et al., 2023; Greene et al., 2024) is substantially higher than the
empirical value observed in the nearby universe (e.g., Kormendy & Ho, 2013).
This implies the existence of an evolutionary stage where growing BHs occupy
$\sim{\rm a~{}few}-10\%$ of the total BH+galaxy system (for theoretical
results see e.g., Inayoshi et al., 2022; Hu et al., 2022; Scoggins et al.,
2023).

In this paper, we study a SED model for LRDs from rest-frame UV to infrared
bands. We hypothesize that the incident radiation from an AGN, characterized
by a typical SED, is embedded in an extended dusty medium with an extinction
law similar to those seen in dense regions such as Orion Nebula or certain AGN
environments. The UV-optical spectrum is described by dust-attenuated AGN
emission, composed of a red optical continuum and a UV excess established
through gray extinction. There is no need for additional stellar emission or
AGN scattered light. In the infrared, the SED is shaped by an extended dust
and gas distribution, which allows relatively cool dust temperatures to
dominate the radiation, thereby shifting the energy peak from near- to mid-
infrared bands. This model can produce an infrared SED flattening that is
consistent with LRD observations through JWST MIRI. We finally provide
discussion and open questions rising from our conclusions regarding (i) the
temporal emergence of LRDs at $z>4$, (ii) the universality of v-shaped SEDs,
and (iii) the cosmological evolution of early BH populations.

##  2 The “v”-shape SED in rest-frame UV/optical

It is evident that dust extinction curves in AGN significantly deviate from
those observed in the Milky Way (e.g., Maiolino et al., 2001b, a, 2004;
Gaskell et al., 2004; Li, 2008; Xie et al., 2017), showing a characteristic
flat wavelength dependence at $\lambda<4000$ Å. This feature is believed to
stem from a deficiency of small-size grains within the AGN environment, which
can be attributed to dust destruction mechanisms such as thermal sublimation
and sputtering (e.g., Laor & Draine, 1993; Tazaki & Ichikawa, 2020), and grain
charging (e.g., Draine & Salpeter, 1979). Recently, Tazaki et al. (2020) have
applied these dust destruction processes specifically at the polar regions of
$1-10~{}{\rm pc}$ directly irradiated by the AGN (e.g., Hönig & Kishimoto,
2017), and show the critical role of Coulomb explosions in flattening the
extinction curve.

Figure 1: Extinction curves as a function of wavelength normalized by the
visual extinction at 5500 $\rm\AA$: the Small Magellanic Cloud (SMC, Gordon et
al., 2003), starburst galaxies (Calzetti et al., 2000), the interstellar dust
in Milky Way with $R_{\mathrm{V}}=3.1$ and 5.0 (Cardelli et al., 1989), high-
redshift galaxies at $6<z<12$ (Markov et al., 2024), the Orion Nebula (Baldwin
et al., 1991), and composite AGN spectra (Gaskell et al., 2004). Extinction
curves in the absence of small-size dust grains, as observed in the Orion
Nebula, are significantly flattened at $\lambda<2500~{}{\rm\AA}$, unlike the
other models commonly used in galaxy SED fitting methodology (e.g. Calzetii’s
law in Labbe et al., 2023).

Figure 1 presents seven different dust extinction laws, normalized at the
visual extinction level at $5500~{}{\rm\AA}$ ($A_{\lambda}/A_{\mathrm{V}}$):
the Small Magellanic Cloud (SMC, Gordon et al., 2003), starburst galaxies
(Calzetti et al., 2000), the interstellar dust in Milky Way with
$R_{\mathrm{V}}=3.1$ and 5.0 (Cardelli et al., 1989), high-redshift galaxies
at $6<z<12$ (Markov et al., 2024), the Orion Nebula (Baldwin et al., 1991),
and composite AGN spectra (Gaskell et al., 2004). The SMC and Calzetti’s laws
are commonly used in galaxy SED fitting methodology (e.g., Carnall et al.,
2018; Boquien et al., 2019). These models suggest that the extinction level
nearly monotonically increases to shorter wavelengths. The extinction curve
for the Milky Way interstellar dust also shows a similar behavior with a bump
at $\lambda\sim 2175~{}{\rm\AA}$, which is considered to originate from small
graphite grains and/or polycyclic aromatic hydrocarbon nanoparticles (Draine &
Lee, 1984). In contrast, the extinction law measured in the Orion Nebula, one
of the brightest nebulae powered by the central massive star forming region,
exhibits a curve that is significantly flattened in the far-UV wavelength
range and shows the absence of the $2175~{}{\rm\AA}$ bump. The relatively grey
extinction curve can be explained by the deficit of small-size dust grains.
Importantly, similar extinction laws are also observed in those of AGNs
(Gaskell et al. 2004; see also Maiolino et al. 2001a, Xie & Ho 2022), which
might preferentially remove small grains with sizes of $a\lesssim
0.1~{}\mu{\rm m}$. Additionally, JWST NIRSpec observations indicate a
flattening trend of dust extinction in star-forming galaxies at $6<z<11$,
potentially due to the formation of large-size dust grains within energetic
supernova ejecta (Markov et al., 2024).

Figure 2: The flux densities of dust-reddened AGNs normalized at
$4000~{}{\rm\AA}$ (solid curves), where the incident quasar SED is adopted
from a SED template model with a removal of host galaxy and dust emission
contributions (dashed curves, Temple et al., 2021). The six different
extinction laws shown in Figure 1 are applied by setting a visual extinction
of $A_{V}=3$ mag. The Markov’s extinction curve is not used because of the
lack of full wavelength coverage. The SED shaped with the Orion Nebula
extinction law reproduces the characteristic v-shape SED of LRDs. The averaged
photometric data of LRDs at $z\sim 5.5$ and $7.5$ are overlaid for comparison
(Barro et al., 2023).

In Figure 2, we show the attenuated flux using the five extinction laws. We
here adopt the composite SED of low-$z$ quasars from Temple et al. (2021) as
the incident flux. This model provides empirical SED templates for unobscured
quasars with $E(B-V)<0.3$ by using a total of $\sim 1.5\times 10^{4}$ SDSS
quasars with redshifts $0<z<5$. The SED covers from the rest-frame UV to NIR
up to 3 $\mu{\rm m}$. In these templates, the contribution from starlight from
the host galaxy is quantified using large quasar samples and subtracted to
isolate the pure quasar SED. For this work, we generate the SED template at
$z\simeq 5$, which excludes both the host galaxy and dust contribution. In
contrast, the widely-used quasar SED template provided by Vanden Berk et al.
(2001) is known to retain significant galaxy contributions, leading to a
substantially redder spectrum at $\lambda\gtrsim 5000~{}{\rm\AA}$ compared to
those in Temple et al. (2021). Although this discrepancy between the two SED
models does not impact our analysis of the rest-frame UV-to-optical
characteristics of LRDs, it introduces a serious inconsistency with rest-frame
NIR observations due to inadequate subtraction of host galaxy contributions
(e.g., MIRI photometry; Pérez-González et al. 2024). Therefore, for the
purpose of accurately modeling LRD SEDs in the absence of host galaxy effects,
we intentionally use the model from Temple et al. (2021) as the incident flux.

For all the cases, we assume a visual extinction $A_{\mathrm{V}}=3$ mag and
normalize the attenuated fluxes at $\lambda=4000~{}{\rm\AA}$. As reference, we
overlay photometric SED data obtained from stacking LRD samples at $z\sim 5.5$
(orange) and $z=7.5$ (blue) (Barro et al., 2023). The attenuated flux density
reproduces a red continuum component at $\lambda\gtrsim 3000-4000~{}{\rm\AA}$
with a slope of $\beta_{\rm opt}\gtrsim 0$, where the slope is defined as
$\beta\equiv{\rm d}\ln f_{\lambda}/{\rm d}\ln\lambda$. However, the SED shape
significantly differs in the rest-UV part, depending on the assumed extinction
curves. Models with extinction laws observed in the Orion Nebula and AGN
composite spectra flatten the SED with a UV slope of $\beta_{\rm UV}\sim-2$.
In contrast, other models result in strongly attenuated flux that fails to
maintain a UV excess. In the two models, the incident UV flux is reduced by
two orders of magnitude, i.e., $1\%$ of the incident flux.

The application of the extinction law derived from AGN composite spectra
(Gaskell et al., 2004) was tested on a LRD (ID: J0647-1045) identified at
$z=4.5321$ by JWST NIRSpec (Killi et al., 2023). For this particular object,
the inferred extinction curve matches that proposed by Gaskell et al. (2004).
Our findings in this work further demonstrate that such gray extinction curves
can explain the characteristic v-shaped SEDs of LRDs, as evidenced by
accumulated large samples, without considering multiple components such as
host galaxies or scattered AGN light.

Previous observational studies have proposed that the observed excess in rest-
frame UV ($\beta_{\rm UV}\sim-2$) may be attributed to the scattering of
unobscured quasar light (Kocevski et al., 2023; Labbe et al., 2023; Greene et
al., 2024). According to this scenario, the fraction of the scattered flux
relative to the primary component heavily depends on both the covering factor
of the scattering medium and the viewing angle, estimated to be about 1-3% to
account for the observed UV characteristics of LRDs. However, this scenario
does not account for the scattering fraction uniformly observed across most
LRDs. Relying on such a coincidental process to explain the universal
characteristics of an abundant AGN population is challenging111 Kokubo &
Harikane (2024) recently reported the analysis of LRD images over a few ($\sim
2-4$) epochs and found no obvious variability features. This analysis
challenges scenarios involving AGNs, which likely show a certain level of
variability. One possible reason why AGN variabilities might be smeared out is
due to scattering at multiple physical scales, by which combining photons with
different path lengths would cancel or dilute variable signatures over time..
In contrast, our model, which applies the dust opacity law observed in nearby
or low-redshift AGNs, offers a more robust scenario for the LRD
characteristics (see also discussion in Section 4.3), and furthermore bring us
new insights on the AGN structure within parsec scales.

Figure 3: The observed SEDs of LRDs across various redshifts of $4\leq z\leq
7$, where the intrinsic AGN luminosity is set to $L_{\rm bol}=10^{46}~{}{\rm
erg~{}s}^{-1}$. We adopt the composite SED template from Temple et al. (2021),
and assume the Orion-like extinction law with $A_{\rm V}=3$ mag. The 5$\sigma$
point-source imaging depths of JWST survey programs are overlaid: COSMOS-Web
(gray, Casey et al., 2023), CEERS (purple, Finkelstein et al., 2023), and
JADES-Medium/Deep (blue solid and dashed, Eisenstein et al., 2023).

Figure 3 presents the SEDs of our LRD model at various redshifts from $z=4$ to
$z=7$. For these models, we assume a bolometric luminosity of $L_{\rm
bol}=10^{46}~{}{\rm erg~{}s}^{-1}$ as the incident radiation from the AGN, and
calculate the monochromatic luminosity at rest-frame $1450~{}{\rm\AA}$ using a
bolometric correction factor of $f_{\rm bol,UV}=4.4$ (Richards et al., 2006).
Each case employs the Orion-Nebula attenuation law with a visual extinction
$A_{\mathrm{V}}=3~{}{\rm mag}$, as shown by the red curve in Figure 2. For
comparison, we overlay the $5~{}\sigma$ point-source imaging depths for each
NIRCam broadband filter used in several JWST survey programs. The depths
achieved by the CEERS and JADES surveys are sufficient to detect the
characteristic SED shape of LRDs from the rest-frame UV to optical
wavelengths. This is consistent with the observational fact that LRDs
identified in these deep surveys exhibit bolometric luminosities similar to
our models, $L_{\rm bol}\sim 10^{45-46}~{}{\rm erg~{}s}^{-1}$ (e.g., Kokorev
et al., 2024). However, the COSMOS-Web survey, which covers wider areas ($\sim
0.5$ deg2, approximately $\sim 20$ times larger than other surveys with $\sim
100$ arcmin2) with relatively shallower depths (as shown in Figure 3),
requires higher bolometric luminosities for the UV component to be detectable
and to capture the characteristic SED shape. Our SED models satisfy the
selection criteria for LRDs; $\beta_{\rm UV}<-0.37$ and $\beta_{\rm opt}>0$
(Kocevski et al. 2024; see also Greene et al. 2024).

##  3 Reprocessed Infrared Radiation from Obscuring Dust

Next, we analyze the SED of LRDs in the IR band, focusing on the re-emission
from dust grains that obscure the AGN light from UV to optical wavelengths.
Since UV photons emitted from the AGN heat the obscuring material, the dust
temperature rises to a sublimation point of $T_{\rm sub}\simeq 1500~{}{\rm
K}$, and decreases with distance from the AGN. In scenarios where dust and gas
are centrally concentrated, as in dusty tori (e.g., Ramos Almeida & Ricci,
2017), the hot dust just behind the sublimation radius significantly
contributes to the dust-grain mass budget and thus preferentially re-emits
thermal heat into the NIR bands, $\lambda\simeq 2.9~{}\mu{\rm m}~{}(T_{\rm
dust}/10^{3}~{}{\rm K})^{-1}$. As a result, the NIR flux density can be
$\gtrsim 2$ dex higher than the attenuated AGN emission at the optical band,
resulting in an extremely red continuum from the optical to infrared regimes,
as observed in dust-obscured AGN at lower redshifts (e.g., Alonso-Herrero et
al., 2006; Dey et al., 2008; Ramos Almeida et al., 2011; Ichikawa et al.,
2015, 2017; Toba et al., 2017; Lyu & Rieke, 2022). However, a good fraction of
high-redshift LRDs observed with MIRI broadband filters (e.g., F770W and
F1800W) do not show such red color but rather exhibit a relatively flat SED
from the rest-frame optical to the NIR wavelengths (Williams et al., 2024;
Wang et al., 2024). This trend is also seen in stacking images of multiple
LRDs (Akins et al., 2024). Nevertheless, there remain LRDs that show extremely
red SEDs from NIR to MIR (Lyu et al., 2024; Pérez-González et al., 2024).

The shape of the IR SED, reprocessed by dust grains heated by an AGN, depends
on the radial distribution of the dust (e.g., Pier & Krolik, 1992; Nenkova et
al., 2008; Nikutta et al., 2021a). When the obscuring material does not form a
concentrated density structure or is dispersed by bipolar outflows, the peak
energy shifts to the MIR since the bulk of the re-emission comes from
relatively cooler dust located farther from the dust sublimation radius (e.g.,
Barvainis, 1987). NIR interferometry observations of nearby AGNs in recent
decades (e.g., Kishimoto et al., 2011a, b) have found that the NIR emitting
region appears to have a radius nearly twice as large as the expected
sublimation radius, as traced by the NIR dust reverberation mapping
observations (Suganuma et al., 2006; Koshida et al., 2014). This suggests that
the dust density profile has a slightly lower level of concentration,
producing a larger emissivity size than the sublimation radius (see also
Nikutta et al., 2021b). Additionally, MIR interferometry observations
resolving the central $\sim$pc scale also showed that the significant fraction
of MIR dust emission originates from the polar region, suggesting the presence
of bipolar dusty outflows (e.g., Jaffe et al., 2004; Hönig et al., 2012, 2013;
Hönig, 2019; Tristram et al., 2014). According to radiation transfer
calculations that explain the IR SED of AGNs, including dust emission towards
the polar region, the gas and dust density profiles with a low level of
concentration ($\rho\propto r^{-\gamma}$; $-0.5<\gamma<0.5$) can reproduce the
observed MIR bump along with a relatively flat NIR SED at $\lambda\simeq
3-5~{}\mu{\rm m}$ (Hönig & Kishimoto, 2017).

Figure 4: The flux densities of LRDs normalized at $4000~{}{\rm\AA}$ for
different values of $n_{0}$ and $\gamma$, where the incident SED is the same
as in Figure 2. The stacked photometric data of LRDs are taken from the JADES
survey (thick curve, Pérez-González et al., 2024) and the COSMOS-Web survey
(circles Akins et al., 2024), and the upper bounds of the mid- to far-IR flux
densities due to non-detection in SCUBA and ALMA bands are taken from Akins et
al. (2024). Left panel: dependence of the SED on the density $n_{0}$. The
highest density case of $n_{0}=10^{3}~{}{\rm cm}^{-3}$ results in overly
bright rest-frame NIR, violating the MIRI F770W constraints. Lower density
cases reproduce SEDs consistent with the observed data points, while the
lowest density of $n_{0}=1~{}{\rm cm}^{-3}$ yields FIR SED peak inconsistent
with the upper bounds from the ALMA observations. Right panel: Variation of
SEDs with different density slopes $\gamma$. As the density is less
concentrated ($\gamma<1$), the energy peak shifts from the NIR to the MIR
regime.

Motivated by these AGN IR observations and models, we focus not only on a
compact dust structure such as a dust torus, but rather a more extended dust
structure described by a radial density profile of $\rho=\rho_{0}(r/r_{\rm
in})^{-\gamma}$, where $\rho_{0}$ is the gas density at the sublimation radius
$r_{\rm in}$, and the power-law index ranges $0<\gamma<2$. With the underlying
density distribution, the temperature of dust grains, $T_{\rm dust}(r)$, can
be determined by balancing dust heating from irradiation with cooling through
thermal emission as

$\int\frac{L_{\nu}e^{-\tau_{\nu}(r)}}{4\pi r^{2}}\kappa_{\nu}^{\rm
abs}d\nu=4\pi\int\kappa_{\nu}^{\rm
abs}B_{\nu}\left[T_{\mathrm{dust}}(r)\right]d\nu,$ (1)

where $L_{\nu}$ is the luminosity density of the central AGN at a frequency of
$\nu$, $\kappa_{\nu}^{\rm abs}$ is the absorption opacity of dust,
$B_{\nu}(T_{\rm dust})$ is the Planck function with a dust temperature $T_{\rm
dust}(r)$ at a distance of $r$ from the center, and $\tau_{\mathrm{\nu}}$ is
the optical depth defined by

$\tau_{\nu}(r)=\int_{r_{\rm in}}^{r}\rho\kappa_{\nu}^{\rm ext}dr^{\prime},$
(2)

where $\kappa_{\nu}^{\rm ext}$ is the sum of the absorption and scattering
opacity.

With a dust temperature profile, the luminosity density of dust re-emission is
calculated by

$L_{\nu,\mathrm{dust}}=4\pi\Omega\int_{r_{\rm in}}^{r_{\rm
out}}\kappa_{\nu}^{\rm
abs}B_{\nu}\left[T_{\mathrm{dust}}(r)\right]\rho{r}^{2}dr,$ (3)

where $0<\Omega\leq 4\pi$ is the solid angle covered by dust and $r_{\rm out}$
is determined by a constraint on the column density so that the visual
extinction reaches $A_{\rm V}\simeq 3$ mag to reproduce the characteristic red
continuum spectra of LRDs in the rest-frame optical bands (see Section 2). In
our calculation, we assume a full-sky covering case (i.e.,
$C\equiv\Omega/4\pi=1$), while this value will be determined when the relative
cosmic abundance between LRDs and unobscured AGNs and the dependence of $C$ on
the AGN properties such as the Eddington ratio are better understood for LRDs
(e.g., Ichikawa et al., 2019; Ricci et al., 2017, 2023). Note that the
equation (3) can be applied in the optically thin limit of dust re-emission in
the dusty surrounding medium. This approximation is valid for a wide range of
parameter combinations of the density value ($n_{0}\equiv\rho_{0}/m_{\rm p}$)
and slope ($\gamma$) of the distribution we explore in our SED modeling for
LRDs (see Section 4.1).

In addition to a diluted density distribution, we consider the feedback effect
on the dusty medium surrounding the central AGN (e.g., Chan & Krolik, 2016;
Namekata & Umemura, 2016; Kudoh et al., 2023; Soliman & Hopkins, 2023).
Radiation pressure expels dust grains between the sublimation radius and the
dust photosphere to UV radiation into a thin shell that accumulates at
$r=r_{\rm ph}$, the radius of the photosphere where the optical depth measured
from the center reaches unity. Within the thin layer at $r=r_{\rm ph}$, the
optical depth transitions from 0 to 1, and the dust temperature varies
depending on the detailed density distribution. For simplicity, we assume a
constant temperature within this layer, corresponding to the inner edge of the
thin shell. This temperature can be determined using Equation (1) with
$r=r_{\rm ph}$ and $\tau_{\nu}(r)=0$, setting an upper bound on the
temperature and emission from the thin shell. In our model, feedback effects
are incorporated by substituting $r_{\rm in}$ with $r_{\rm ph}$ in Equation
(3). While detailed radiation hydrodynamic simulations can pinpoint the exact
location of the photosphere, we evaluate $r_{\rm ph}$ using the original
unperturbed density profile, because the optical depth value within a radius
of $r$ remains unchanged if mass is conserved. As a result of the feedback,
hot dust just behind the sublimation radius $r_{\rm in}$ is pushed to larger
radii and becomes cooler, shifting the NIR emission to MIR. Despite these
change within $r_{\rm ph}$, the temperature profile and emission of the region
outside $r_{\rm ph}$ remain unaffected.

Figure 4 presents multi-wavelength SEDs of a LRD from the rest-frame UV to far
infrared with various combinations of ($n_{0}$, $\gamma$). This wavelength
coverage corresponds to the observed wavelengths frame from JWST NIRCam/MIRI
to ALMA bands for sources at $z=6$. For comparison, we overlay the stacking
photometric data of observed LRDs reported by the JADES survey (thick curve,
Pérez-González et al., 2024) and the COSMOS-Web survey (circles, Akins et al.,
2024). Furthermore, the upper bounds of the mid- to far-IR flux densities due
to non-detection in SCUBA and ALMA bands are shown (triangles, Akins et al.
2024; see also Labbe et al. 2023).

In the left panel of Figure 4, we show the dependence of the SED shape on the
gas density ($1\leq n_{0}/{\rm cm}^{-3}\leq 10^{3}$) at the sublimation
radius, while the slope of its radial profile is fixed to $\gamma=0$. When the
density is as high as $n_{0}=10^{3}~{}{\rm cm}^{-3}$, the outer edge of the
obscuring material is located at $r_{\rm out}\gtrsim 30~{}r_{\rm in}$, where
the dust temperature ranges from $T_{\rm dust}\simeq 1500~{}{\rm K}$ to $\sim
200~{}{\rm K}$. Therefore, the reprocessed IR spectrum begins to rise from
$\lambda\simeq 1~{}\mu{\rm m}$ and exhibits a near- to mid-infrared bump that
covers $\lambda\simeq 2-15~{}\mu{\rm m}$. In this case, the bright NIR bump
makes the SED shape substantially redder from $\lambda\simeq 6000~{}{\rm\AA}$
to $\sim 2~{}\mu{\rm m}$, consistently observed in AGNs with hot dust tori
(Polletta et al., 2007). However, this extremely red SED color does not match
the JWST MIRI photometric data for LRDs, which suggest a modest red color or
even a flat spectral shape. We note that such a contradiction with the MIRI
data appears as long as $n_{0}\gtrsim 10^{3}~{}{\rm cm}^{-3}$.

Figure 5:  Opt-NIR slope $\beta_{\mathrm{opt/NIR}}$ as the function of density
power-law index $\gamma$ and gas density $n_{0}$. High density cases lead to
bright rest-frame NIR emission and thus a steeper slope. The solid curve shows
the boundary of $\beta_{\mathrm{opt/NIR}}=-1$, below which the NIR SED in our
model is consistent wit the requirement from JWST/MIRI observations for LRDs.
The dashed curves present constant $30\leq T_{\mathrm{out}}/{\rm K}\leq 120$
and the temperature of Cosmic Microwave Background (CMB) at $z=6$. Below
$T_{\rm out}\simeq 30~{}{\rm K}$, the SED becomes inconsistent with the non-
detection of LRDs with ALMA observations. The allowed parameter space is
approximated as $10^{2.9\gamma+0.83}\lesssim n_{0}/{\rm cm}^{-3}\lesssim
10^{3}$.

On the other hand, as the density $n_{0}$ decreases and the outer radius
expands to reach a sufficient attenuation level, the peak energy of the IR SED
decreases because relatively cooler dust grains dominate both the mass and IR
emissivity. The energy transfer from NIR to MIR significantly suppresses the
NIR bump and eases the mismatch between our SED model and the JWST MIRI
photometric data. For the cases with $n_{0}=10$ and $100~{}{\rm cm}^{-3}$, the
SED shapes from the optical band to $\lambda=2-3~{}\mu{\rm m}$ agree to the
observed spectra of LRDs (or their stacked photometry based on the F444W and
F770W filters for $z\simeq 6$ LRDs), and also the peak energy is shifted to
the FIR band at $\lambda\sim 20-50~{}\mu{\rm m}$. As the density decreases to
$n_{0}=1~{}{\rm cm}^{-3}$, the IR SED peaks at $\lambda\sim 100~{}\mu{\rm m}$
and exceeds the upper bounds obtained by ALMA observations.

We also observe that our SED model with $\gamma=0$ and $n_{0}<10^{3}~{}{\rm
cm}^{-3}$ aligns with the upper bounds of the F1800W photometry. It is
important to note that the F1800W upper limit is derived from stacking images
covering less than 20% of the LRD samples observed with the F770W MIRI filter
(Akins et al., 2024). Further investigations into the rest-frame MIR band,
such as through F1800W, are crucial for enhancing our understanding of the
nature of LRDs (see the PRIMER survey; Dunlop et al. 2021).

The right panel of Figure 4 presents cases for a centrally concentrated
density distribution with $\gamma=2$ and $\gamma=1$, alongside the case of
$\gamma=0$ and $n_{0}=100~{}{\rm cm}^{-3}$ (see the left panel). At the
sublimation radius, the densities are set to $n_{0}=4\times 10^{4}~{}{\rm
cm}^{-3}$ and $4\times 10^{3}~{}{\rm cm}^{-3}$ for $\gamma=2$ and $1$,
respectively, ensuring sufficient extinction ($A_{\mathrm{V}}\simeq 3$ mag).
When the density is centrally concentrated with $\gamma=1-2$, a significant
fraction of the obscuring material is efficiently heated by the AGN, and thus
the hot dust just behind the sublimation radius is responsible for the IR
emissivity. As a result, the spectrum appears redder in the NIR wavelengths,
similar to the high-density case shown in the left panel.

##  4 Discussion

###  4.1 Physical properties of gas obscuring LRDs

In this section, we explore the parameter regions $(n_{0},\gamma)$ that
satisfy the conditions required for capturing the characteristic SED shape of
LRDs. We consider the following three key conditions.

  1. 1.

Dust attenuation: To reproduce the rest-frame UV-to-optical SED of LRDs, the
dust attenuation is set to $A_{\rm V}\simeq 3$ mag. This level of attenuation
corresponds to a hydrogen column density of $N_{\rm H}\simeq 7.5\times
10^{22}~{}{\rm cm}^{-2}~{}(Z/0.1~{}Z_{\odot})^{-1}$ (Baldwin et al., 1991).

  2. 2.

Optical-NIR SED slope: The spectral slope between the rest-frame optical and
NIR regions needs to be consistent with the JWST NIRCam/MIRI photometry,
showing that the F444W$-$F770W colors for $z\simeq 5-6$ LRDs are $\simeq 0.5$
mag in stacking images. To quantify the optical-to-NIR slope, we analyze flux
densities at $\lambda=0.65$ and $1.2~{}\mu{\rm m}$ from our SED model. An
crucial consideration is the influence of ${\rm H}\alpha$ line emission, which
can cause an excess in the NIRCam filter-convolved flux density at
$\lambda\simeq 0.65~{}\mu{\rm m}$. For LRD samples identified by slitless
spectroscopy (Matthee et al., 2023), we adopt an average rest-frame equivalent
width (EW) of ${\rm EW}_{0,\rm H\alpha}=500~{}{\rm\AA}$. The prominent line
emission enhances the filter-convolved flux densities within a filter with a
bandwidth of $\Delta\lambda_{\rm filter}$ by a factor of $\ell\equiv[{\rm
EW}_{0,\rm H\alpha}(1+z)+\Delta\lambda_{\rm filter}]/\Delta\lambda_{\rm
filter}$. Assuming that the ${\rm H}\alpha$ line emission falls within the
F444W, this excess factor of $\ell\simeq 1.35$ applies for a $z=6$,
corresponding to an additional $0.33$ mag. Therefore, we define the optical-
to-NIR slope as

$\beta_{\rm opt/NIR}=\frac{\log_{10}[f_{\nu_{2}}/(\ell
f_{\nu_{1}})]}{\log_{10}(\lambda_{2}/\lambda_{1})}-2,$ (4)

where $\lambda_{1}=0.65~{}\mu{\rm m}$, $\lambda_{2}=1.2~{}\mu{\rm m}$ to
measure the continuum flux $f_{\nu}$, and $\nu_{1(2)}=c/\lambda_{1(2)}$. Here
we impose $\beta_{\rm opt/NIR}<-1.0$ for this study.

  3. 3.

ALMA non-detection: Given the non-detection of LRDs in ALMA observations, we
impose a lower limit on the dust energy that re-emits the AGN light. This
requires that the minimum dust temperature, which is analytically calculated
in our SED model, exceeds a temperature threshold of $T_{\rm min}=30~{}{\rm
K}$.

Figure 6: Observed SEDs of LRDs in our model assuming an AGN bolometric
luminosity of $L_{\mathrm{bol}}=5\times 10^{45}$ erg s-1 at $z=6$. The colored
curves presents variations in SEDs with gas density $n_{0}$, ranging from
$150~{}{\rm cm}^{-3}$ (dark color) to $10~{}{\rm cm}^{-3}$ (light color),
while maintaining a constant radial distribution index $\gamma=0$. The stacked
photometric data for LRDs and upper bounds talen from Akins et al. (2024) are
overlaid. For comparison, we also show the SED templates of AGN hot dust tori
(Torus model), Apr 220, and Mrk 231 provided by (light blue curves, Polletta
et al., 2007), for which the same bolometric luminosity is assumed. At the
bottom of the plot, we show the wavelength coverage of related telescopes.

In the analysis described in Section 3, we find that the SED peaks in the NIR
regime and fails to meet the slope criterion of $\beta_{\rm opt/NIR}<-1$ when
a centrally concentrated density profile with $\gamma>1$ is assumed (see the
right panel of Figure 4). We thus focus on the parameter region with
$\gamma\leq 1$ in the following discussion.

Figure 5 summarizes the characteristic properties of the SED across various
combinations of $n_{0}$ and $\gamma$. The color contours indicate the optical-
to-NIR slope as defined by Equation (4). Employing the criterion of
$\beta_{\rm opt/NIR}=-1$, the parameter space is divided into two regions. In
the high density region ($n_{0}\gtrsim 10^{3}~{}{\rm cm}^{-3}$), the NIR SED
is influenced by re-emission from hot dust grains heated by the AGN, resulting
in a steeper optical-to-NIR slope that is inconsistent with the observational
requirement. We also present curves showing a constant dust temperature
$T_{\rm out}$ at the outer edge of the obscuring material at $r=r_{\rm out}$.
In cases with lower density, the outer radius further extends to achieve
$A_{\mathrm{V}}\simeq 3$ mag and $T_{\rm out}$ also decreases. By imposing
$T_{\rm out}\geq T_{\rm min}=30~{}{\rm K}$, we ensure that the FIR part of the
SED aligns with non-detection of LRDs with ALMA.

In summary, the region between the two curves of $\beta_{\rm opt/NIR}=-1$ and
$T_{\rm out}=30~{}{\rm K}$ defines the allowed parameter space of
$(n_{0},\gamma)$. Here, the density range is approximated as

$10^{2.9\gamma+0.83}\lesssim n_{0}/({\rm cm}^{-3})\lesssim 10^{3},$ (5)

which is valid $\gamma<0.8$. We also note that the density values are the
lower side of the typical electron density of the low redshift narrow line
region (NLR) with $n_{\rm e}\sim 10-10^{6}~{}{\rm cm}^{-3}$ (e.g., Netzer,
1990; Nesvadba et al., 2006; Liu et al., 2013; Harrison et al., 2014).
Recently, Integral Field Unit (IFU) observations that spatially resolves NLRs
of low-redshift AGNs show that the electron density at the inner $\sim
100~{}{\rm pc}$ scales is found to be $n_{\rm e}\sim 10^{2}-10^{3}~{}{\rm
cm}^{-3}$, and decreases toward larger radii (Kakkad et al., 2018). For high-
redshift star forming galaxies at $z\sim 6$ observed by JWST/NIRSpec, the
electron density in the nebular region is estimated as $n_{\rm e}\gtrsim
10^{2}~{}{\rm cm}^{-3}$ (e.g. Isobe et al., 2023; Abdurro’uf et al., 2024).

Finally, Figure 6 shows the compiled SED models that agree to the MIRI (F770W
and F1800W) and ALMA constraints. We consider a LRD at $z=6$ with a bolometric
luminosity of $L_{\mathrm{bol}}=5\times 10^{45}$ erg s-1. The color code
presents variations of gas density $n_{0}$, while maintaining a constant
radial distribution index $\gamma=0$. For comparison, we also show the SED
templates of AGN hot dust tori (Torus model), Apr 220, and Mrk 231 provided by
(light blue curves Polletta et al., 2007), for which the same bolometric
luminosity is assumed. The SED shapes of these compared templates violate
observational constraints either at the NIR or FIR regimes.

With the median value of the density in the log space ($n_{0}\sim 100~{}{\rm
cm}^{-3}$), the dust mass enclosed within $r\leq r_{\rm out}$ is estimated as
$M_{\rm dust}\simeq 1.2~{}(9.6)\times 10^{4}~{}M_{\odot}$ for $\gamma=0$ and
$0.2$, respectively, assuming a dust-to-gas mass ratio of
$0.001(Z/0.1~{}Z_{\odot})$. These results agree well to the mean value of the
dust mass $\langle M_{\rm dust}\rangle\simeq 1.6^{+4.8}_{-0.9}\times
10^{4}~{}M_{\odot}$ estimated using NIRCam images of LRDs (Casey et al.,
2024). Adopting the dust-to-stellar mass ratio of $f_{\star,\rm dust}\sim
10^{-4}-10^{-3}$ in a galaxy younger than one billion years (e.g., Valiante et
al., 2009), the expected mass of stars that form these dust grains would be
the order of $\sim 10^{7-9}~{}M_{\odot}$. This estimate can be used in
comparison to other independent stellar-mass measurements, e.g., by image
decompositions and SED fitting.

###  4.2 Host galaxy contributions?

Previous observations of LRDs have provided limited information about their
host galaxies. While some LRDs exhibit extended morphology in short-wavelength
images, achieving a robust confirmation of host galaxy features remains
challenging. Detailed measurements of JWST/PRISM spectra have successfully
identified a continuum break feature at $\lambda\simeq 4000~{}{\rm\AA}$ in two
LRDs: MSAID13123 at $z=7.04$ in Greene et al. (2024) and RUBIES-BLAGN-1 at
$z=3.1$ in Wang et al. (2024). These findings suggest the presence of an
underlying host galaxy. However, if stellar light accounted for the continuum
flux across the entire optical band, the required stellar mass would approach
$M_{\star}\gtrsim 10^{11}~{}M_{\odot}$, comparable to that of the Milky Way.
When combining such a high mass and their high abundance in the cosmic volume,
this scenario implies an extraordinary stellar mass density that contradicts
the standard cosmological framework (see also Inayoshi & Ichikawa, 2024; Akins
et al., 2024). As a possible solution, Wang et al. (2024) suggest that the
stellar contribution on the SED is limited to the short-wavelength range at
$\lambda\lesssim 4000~{}{\rm\AA}$ and the AGN is responsible for the red
optical continuum at the longer wavelengths. Under this assumption, the SED
fitting indicates a minimal stellar mass of $M_{\star}\sim
10^{9}~{}M_{\odot}$, resolving the tension with the cosmological prediction.

Here, we emphasize that our approach does not rule out the existence of host
galaxies in LRDs, but demonstrates that the SED model focusing only on AGNs
successfully captures the characteristic spectral properties of LRDs: (1) the
v-shaped SED in rest-frame UV to optical bands and (2) modest NIR emission. In
our AGN-only spectral model, as a natural consequence, the red optical
continuum is attributed to the attenuated AGN light. This interpretation is
indeed consistent with observations that the observed continuum and ${\rm
H}\alpha$ line fluxes for LRDs follow the empirical correlation observed in
nearby AGNs (Greene & Ho, 2005). In other words, if stellar light dominates
the observed continuum flux of LRDs, the EW of broad ${\rm H}\alpha$ emission
line would be unusually high, compared to normal AGNs. This consistency is
crucial, as it strengthens the interpretation of BH mass measurements using
the single-epoch method.

###  4.3 Why do LRDs emerge at high redshifts of $z>4$
with a universal SED shape?

There is a critical question: Why do LRDs emerge at $z>4$ with a uniform
extinction level of $A_{\rm V}=3$ mag? There is no clear answer to this
question at this time, but some speculations are given below.

The key assumption in our model is the relaxation of density concentration
around the AGN. This idea implies that LRDs form in the early stages of galaxy
and AGN formation, characterized by dynamically evolving systems (e.g., clumpy
gaseous structures, irregularly shaped and misaligned galactic disks,
significant perturbations from mergers) before a well-organized dust torus
structure is established (e.g., Ceverino et al., 2010, 2015; Hopkins et al.,
2024). In the nascent stages of galaxy formation, the unsettled nature of the
gas distribution permits a variety of attenuation levels and column densities.
At higher redshifts, if the column density is typically higher and spans a
broader range, numerous objects will exhibit $A_{\rm V}\gg 1$. However, in the
context of our proposed model, only those with moderate extinction levels with
$A_{\rm V}\sim 1-3$ would reveal the characteristic SED shape, as seen in
observed LRD samples (e.g., Barro et al., 2023; Matthee et al., 2023; Greene
et al., 2024; Kocevski et al., 2024). With $A_{\rm V}\gtrsim 5$, which
corresponds to $A_{3000}=1.6~{}A_{\rm V}\gtrsim 8$ at
$\lambda=3000~{}{\rm\AA}$ with the Orion-Nebula extinction curve, it is
challenging to observe the rest-frame UV part of LRDs even with the JWST
sensitivity. In this case, these objects appear to have only a red continuum
without an UV excess, which are categorized as hot-dust-obscured AGNs observed
at the lower-redshift universe (Alonso-Herrero et al., 2006; Ramos Almeida et
al., 2011; Ichikawa et al., 2014; Hickox & Alexander, 2018; Lyu & Rieke,
2022). Our scenario seems more relevant in earlier epochs before a robust dust
torus structure forms in its galactic nucleus. When a centrally concentrated
dust torus is established, as seen in nearby AGNs, there is a dichotomy
between type 1 and type 2 AGN populations due to the geometric arrangement
(e.g., Ramos Almeida & Ricci, 2017; Hickox & Alexander, 2018) rather than a
modest variation due to $A_{\rm V}$.

Another critical element in our model is the flattening of dust extinction
curves at shorter wavelengths, $\lambda\lesssim 3000~{}{\rm\AA}$. This
characteristic may help explain why dust-reddened objects with v-shaped SEDs
predominantly appear at high redshifts. In the early stages of protogalaxies,
particularly those younger than one billion years ($z\gtrsim 6-7$), dust
production is primarily driven by massive stars that expel heavy elements into
the interstellar medium (ISM) through core-collapse supernovae (Todini &
Ferrara, 2001; Valiante et al., 2009). Due to effective dust destruction
mechanisms that target smaller grains, the dust-size distribution tends to
favor larger particles (Nozawa et al., 2006; Asano et al., 2013). These larger
grains, once released into the ISM, undergo shattering processes and break
into smaller fragments in a typical timescale of $\gtrsim 0.5-1$ Gyr. As a
result, the extinction curve at shorter wavelengths begins to evolve towards
those observed in the ISM of the Milky Way. The temporal evolution in dust
characteristics would be critical for shaping the unique spectral features of
LRDs only at high redshifts.

###  4.4 Implications to the cosmological BH evolution

In this paper, we propose that LRDs originate from rapidly growing BH embedded
within dusty flows with an extended density distribution. In this framework,
we reach two intriguing conclusions regarding the cosmological BH evolution.

First, the luminosity of LRDs is primarily powered by AGNs that host accreting
BHs. These LRDs are more abundant than X-ray selected AGNs, suggesting an
increase in the cosmic growth rate of BHs beyond $z\simeq 6$ (Kokorev et al.,
2024; Kocevski et al., 2024; Akins et al., 2024). The BH accretion rate
density derived from the bolometric luminosity functions of LRDs is
significantly higher than that estimated from other AGN surveys. To reconcile
this high value with the BH mass density at $z\sim 5$, a radiative efficiency
of $\gtrsim 20\%$ – twice the canonical 10% value – is required with a
$3\sigma$ confidence level (Inayoshi & Ichikawa, 2024). This high radiative
efficiency implies that the BHs are rapidly spinning in the LRD stage through
violent mass accretion during the early galaxy assembly.

Second, when the extinction level is accurately determined from the reddened
optical continuum shape, the AGN continuum and broad-line luminosity can be
more precisely estimated. These luminosity values are crucial for measuring BH
mass using the single-epoch method. For instance, Matthee et al. (2023)
calculated BH masses for LRDs using the ${\rm H}\alpha$ luminosity as an
indicator of AGN continuum luminosity, but without accounting for dust
reddening correction. According to the empirical relationship in the single-
epoch method, we derive that the Eddington ratio follows $\lambda_{\rm
Edd}\propto L_{\rm H\alpha}^{0.394}$ (see Greene & Ho, 2005; Matthee et al.,
2023). Therefore, by applying a dust correction, the Eddington ratio increases
by a factor of $\mathscr{L}\sim 10^{0.14A_{\rm V}}\sim 2.6$ for $A_{\rm V}=3$.
With this correction, the Eddington ratios of BHs reported by Matthee et al.
(2023) approach unity, suggesting that a significant fraction of LRDs may be
super-Eddington accretors. This result highlights the importance of accurate
dust correction for high-redshift AGNs, particularly LRDs, to understand the
true nature of their accretion dynamics (Lupi et al. 2024; see also Du et al.
2014).

##  5 Conclusions

JWST observations have uncovered the presence of very compact and red objects
at the high-redshift universe, referred to as little red dots (LRDs).
Spectroscopic observations have confirmed the broad components of emission
lines and the presence of AGNs within LRDs. Despite the breakthrough
discovery, the characteristic v-shape spectral feature observed through
JWST/NIRCam challenges us to understand the origin of these LRDs and the
contribution from AGNs powered by accreting BHs.

In this work, we propose a SED model for LRDs spanning from rest-frame UV to
infrared bands. We hypothesize that the incident radiation from an AGN,
characterized by a typical SED, is embedded in an extended dusty medium, which
has an extinction law similar to those seen in dense regions such as Orion
Nebula or certain AGN environments. The UV-optical spectrum is described by
dust-attenuated AGN emission with $A_{\rm V}\simeq 3$ mag, featuring a red
optical continuum at $\lambda>4000~{}{\rm\AA}$ and a flat UV spectral shape
established through a gray extinction curve at $\lambda<3000~{}{\rm\AA}$, due
to the absence of small-size grains. In the infrared, the SED is shaped by an
extended dust and gas distribution ($\gamma<1$; $\rho\propto r^{-\gamma}$)
with a characteristic gas density of $\simeq 10-10^{3}~{}{\rm cm}^{-3}$ at the
dust sublimation radius, which allows relatively cool dust temperatures to
dominate the radiation, thereby shifting the energy peak from near- to mid-
infrared bands. This model, unlike the typical AGN hot torus models, can
produce a infrared SED flattening that is consistent with LRD observations
through JWST MIRI. In this scenario, there is no need for additional stellar
emission or AGN scattered light.

We greatly thank Changhao Chen, Seiji Fujimoto, Jinyi Shangguan, and Mingyang
Zhuang for constructive discussions. We acknowledge support from the National
Natural Science Foundation of China (12073003, 12003003, 11721303, 11991052,
11950410493), and the China Manned Space Project (CMS-CSST-2021-A04 and CMS-
CSST-2021-A06). This work is also supported by Japan Society for the Promotion
of Science (JSPS) KAKENHI (20H01939; K. Ichikawa). LCH was supported by the
National Science Foundation of China (11991052, 12233001), the National Key
R&D Program of China (2022YFF0503401), and the China Manned Space Project
(CMS-CSST-2021-A04, CMS-CSST-2021-A06). This work made use of
Astropy:222http://www.astropy.org a community-developed core Python package
and an ecosystem of tools and resources for astronomy (Astropy Collaboration
et al., 2022). This work use the dust extinction data from version 23.01 of
CLOUDY, last described by Gunasekera et al. (2023).

## References

  * Abdurro’uf et al. (2024) Abdurro’uf, Larson, R. L., Coe, D., et al. 2024, arXiv e-prints, arXiv:2404.16201, doi: 10.48550/arXiv.2404.16201
  * Akins et al. (2024) Akins, H. B., Casey, C. M., Lambrides, E., et al. 2024, arXiv e-prints, arXiv:2406.10341, doi: 10.48550/arXiv.2406.10341
  * Alonso-Herrero et al. (2006) Alonso-Herrero, A., Pérez-González, P. G., Alexander, D. M., et al. 2006, ApJ, 640, 167, doi: 10.1086/499800
  * Asano et al. (2013) Asano, R. S., Takeuchi, T. T., Hirashita, H., & Inoue, A. K. 2013, Earth, Planets and Space, 65, 213, doi: 10.5047/eps.2012.04.014
  * Astropy Collaboration et al. (2022) Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, ApJ, 935, 167, doi: 10.3847/1538-4357/ac7c74
  * Baldwin et al. (1991) Baldwin, J. A., Ferland, G. J., Martin, P. G., et al. 1991, ApJ, 374, 580, doi: 10.1086/170146
  * Barro et al. (2023) Barro, G., Perez-Gonzalez, P. G., Kocevski, D. D., et al. 2023, arXiv e-prints, arXiv:2305.14418, doi: 10.48550/arXiv.2305.14418
  * Barvainis (1987) Barvainis, R. 1987, ApJ, 320, 537, doi: 10.1086/165571
  * Boquien et al. (2019) Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A, 622, A103, doi: 10.1051/0004-6361/201834156
  * Calzetti et al. (2000) Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682, doi: 10.1086/308692
  * Cardelli et al. (1989) Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1989, ApJ, 345, 245, doi: 10.1086/167900
  * Carnall et al. (2018) Carnall, A. C., McLure, R. J., Dunlop, J. S., & Davé, R. 2018, MNRAS, 480, 4379, doi: 10.1093/mnras/sty2169
  * Casey et al. (2024) Casey, C. M., Akins, H. B., Kokorev, V., et al. 2024, arXiv e-prints, arXiv:2407.05094, doi: 10.48550/arXiv.2407.05094
  * Casey et al. (2023) Casey, C. M., Kartaltepe, J. S., Drakos, N. E., et al. 2023, ApJ, 954, 31, doi: 10.3847/1538-4357/acc2bc
  * Ceverino et al. (2010) Ceverino, D., Dekel, A., & Bournaud, F. 2010, MNRAS, 404, 2151, doi: 10.1111/j.1365-2966.2010.16433.x
  * Ceverino et al. (2015) Ceverino, D., Dekel, A., Tweed, D., & Primack, J. 2015, MNRAS, 447, 3291, doi: 10.1093/mnras/stu2694
  * Chan & Krolik (2016) Chan, C.-H., & Krolik, J. H. 2016, ApJ, 825, 67, doi: 10.3847/0004-637X/825/1/67
  * Dayal et al. (2024) Dayal, P., Volonteri, M., Greene, J. E., et al. 2024, arXiv e-prints, arXiv:2401.11242, doi: 10.48550/arXiv.2401.11242
  * Dey et al. (2008) Dey, A., Soifer, B. T., Desai, V., et al. 2008, ApJ, 677, 943, doi: 10.1086/529516
  * Draine & Lee (1984) Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89, doi: 10.1086/162480
  * Draine & Salpeter (1979) Draine, B. T., & Salpeter, E. E. 1979, ApJ, 231, 77, doi: 10.1086/157165
  * Du et al. (2014) Du, P., Hu, C., Lu, K.-X., et al. 2014, ApJ, 782, 45, doi: 10.1088/0004-637X/782/1/45
  * Dunlop et al. (2021) Dunlop, J. S., Abraham, R. G., Ashby, M. L. N., et al. 2021, PRIMER: Public Release IMaging for Extragalactic Research, JWST Proposal. Cycle 1, ID. #1837
  * Eisenstein et al. (2023) Eisenstein, D. J., Willott, C., Alberts, S., et al. 2023, arXiv e-prints, arXiv:2306.02465, doi: 10.48550/arXiv.2306.02465
  * Finkelstein et al. (2023) Finkelstein, S. L., Bagley, M. B., Ferguson, H. C., et al. 2023, ApJ, 946, L13, doi: 10.3847/2041-8213/acade4
  * Furtak et al. (2023) Furtak, L. J., Zitrin, A., Plat, A., et al. 2023, ApJ, 952, 142, doi: 10.3847/1538-4357/acdc9d
  * Gaskell et al. (2004) Gaskell, C. M., Goosmann, R. W., Antonucci, R. R. J., & Whysong, D. H. 2004, ApJ, 616, 147, doi: 10.1086/423885
  * Gordon et al. (2003) Gordon, K. D., Clayton, G. C., Misselt, K. A., Landolt, A. U., & Wolff, M. J. 2003, ApJ, 594, 279, doi: 10.1086/376774
  * Greene & Ho (2005) Greene, J. E., & Ho, L. C. 2005, ApJ, 630, 122, doi: 10.1086/431897
  * Greene et al. (2024) Greene, J. E., Labbe, I., Goulding, A. D., et al. 2024, ApJ, 964, 39, doi: 10.3847/1538-4357/ad1e5f
  * Gunasekera et al. (2023) Gunasekera, C. M., van Hoof, P. A. M., Chatzikos, M., & Ferland, G. J. 2023, Research Notes of the American Astronomical Society, 7, 246, doi: 10.3847/2515-5172/ad0e75
  * Harrison et al. (2014) Harrison, C. M., Alexander, D. M., Mullaney, J. R., & Swinbank, A. M. 2014, MNRAS, 441, 3306, doi: 10.1093/mnras/stu515
  * He et al. (2023) He, W., Akiyama, M., Enoki, M., et al. 2023, arXiv e-prints, arXiv:2311.08922, doi: 10.48550/arXiv.2311.08922
  * Hickox & Alexander (2018) Hickox, R. C., & Alexander, D. M. 2018, ARA&A, 56, 625, doi: 10.1146/annurev-astro-081817-051803
  * Hönig (2019) Hönig, S. F. 2019, ApJ, 884, 171, doi: 10.3847/1538-4357/ab4591
  * Hönig & Kishimoto (2017) Hönig, S. F., & Kishimoto, M. 2017, ApJ, 838, L20, doi: 10.3847/2041-8213/aa6838
  * Hönig et al. (2012) Hönig, S. F., Kishimoto, M., Antonucci, R., et al. 2012, ApJ, 755, 149, doi: 10.1088/0004-637X/755/2/149
  * Hönig et al. (2013) Hönig, S. F., Kishimoto, M., Tristram, K. R. W., et al. 2013, ApJ, 771, 87, doi: 10.1088/0004-637X/771/2/87
  * Hopkins et al. (2024) Hopkins, P. F., Grudic, M. Y., Su, K.-Y., et al. 2024, The Open Journal of Astrophysics, 7, 18, doi: 10.21105/astro.2309.13115
  * Hu et al. (2022) Hu, H., Inayoshi, K., Haiman, Z., et al. 2022, ApJ, 935, 140, doi: 10.3847/1538-4357/ac7daa
  * Ichikawa et al. (2014) Ichikawa, K., Imanishi, M., Ueda, Y., et al. 2014, ApJ, 794, 139, doi: 10.1088/0004-637X/794/2/139
  * Ichikawa et al. (2017) Ichikawa, K., Ricci, C., Ueda, Y., et al. 2017, ApJ, 835, 74, doi: 10.3847/1538-4357/835/1/74
  * Ichikawa et al. (2015) Ichikawa, K., Packham, C., Ramos Almeida, C., et al. 2015, ApJ, 803, 57, doi: 10.1088/0004-637X/803/2/57
  * Ichikawa et al. (2019) Ichikawa, K., Ricci, C., Ueda, Y., et al. 2019, ApJ, 870, 31, doi: 10.3847/1538-4357/aaef8f
  * Inayoshi & Ichikawa (2024) Inayoshi, K., & Ichikawa, K. 2024, arXiv e-prints, arXiv:2402.14706, doi: 10.48550/arXiv.2402.14706
  * Inayoshi et al. (2022) Inayoshi, K., Nakatani, R., Toyouchi, D., et al. 2022, ApJ, 927, 237, doi: 10.3847/1538-4357/ac4751
  * Isobe et al. (2023) Isobe, Y., Ouchi, M., Nakajima, K., et al. 2023, ApJ, 956, 139, doi: 10.3847/1538-4357/acf376
  * Jaffe et al. (2004) Jaffe, W., Meisenheimer, K., Röttgering, H. J. A., et al. 2004, Nature, 429, 47, doi: 10.1038/nature02531
  * Kakkad et al. (2018) Kakkad, D., Groves, B., Dopita, M., et al. 2018, A&A, 618, A6, doi: 10.1051/0004-6361/201832790
  * Killi et al. (2023) Killi, M., Watson, D., Brammer, G., et al. 2023, arXiv e-prints, arXiv:2312.03065, doi: 10.48550/arXiv.2312.03065
  * Kishimoto et al. (2011a) Kishimoto, M., Hönig, S. F., Antonucci, R., et al. 2011a, A&A, 527, A121, doi: 10.1051/0004-6361/201016054
  * Kishimoto et al. (2011b) —. 2011b, A&A, 536, A78, doi: 10.1051/0004-6361/201117367
  * Kocevski et al. (2023) Kocevski, D. D., Onoue, M., Inayoshi, K., et al. 2023, arXiv e-prints, arXiv:2302.00012.  https://arxiv.org/abs/2302.00012
  * Kocevski et al. (2024) Kocevski, D. D., Finkelstein, S. L., Barro, G., et al. 2024, arXiv e-prints, arXiv:2404.03576, doi: 10.48550/arXiv.2404.03576
  * Kokorev et al. (2023) Kokorev, V., Fujimoto, S., Labbe, I., et al. 2023, ApJ, 957, L7, doi: 10.3847/2041-8213/ad037a
  * Kokorev et al. (2024) Kokorev, V., Caputi, K. I., Greene, J. E., et al. 2024, arXiv e-prints, arXiv:2401.09981, doi: 10.48550/arXiv.2401.09981
  * Kokubo & Harikane (2024) Kokubo, M., & Harikane, Y. 2024, arXiv e-prints, arXiv:2407.04777, doi: 10.48550/arXiv.2407.04777
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511, doi: 10.1146/annurev-astro-082708-101811
  * Koshida et al. (2014) Koshida, S., Minezaki, T., Yoshii, Y., et al. 2014, ApJ, 788, 159, doi: 10.1088/0004-637X/788/2/159
  * Kudoh et al. (2023) Kudoh, Y., Wada, K., Kawakatu, N., & Nomura, M. 2023, ApJ, 950, 72, doi: 10.3847/1538-4357/accc2b
  * Labbe et al. (2023) Labbe, I., Greene, J. E., Bezanson, R., et al. 2023, arXiv e-prints, arXiv:2306.07320, doi: 10.48550/arXiv.2306.07320
  * Laor & Draine (1993) Laor, A., & Draine, B. T. 1993, ApJ, 402, 441, doi: 10.1086/172149
  * Li (2008) Li, A. 2008, arXiv e-prints, arXiv:0808.4117, doi: 10.48550/arXiv.0808.4117
  * Li et al. (2023) Li, W., Inayoshi, K., Onoue, M., et al. 2023, arXiv e-prints, arXiv:2306.06172, doi: 10.48550/arXiv.2306.06172
  * Liu et al. (2013) Liu, G., Zakamska, N. L., Greene, J. E., Nesvadba, N. P. H., & Liu, X. 2013, MNRAS, 436, 2576, doi: 10.1093/mnras/stt1755
  * Lupi et al. (2024) Lupi, A., Trinca, A., Volonteri, M., Dotti, M., & Mazzucchelli, C. 2024, arXiv e-prints, arXiv:2406.17847, doi: 10.48550/arXiv.2406.17847
  * Lyu & Rieke (2022) Lyu, J., & Rieke, G. 2022, Universe, 8, 304, doi: 10.3390/universe8060304
  * Lyu et al. (2024) Lyu, J., Alberts, S., Rieke, G. H., et al. 2024, ApJ, 966, 229, doi: 10.3847/1538-4357/ad3643
  * Madau et al. (2024) Madau, P., Giallongo, E., Grazian, A., & Haardt, F. 2024, arXiv e-prints, arXiv:2406.18697, doi: 10.48550/arXiv.2406.18697
  * Maiolino et al. (2001a) Maiolino, R., Marconi, A., & Oliva, E. 2001a, A&A, 365, 37, doi: 10.1051/0004-6361:20000012
  * Maiolino et al. (2001b) Maiolino, R., Marconi, A., Salvati, M., et al. 2001b, A&A, 365, 28, doi: 10.1051/0004-6361:20000177
  * Maiolino et al. (2004) Maiolino, R., Schneider, R., Oliva, E., et al. 2004, Nature, 431, 533, doi: 10.1038/nature02930
  * Markov et al. (2024) Markov, V., Gallerani, S., Ferrara, A., et al. 2024, arXiv e-prints, arXiv:2402.05996, doi: 10.48550/arXiv.2402.05996
  * Matsuoka et al. (2018) Matsuoka, Y., Strauss, M. A., Kashikawa, N., et al. 2018, ApJ, 869, 150, doi: 10.3847/1538-4357/aaee7a
  * Matsuoka et al. (2023) Matsuoka, Y., Onoue, M., Iwasawa, K., et al. 2023, arXiv e-prints, arXiv:2305.11225, doi: 10.48550/arXiv.2305.11225
  * Matthee et al. (2023) Matthee, J., Naidu, R. P., Brammer, G., et al. 2023, arXiv e-prints, arXiv:2306.05448, doi: 10.48550/arXiv.2306.05448
  * Namekata & Umemura (2016) Namekata, D., & Umemura, M. 2016, MNRAS, 460, 980, doi: 10.1093/mnras/stw862
  * Nenkova et al. (2008) Nenkova, M., Sirocky, M. M., Ivezić, Ž., & Elitzur, M. 2008, ApJ, 685, 147, doi: 10.1086/590482
  * Nesvadba et al. (2006) Nesvadba, N. P. H., Lehnert, M. D., Eisenhauer, F., et al. 2006, ApJ, 650, 693, doi: 10.1086/507266
  * Netzer (1990) Netzer, H. 1990, in Active Galactic Nuclei, ed. R. D. Blandford, H. Netzer, L. Woltjer, T. J. L. Courvoisier, & M. Mayor, 57–160
  * Niida et al. (2020) Niida, M., Nagao, T., Ikeda, H., et al. 2020, ApJ, 904, 89, doi: 10.3847/1538-4357/abbe11
  * Nikutta et al. (2021a) Nikutta, R., Lopez-Rodriguez, E., Ichikawa, K., et al. 2021a, ApJ, 919, 136, doi: 10.3847/1538-4357/ac06a6
  * Nikutta et al. (2021b) —. 2021b, ApJ, 923, 127, doi: 10.3847/1538-4357/ac2949
  * Nozawa et al. (2006) Nozawa, T., Kozasa, T., & Habe, A. 2006, ApJ, 648, 435, doi: 10.1086/505639
  * Onoue et al. (2023) Onoue, M., Inayoshi, K., Ding, X., et al. 2023, ApJ, 942, L17, doi: 10.3847/2041-8213/aca9d3
  * Pacucci et al. (2023) Pacucci, F., Nguyen, B., Carniani, S., Maiolino, R., & Fan, X. 2023, ApJ, 957, L3, doi: 10.3847/2041-8213/ad0158
  * Pérez-González et al. (2024) Pérez-González, P. G., Barro, G., Rieke, G. H., et al. 2024, arXiv e-prints, arXiv:2401.08782, doi: 10.48550/arXiv.2401.08782
  * Pier & Krolik (1992) Pier, E. A., & Krolik, J. H. 1992, ApJ, 401, 99, doi: 10.1086/172042
  * Polletta et al. (2007) Polletta, M., Tajer, M., Maraschi, L., et al. 2007, ApJ, 663, 81, doi: 10.1086/518113
  * Ramos Almeida & Ricci (2017) Ramos Almeida, C., & Ricci, C. 2017, Nature Astronomy, 1, 679, doi: 10.1038/s41550-017-0232-z
  * Ramos Almeida et al. (2011) Ramos Almeida, C., Levenson, N. A., Alonso-Herrero, A., et al. 2011, ApJ, 731, 92, doi: 10.1088/0004-637X/731/2/92
  * Ricci et al. (2017) Ricci, C., Trakhtenbrot, B., Koss, M. J., et al. 2017, Nature, 549, 488, doi: 10.1038/nature23906
  * Ricci et al. (2023) Ricci, C., Ichikawa, K., Stalevski, M., et al. 2023, ApJ, 959, 27, doi: 10.3847/1538-4357/ad0733
  * Richards et al. (2006) Richards, G. T., Strauss, M. A., Fan, X., et al. 2006, AJ, 131, 2766, doi: 10.1086/503559
  * Scoggins et al. (2023) Scoggins, M. T., Haiman, Z., & Wise, J. H. 2023, MNRAS, 519, 2155, doi: 10.1093/mnras/stac3715
  * Soliman & Hopkins (2023) Soliman, N. H., & Hopkins, P. F. 2023, MNRAS, 525, 2668, doi: 10.1093/mnras/stad2460
  * Suganuma et al. (2006) Suganuma, M., Yoshii, Y., Kobayashi, Y., et al. 2006, ApJ, 639, 46, doi: 10.1086/499326
  * Tazaki & Ichikawa (2020) Tazaki, R., & Ichikawa, K. 2020, ApJ, 892, 149, doi: 10.3847/1538-4357/ab72f6
  * Tazaki et al. (2020) Tazaki, R., Ichikawa, K., & Kokubo, M. 2020, ApJ, 892, 84, doi: 10.3847/1538-4357/ab7822
  * Temple et al. (2021) Temple, M. J., Hewett, P. C., & Banerji, M. 2021, MNRAS, 508, 737, doi: 10.1093/mnras/stab2586
  * Toba et al. (2017) Toba, Y., Bae, H.-J., Nagao, T., et al. 2017, ApJ, 850, 140, doi: 10.3847/1538-4357/aa918a
  * Todini & Ferrara (2001) Todini, P., & Ferrara, A. 2001, MNRAS, 325, 726, doi: 10.1046/j.1365-8711.2001.04486.x
  * Tristram et al. (2014) Tristram, K. R. W., Burtscher, L., Jaffe, W., et al. 2014, A&A, 563, A82, doi: 10.1051/0004-6361/201322698
  * Valiante et al. (2009) Valiante, R., Schneider, R., Bianchi, S., & Andersen, A. C. 2009, MNRAS, 397, 1661, doi: 10.1111/j.1365-2966.2009.15076.x
  * Vanden Berk et al. (2001) Vanden Berk, D. E., Richards, G. T., Bauer, A., et al. 2001, AJ, 122, 549, doi: 10.1086/321167
  * Wang et al. (2024) Wang, B., de Graaff, A., Davies, R. L., et al. 2024, arXiv e-prints, arXiv:2403.02304, doi: 10.48550/arXiv.2403.02304
  * Williams et al. (2024) Williams, C. C., Alberts, S., Ji, Z., et al. 2024, ApJ, 968, 34, doi: 10.3847/1538-4357/ad3f17
  * Xie & Ho (2022) Xie, Y., & Ho, L. C. 2022, ApJ, 925, 218, doi: 10.3847/1538-4357/ac32e2
  * Xie et al. (2017) Xie, Y., Li, A., & Hao, L. 2017, ApJS, 228, 6, doi: 10.3847/1538-4365/228/1/6