# Revisiting the dust torus size – luminosity relation based on a uniform
reverberation mapping analysis

Amit Kumar Mandal Department of Physics $\&$ Astronomy, Seoul National
University
Seoul 08826, Republic of Korea<EMAIL_ADDRESS>Jong-Hak Woo Department of
Physics $\&$ Astronomy, Seoul National University
Seoul 08826, Republic of Korea<EMAIL_ADDRESS>Shu Wang Department of
Physics $\&$ Astronomy, Seoul National University
Seoul 08826, Republic of Korea<EMAIL_ADDRESS>Suvendu Rakshit Aryabhatta
Research Institute of Observational Sciences
Manora Peak, Nainital-263001, Uttarakhand, India  Hojin Cho Department of
Physics $\&$ Astronomy, Seoul National University
Seoul 08826, Republic of Korea<EMAIL_ADDRESS>Donghoon Son Department of
Physics $\&$ Astronomy, Seoul National University
Seoul 08826, Republic of Korea<EMAIL_ADDRESS>C. S. Stalin Indian Institute
of Astrophysics
Block II, Koramangala, Bangalore, 560 034, India

###### Abstract

We investigate the torus size – luminosity relation of Type 1 AGNs based on
the reverberation-mapping analysis using the light curves of the optical
continuum and the IR continuum obtained with the W1 and W2- bands of the Wide-
field Infrared Survey Explorer (WISE) survey. The final sample consists of 446
and 416 AGNs, respectively, for W1 and W2- band light curves, covering a large
dynamic range of bolometric luminosity from $10^{43.4}$ to $10^{47.6}$
$erg\,s^{-1}$, which show reliable lag measurements based on our quality
assessment analysis. After correcting for the accretion disk contamination in
the observed IR flux, we constrain the torus size ($R_{dust}$) and AGN
bolometric luminosity ($L_{bol}$) relationship with the best-fit slope of 0.39
(0.33) for the W1- (W2) band, which is shallower than expected from the dust
radiation equilibrium model. By combining the previous K- band lag
measurements, we find that the measured torus size depends on the observed
wavelength of the dust radiation, as $R_{dust,K}:R_{dust,W1}:R_{dust,W2}$ =
1.0:1.5:1.8 ($R_{dust}\,\propto\,\lambda^{0.80}$) at $L_{bol}$ =
$10^{46}\,erg\,s^{-1}$, confirming a stratified structure of the torus, where
wavelength-dependent emissions originate from distinct regions of the torus.
By investigating the deviation from the best-fit torus size – luminosity
relation, we find a moderate correlation between the offset from the
$R_{dust}$–$L_{bol}$ relation and Eddington ratio. This suggests a possible
influence of the Eddington ratio on the observed flattening of the
$R_{dust}$–$L_{bol}$ relationship.

AGN — quasars — galaxies — reverberation mapping — torus — BLR

##  1 Introduction

Active galactic nuclei (AGNs) are amongst the most luminous objects in the
Universe, powered by the accretion of matter onto a supermassive black hole
(SMBH) and emit radiation across the entire electromagnetic spectrum (Rees,
1984). The unification model of AGN posits the presence of a dust torus
surrounding the broad emission line region (BLR). Based on different
obscuration by the torus, AGNs are classified into two main categories, namely
Type 1 and Type 2, characterized by the presence or absence of broad emission
lines in their optical/infrared spectrum, respectively (Seyfert, 1943; Urry &
Padovani, 1995).

In AGNs with UV luminosities ranging from $10^{42}$ to $10^{44}$
$erg\,s^{-1}$, the BLR typically resides at a distance of $\sim$ 0.003 to 0.02
pc from the accretion disk (AD), while the surrounding dust torus can extend
from 0.01 to 0.1 pc. Consequently, the central regions of AGNs are compact and
pose challenges for direct imaging techniques. Nevertheless, the advancements
in near-infrared (NIR) and MIR interferometry have enabled the measurement of
torus sizes in more than 40 sources, primarily focusing on nearby and bright
AGNs (Kishimoto et al., 2009, 2011a, 2011b; Burtscher et al., 2013; Gravity
Collaboration et al., 2020).

An alternative way to resolve the central part of AGNs is the technique of
reverberation mapping (RM; Blandford & McKee, 1982; Peterson, 1993). The
broad-band spectral energy distribution (SED) of an AGN presents an excess
emission in the IR, which is attributed to the thermal emission from the dusty
torus (Salpeter, 1964; Sanders et al., 1989; Elvis et al., 1994). By measuring
the delayed response ($\tau$) of the reprocessed IR continuum relative to the
ionizing UV/optical continuum from the AD, the torus size can be obtained as
$R_{torus}\,<\,c\times\tau$, where c is the speed of light.

To measure the inner extent of the torus, the dust RM (DRM) between the
optical and the NIR continuum in the K- band has been applied to $\sim$ 40
AGNs over a relatively narrow range of luminosity, i.e., V- band luminosity
$L_{V}\,<\,10^{45}\,erg\,s^{-1}$ (Suganuma et al., 2006; Koshida et al., 2014;
Mandal et al., 2018; Minezaki et al., 2019; Mandal et al., 2021). Koshida et
al. (2014) found that the dust torus radius ($R_{dust}$) correlated with the
optical luminosity ($L$) as $R_{dust}\,\propto\,L^{\alpha}$ with the exponent
$\alpha$ = 0.5, which is consistent with the expectation when the dust
particles are in radiative equilibrium with the dust sublimation temperature
($\mathrm{T_{sub}=1500K}$; Barvainis, 1987; Kishimoto et al., 2007). However,
Minezaki et al. (2019) reported a somewhat shallower slope with
$\alpha=0.424$, proposing that the anisotropic illumination of the AD leads to
a smaller dust-sublimation radius in the equatorial plane compared to the
polar directions. They also suggested that the dust torus does not respond
instantaneously to the changes in the luminosity of AD, leading to larger
torus sizes being more frequently observed in less luminous AGNs, which show
relatively small timescales of flux variation as the AD of low-luminosity AGNs
has dimmed (Koshida et al., 2009). Consequently, the average inner radius of
the dust torus tends to be larger for less luminous AGNs. In contrast, the
torus radius from the K-band interferometric measurements is roughy a factor
of two larger than the dust reverberation radius measured in the same band.
This difference could be attributed to the distinction between the flux-
weighted and response-weighted radii of the innermost dust torus obtained from
interferometry and DRM, respectively (Kishimoto et al., 2011b; Koshida et al.,
2014). These results indicate that the torus size and luminosity relation is
complex, requiring more detailed investigations.

Table 1: Details of the surveys used in this work

Name | Filter | Time Span | Cadence | $\mathrm{N_{epoch}}$ | Coverage | mag
---|---|---|---|---|---|---
|  |  | days ($\\#$/season) |  | $\mathrm{deg^{2}}$ |
(1) | (2) | (3) | (4) | (5) | (6) | (7)
|  |  | Sample 1 |  |  |
CRTS | unfiltered | 2005–2013 |  $\sim$ 20 | 15–120 | 33,000 $\mathrm{deg^{2}}$ | 19–21
PTF |  $g^{\prime}$, R | 2009–2014 |  $\sim$ 5–15 | 5–100 | 11,233 $\mathrm{deg^{2}}$ | 21.3($g^{\prime}$), 20.6 (R)
ASAS–SN | V | 2012–2019 |  $\sim$ 3–4 | 150–350 | All sky | 17
ZTF | ZTF-$g$, ZTF-$r$ | 2018–2021 |  $\sim$ 3–4 | 20–350 | 25,000–30,000 $\mathrm{deg^{2}}$ | 20.8($g$), 20.6($r$)
WISE/NEOWISE | W1, W2 | 2010–2020 |  $\sim$ 180 | 15–17 | All sky | 16.6 (W1), 16.0 (W2)
|  |  | Sample 2 |  |  |
SDSS(S82) | $g$ | 1998–2007 |  $\sim$ 5 |  $\sim$ 60 | 300 $\mathrm{deg^{2}}$ | 22.2
CRTS | unfiltered | 2005–2013 |  $\sim$ 21 |  $\sim$ 40 | 33,000 $\mathrm{deg^{2}}$ | 19–21
PTF | $g^{\prime}$ | 2009–2014 |  $\sim$ 5 |  $\sim$ 4 | 11,233 $\mathrm{deg^{2}}$ | 21.3
PS1 | $g(r,i)$ | 2011–2014 | 2/season |  $\sim$ 10 | 3$\pi$ | 22.0
ASAS–SN | V | 2012–2019 |  $\sim$ 3–4 |  $\sim$ 200 | All sky | 17
DES | $g(r,i)$ | 2013–2018 | 1–4/season |  $\sim$ 10 | 5100 $\mathrm{deg^{2}}$ | 23.6
ZTF | ZTF-$g$ | 2018 |  $\sim$ 3–4 |  $\sim$ 30 | 25,000–30,000 $\mathrm{deg^{2}}$ | 20.8
WISE/NEOWISE | W1, W2 | 2010–2020 |  $\sim$ 180 | 12–16 | All sky | 16.6 (W1), 16.0 (W2)

Note: Columns are (1) name of the survey, (2) name of the filter, (3)
observing time period of the light curves, (4) cadence of the light curves,
(5) number of data points present in the light curves, (6) sky coverage of the
corresponding survey, and (7) limiting magnitude in mag unit.

Table 2: Details of the initially selected samples with optical and IR (W1,
W2- band) light curve data

Sample | $\mathrm{N_{targets}}$ | redshift | $L_{bol}$ | WISE-band
---|---|---|---|---
|  | $z$ | $erg\,s^{-1}$ |
(1) | (2) | (3) | (4) | (5)
Sample 1 | 142 | 0.02 – 0.34 |  $10^{43.4}$ – $10^{46.5}$ | W1, W2
Sample 2 | 587 | 0.3 – 2.0 |  $10^{44.5}$ – $10^{47.6}$ | W1, W2

Note: Columns are (1) name of the sample, (2) initial number of AGNs present
in the sample, (3) redshift range of the sample, (4) AGN bolometric luminosity
range of the sample, and (5) name of the WISE bands used.

By expanding the K-band studies, Lyu et al. (2019) and Yang et al. (2020)
performed RM studies by utilizing various optical surveys and IR imaging data
from the all-sky Wide-field Infrared Survey Explorer (WISE) (Wright et al.,
2010; Mainzer et al., 2014). They focused on high-luminosity AGNs with
bolometric luminosities ($L_{bol}$) exceeding $10^{45}\,erg\,s^{-1}$,
effectively excluding the low-luminosity counterparts. Using the WISE W1 and
W2-band lag analysis, Lyu et al. (2019) derived reliable torus sizes for 67
among 87 Palomar–Green (PG) quasars at redshifts $z\leq 0.5$, reporting that
the torus size correlates with AGN bolometric luminosity with a slope of 0.47
(0.45) for W1- (W2) band lags. Using a large sample of 587 AGNs within the
SDSS stripe 82 region, Yang et al. (2020) performed a DRM study, although they
did not report an independent torus size – luminosity relationship. More
recently, Chen et al. (2023) conducted DRM analysis for 78 AGNs, for which
$H\beta$ lag measurements are available, reporting a much shallower slope,
$\alpha=0.37$ between the torus size and bolometric luminosity. They also
uncovered a correlation between the deviations from the best-fit torus size –
luminosity relation and the accretion rate (Du et al., 2014), concluding AGNs
with higher mass accretion tend to have a more shortened torus size because of
the self-shadowing effect of the slim disk (Kawaguchi & Mori, 2011; Wang et
al., 2014a).

While the thermal radiation from the dust torus dominates the IR flux of Type
1 AGNs, it includes emission from the AD. Consequently, the presence of AD
contamination in the observed IR light curves leads to derived dust lags
appearing shorter than the actual lag (Koshida et al., 2014; Mandal et al.,
2021), because the continuum emission from the AD in the IR is nearly
synchronous with the optical emission compared to the lag of the torus IR
emission (Tomita et al., 2006; Kishimoto et al., 2008; Lira et al., 2011).
Therefore, it is necessary to correct for AD contamination in the observed IR
emission in order to reliably measure the dust lag. Notably, while Lyu et al.
(2019) performed AD contamination correction in the observed IR fluxes, Chen
et al. (2023) and Yang et al. (2020) did not consider this correction,
potentially leading to an underestimation of their obtained torus sizes.

In a different approach, Li & Shen (2023) employed an ensemble of structure
functions (SFs) based on optical and IR light curves in W1- band of 587 AGNs
in the SDSS Stripe 82 region. They determined the torus size – luminosity
relation with an exponent of $\alpha=0.51$ by binning the luminosity range
into nine sub-samples and calculating the average effective torus size from
the ensemble SF analysis for each sub-sample. Unlike the RM technique, their
method relies on fitting parameters, such as the inclination angle
($\theta_{inc}$), inner radius ($R_{in}$), outer-to-inner radius ratio ($Y$),
torus half-opening angle ($\theta_{opn}$), and radial power-law index ($p$).

The previous studies of DRM showed a notable diversity in analysis method and
the slope of the correlation between torus size and luminosity varies among
these studies (See Lyu et al., 2019; Yang et al., 2020; Li & Shen, 2023; Chen
et al., 2023). Thus, it is necessary to further investigate the torus size –
luminosity relation using a uniform and consistent analysis method over a
large dynamic range of AGN luminosity. Such investigations would also
contribute to establishing AGNs as standard candles for constraining
cosmological parameters (Yoshii et al., 2014; Hönig et al., 2017).

In this paper, we revisit the torus size – luminosity ($R_{dust}$–$L_{AGN}$)
relation with a uniform analysis of the torus lag measurements, using of a
large sample of AGNs. First, our sample covers a large dynamic range in
luminosities ($10^{43.4}$ to $10^{47.6}$ $erg\,s^{-1}$) with a redshift
ranging from 0.02 to 2. Second, we employ a uniform method for cross
correlation analysis and perform quality assessment, ensuring consistent and
reliable lag measurements. Third, we investigate the effect of the AD
contamination on the observed IR fluxes in W1 and W2- bands, and correct for
the contamination in order to determine the $R_{dust}$–$L_{AGN}$ relationship.
The presentation of the paper is as follows. In Section 2, we describe the
sample and data acquired from different surveys and literature. We discuss the
construction of the optical and IR light curves in Section 3. In Section 4, we
describe the time series analysis. The results are presented in Section 5. We
discuss the wavelength dependent lag and scatter in our $R_{dust}$–$L_{AGN}$
relation in Section 6 and a summary in Section 7.

.

Figure 1: Distribution of redshift (left) and bolometric luminosity (right)
for different samples of AGNs with IR lag measurements.

##  2 Sample and Data

In order to investigate the relationship between the torus size and AGN
luminosity over a wide range of luminosities, we conducted a selection process
of Type 1 AGNs with available optical and IR light curve data. We selected
targets from two distinct categories based on their redshifts and
luminosities.

###  2.1 Sample 1 $(z<0.35)$

To explore the low luminosity end of the $R_{dust}$–$L_{AGN}$ relationship, we
selected targets from the catalog of Liu et al. (2019). This catalog comprises
a uniformly defined sample of 14,584 Type 1 AGNs compiled from SDSS DR7 at $z$
$<$ 0.35 and covering a range of AGN bolometric luminosities between
$10^{41.5}$ and $10^{46.6}$ $erg\,s^{-1}$. Since the AGNs in this catalog are
selected based on the robustness of the broad components of the Balmer lines
instead of imposing a lower limit on the FWHM of emission line, the catalog
contains many low BH mass ($M_{BH}$) AGNs with lower luminosity than other
existing Type 1 AGN catalogs. Thus, we can explore the low luminosity end of
the $R_{dust}$–$L_{AGN}$ relation using this catalog. The authors carefully
removed the stellar continuum by fitting their SDSS spectra. We used other key
parameters (i.e., optical luminosity at 5100 Å, $H\beta$, [O iii] and Fe ii
luminosities, FWHM of $H\beta$, Eddington ratio etc.) from the catalog.

###  2.2 Sample 2 $(z>0.3)$

To cover a moderate to high luminosity range in the $R_{dust}$–$L_{AGN}$
relationship, we used the SDSS Stripe 82 (S82) region (MacLeod et al., 2012),
where a total of approximately 9258 spectroscopically confirmed broad-line
AGNs are available. These AGNs span a wide range of redshifts and
luminosities, useful for investigating luminous AGN population in the
$R_{dust}$–$L_{AGN}$ relation. Yang et al. (2020) performed a DRM study of S82
AGNs using $\sim$ 20 yr ground-based optical light curves and 10 yr IR light
curves in the W1- band from the WISE. They constructed optical light curves
using data from the SDSS, Pan-STARRS (PS1; Chambers et al., 2016), Dark Energy
Survey (DES; Abbott et al., 2018), Catalina Real-Time Transient Survey (CRTS;
Drake et al., 2009), Palomar Transient Factory (PTF; Law et al., 2009), Zwicky
Transient Facility (ZTF; Bellm et al., 2019) and All-Sky Automated Survey for
Supernovae (ASAS-SN; Shappee et al., 2014). They were able to find lag between
the optical and W1- bands for 587 AGNs with redshift in the range 0.3 $<$ $z$
$<$ 2, and luminosity in the range $10^{44.5}$ to $10^{47.6}$ $erg\,s^{-1}$.
They employed JAVELIN (Zu et al., 2011) to search for the lags between the
optical and W1- band light curves. Thus, Sample 2 comprises a total of 587
AGNs from S82 region. We used the same object ID number as assigned by Yang et
al. (2020) for Sample 2.

###  2.3 Light curves data

We gathered optical and IR data from both ground-based and space-based
telescopes, as well as from the relevant literature for our target AGNs.

####  2.3.1 Optical data from ground based surveys

For Sample 1, we compiled the available optical photometric data from
different ground-based transient surveys. Initially, we conducted a search for
data from the ASAS-SN survey, as it overlaps with the major part of the IR
light curves from WISE and NEOWISE. Since 2012, ASAS-SN monitored nearly the
entire visible sky to a V- band depth of 17 mag with a cadence of 2–3 days
using twenty telescopes each of 14 cm diameter. Therefore, we selected AGNs
with $g<$ 17.5 mag from our sample of 14,584 Type 1 AGNs, which resulted in
732 targets. We searched light curve data from the publicly available ASAS-SN
light curve servers, i.e., ASAS-SN Photometry Database 111https://asas-
sn.osu.edu/photometry (Shappee et al., 2014; Jayasinghe et al., 2019). Among
the 732 targets, we found ASAS-SN light curve data for 142 AGNs. Thus in
Sample 1, we arrived at an initial number of 142 sources.

We also collected V- band photometric light curve data for the 142 targets
from the CRTS. This survey uses the 0.7 m Catalina Schmidt Telescope in
Tucson, Arizona, the 1.5 m telescope on Mt. Lemmon, Arizona and the 0.5 m
Siding Springs Survey (SSS) telescope in Australia to search for optical
transients with timescales of minutes to yr in the northern sky. CRTS uses
SExtractor (Bertin & Arnouts, 1996) to process the images following standard
techniques. The limiting V- band magnitude for the 1.5 m and 0.7 m are 21.5
mag and 19.5 mag, respectively. The photometry for sources brighter than 13
V-mag is problematic because of saturation and non-linearity response of the
CCD.

In addition, we obtained photometric light curve data from the PTF Data
Release 3 (DR3) when available. This survey used the 48 inch Samuel Oschin
Telescope at the Palomar Observatory in northern San Diego County, California.
The photometry was performed in two broadband filters, SDSS $g$′ and Mould
$R$. The survey provides a depth of $m_{g}$ $\sim$ 21.3 mag and $m_{R}$ $\sim$
20.6 mag with an exposure of 60 s. The PTF data were used to fill the gap
between the CRTS and ASAS-SN light curves.

Finally, we gathered light curve data from the ZTF Data Release 14 (DR14). ZTF
is a wide-field time-domain survey, that uses the 1.2 m Schmidt telescope at
the Palomar Observatory with 47 deg2 field of view. It provides photometric
light curve data with a cadence of 3–4 days in the $gri$ bands since 2018
March (MJD $\sim$ 58194). It scans the entire northern visible sky to median
depths of $g$ $\sim$ 20.8 mag and $r$ $\sim$ 20.6 mag (AB, 5$\sigma$) with a
30 s exposure (Masci et al., 2019). ZTF uses point-source-function-fit
photometry to build the light curves. We selected ZTF data with the condition
catflags=0 to avoid the data contaminated by clouds, moonlight or bad seeing.

In summary, we collected optical photometric light curve data for 142 Type 1
AGNs in Sample 1 with a baseline from 2005 to 2021. This sample consists of
targets with bolometric luminosity within the range of $\mathrm{10^{43.4}}$ to
$\mathrm{10^{46.5}}$ $\mathrm{erg\,s^{-1}}$ at $z<0.3$.

We utilized the publicly available light curve data in the optical band from
Yang et al. (2020) for 587 Type 1 AGNs in Sample 2 for our own analysis.

Figure 2: From top left to bottom right: optical and IR W1 and W2- bands light
curves of the targets (a) OB21, (b) OB44, (c) OB64, and (d) OB76 from Sample
1, as well as (e) ID24646, and (f) ID2497212 from Sample 2. The final
calibrated optical light curves including data from CRTS, ASAS-SN, PTF and ZTF
from Sample 1 are shown in the top panel, while for Sample 2 finally
calibrated DES $g$\- band optical light curves are presented. The bottom two
panels represent the IR light curves in W1 and W2- bands. The W1 and W2- bands
light curves shown in the figure are not corrected for the AD contamination.
The shifted (by the lag from ICCF) and scaled version of the optical light
curves (yellow circles) are also shown with the observed IR light curves,
where the shade represents the $1\sigma$ uncertainty obtained in the scaling
process.

####  2.3.2 IR data from WISE

WISE surveys the full sky at wavelengths centered at 3.4 (W1), 4.6 (W2), 12
(W3) and 22 $\mu m$ (W4) using a 40 cm diameter infrared telescope since 2010.
We collected profile-fit photometric light curve data of our targets in the W1
and W2– bands from the Near-Earth Object Wide-field Infrared Survey Explorer
Reactivation Database (NEOWISE-R), WISE All-Sky and WISE Post-Cryo Database
from the NASA/IPAC Infrared Science Archive
222https://irsa.ipac.caltech.edu/cgi-bin/Gator/nph-
scan?mission=irsa&submit=Select&projshort=WISE. The retrieved WISE data covers
a time period from 2010 to 2020. There is a gap of 3 yr between 2011 and 2014.
We only selected good-quality photometric measurements following the
suggestions given at
https://wise2.ipac.caltech.edu/docs/release/neowise/expsup/sec2_3.html.

For Sample 2, we acquired W1- band light curve data from Yang et al. (2020),
while the light curve data in the W2- band was retrieved from the WISE
database.

The survey details utilized in this study are summarized in Table 1, while an
overview of the initially selected sample properties can be found in Table 2.
Figure 1 presents a comparison between the redshift and luminosity
distributions of our initially selected sample and those reported in the
literature for DRM. Sample 1 covers a range of luminosity with redshifts below
0.33, while Sample 2 mainly represents high-luminosity and high-redshift AGNs.
Therefore, our total sample encompasses the range of luminosity and redshift
that has been used to measure the torus size in AGNs to date.

##  3 Light Curves Construction

###  3.1 Photometric Calibration of the optical light curves

The optical light curves originate from different surveys, each observed with
distinct telescopes and filter properties. For instance, the CRTS observations
are unfiltered. We converted the CRTS magnitudes to V- band fluxes following
Hovatta et al. (2014). While ASAS-SN provides light curve data in V- band, PTF
have data either or both in SDSS-$g$ and Mould R- bands. The ASAS-SN and PTF
magnitudes were converted to fluxes using the conversion factors given in
Bessell et al. (1998) and Fukugita et al. (1996) in their respective filters.
We also utilized ZTF-$g$ band data to construct our optical light curves and
used ZTF-$r$ when ZTF-$g$ is not available. The ZTF magnitudes were converted
into fluxes using the zero-points provided by the Spanish Virtual Observatory
(SVO) filter profile service 333http://svo2.cab.inta-
csic.es/theory/fps/index.php?mode=browse&gname=Palomar&gname2=ZTF&asttype=. We
took the average of the fluxes if more than one observation were taken in a
single night.

Various optical monitoring surveys utilize different filters and adopt various
aperture schemes in photometry, generating systematic offset among the
continuum light curves. Thus, it is necessary to inter-calibrate heterogeneous
data sets and include systematic errors. We inter-calibrated the light curves
from different surveys with various optical bands using the python-based
software PyCALI444https://github.com/LiyrAstroph/PyCALI(Li et al., 2014).
PyCALI uses Bayesian statistics to determine the best scaling factors for
inter-calibration by modeling each light curve using damped random walk model
(DRW; Kelly et al., 2009; MacLeod et al., 2010) for describing AGN
variability. In this process, we used ASAS-SN light curve as the reference to
scale all the other light curves, because ASAS-SN covers a significant portion
of the overlapping region between the optical and WISE light curves. In this
inter-calibration we ignored the small optical inter-band continuum lags
(accretion disk RM; Jiang et al., 2017; Mudd et al., 2018; Guo et al., 2022),
which are expected to be a few days ($<$ 10 days; Jha et al., 2022),
negligible compared to IR time scale. In our subsequent analysis, we adopted
the total error computed by PyCALI, which accounts for measurement and
systematic errors. For Sample 2, we used DES $g$\- band light curves that were
finally calibrated using various optical data sources from Yang et al. (2020).

###  3.2 WISE Light Curves

The WISE magnitudes are obtained using profile-fitting photometry and circular
aperture photometry
555https://wise2.ipac.caltech.edu/docs/release/allsky/expsup/sec4_4c.html.
WISE Vega magnitudes were then converted to AB magnitudes using
$m_{AB}\,=\,m_{Vega}\,+\,\Delta m$, where $\Delta m$ = 2.699, 3.339, 5.174,
and 6.620 for the W1, W2, W3, and W4- bands, respectively (Jarrett et al.,
2011), and finally into flux densities following Jarrett et al. (2011). WISE
observed each object for 15–17 epochs, with 10–20 frames within a 36 hour
window in each epoch during the period 2010–2020. Thus, we calculated the mean
value of the IR fluxes within a 6-month window for both W1 and W2- band light
curves.

The light curves were not corrected for the host galaxy contribution
considering a constant host star-light contamination, which makes a negligible
effect on the lag measurements. We obtained the errors in the flux
measurements through the standard propagation of errors. In Figure 2, we show
light curves for 6 targets in the optical and IR bands. The displayed light
curves reveal a correlated flux variations between the optical and IR- (W1,
W2) bands for the selected AGNs. A few more light curves for objects with
various lags and luminosities are shown in Figures 16 and 17 for Sample 1 and
Sample 2, respectively. The light curve data are presented in Table 3.

##  4 Analysis

The final calibrated optical and WISE light curves, with a typical cadence of
5 days in the optical and 180 days in the WISE light curves, span over 16–20,
and 10 yr, respectively. In this section we describe our lag measurements
based on these light curves, quality assessment, and correction for $z$ and AD
contamination in WISE light curves.

= -15mm  ID no MJD $F$ $F_{err}$ Survey mJy mJy (1) (2) (3) (4) (5) OB1
53627.295443 1.2872 0.1976 CRTS OB1 53637.455742 1.2723 0.1978 CRTS OB1
53646.362378 1.2807 0.1978 CRTS … … … … … … … … … …

Table 3: Optical and IR photometric light curves. The full table is available
for Sample 1 and Sample 2 separately in a machine-readable form online.

Note: Columns are (1) object ID number, (2) Modified Julian Date, (3) flux,
(4) error in flux, and (5) name of the imaging survey. Data table for the
observed and AD contamination corrected IR light curves for Sample 1 and
Sample 2 will be provided separately with the same format.

Figure 3: From top left to bottom right: the ICCFs and JAVELIN lag analysis
between optical and W1- band light curves (not corrected for AD contamination)
for the targets (a) OB21, (b) OB44, (c) OB64, and (d) OB76, respectively, from
Sample 1. Black solid line represents the ICCF between the optical and W1-
band light curves and dashed green line is the ACF of the optical light curve.
The CCCD obtained from ICCF is shown by the blue histogram, whereas the red
histogram represents the lag probability distribution of the time delay from
JAVELIN. The last two in the bottom panel are the ICCFs for the targets (e)
ID24646, and (f) ID2497212 from Sample 2, where the vertical red dashed lines
represent the JAVELIN lag measurements from Yang et al. (2020) with 1$\sigma$
uncertainties shown by the red dotted lines. The observed-frame time lags
obtained from ICCF with their uncertainties are also mentioned.

###  4.1 Lag measurements between the optical and IR light curves

All the AGNs selected for this study exhibit synchronized flux variations over
time in both the optical and IR W1 and W2- bands, with the optical
observations consistently leading the IR measurements. To determine the time
lag between the optical and IR continua, we employed the Interpolated Cross-
Correlation Function (ICCF; Gaskell & Sparke, 1986; Gaskell & Peterson, 1987)
as the primary method. It determines time lags by linearly interpolating two
light curves and calculating the cross-correlation coefficient (r). We
evaluated CCF by setting the lag search range to $\sim$ 70 % of the baseline
defined by the optical light curve, which is sufficient for robust time lag
measurements (Grier et al., 2017, 2019; Yang et al., 2020). We applied the
ICCF method to measure time lags in 142 AGNs from Sample 1 and 587 AGNs from
Sample 2, respectively. We first located the most significant peak within the
search range and the lag was estimated from the centroid ($\tau_{cent}$) of
the CCF around the most significant peak, defined as (Peterson et al., 1998)

$\tau_{cent}=\frac{\Sigma_{j}\,\tau_{j}CCF_{j}}{\Sigma_{j}\,CCF_{j}}$ (1)

where, $\tau_{cent}$ was obtained by considering all the points within 80 % of
the maximum of the CCF ($r_{max}$). The uncertainties in the measured lags
were estimated using a model-independent Monte Carlo simulation based on flux
randomization (FR) and random subset selection (RSS) as described in Peterson
et al. (1998); Wandel et al. (1999); and Peterson et al. (2004). This process
was repeated for 5000 iterations, and the centroid of the lag was determined
in each iteration. The cross-correlation centroid distribution (CCCD) was then
constructed. As the CCCD was not Gaussian, the lag was evaluated from the
median of the distribution and the lower and upper uncertainties in the
measured lag are that values at the 15.9 and 84.1 percentile of the
distribution. This is equivalent to $1\sigma$ error for a Gaussian
distribution.

We also used JAVELIN to find the time delays between the optical and IR light
curves and validated the ICCF measurements. In JAVELIN the driving optical
continuum was first modeled using the DRW model with two free parameters;
amplitude ($\sigma_{d}$) and time scale of variability ($\tau_{d}$). To
recover the IR continuum light curve a top-hat transfer function was convolved
with the driving optical continuum light curve. We used the same search range
as used for ICCF to find the lag, and adopted Markov Chain Monte Carlo
approach to find the best-fitting model maximizing the likelihood. However, we
observed that JAVELIN tends to produce strong aliases, particularly near the
boundaries of the lag search window when a large lag search range was applied,
leading to multiple peaks in the lag posterior distribution. This is a feature
of aliasing in the light curves resulting from sparsely sampled multi-year
data with seasonal gaps (Grier et al., 2019; Homayouni et al., 2020). Thus, we
excluded the false peaks in the posterior distribution that appeared near the
boundary of the lag-search window. The final lag was obtained from the median
of the lag probability distribution around the most significant peak of the
distribution guided by the ICCF peak.

Figure 3 shows examples of the results from ICCF and JAVELIN for the targets
shown in Figure 2. For a few targets, along with a primary peak at smaller
lag, there is a secondary peak at a larger lag (for examples, see the CCFs of
OB44 and OB64 in Figure 3). While the primary peak represents the actual lag
between the optical and W1- band light curves, the other peak arises mainly
because of the shape and aliasing in the light curves (Guise et al., 2022).
For example, in Figure 3, the CCF of the target OB44 shows a primary peak at
around 110 days, and a secondary peak around 2200 days. The corresponding Auto
Correlation Function (ACF) of the optical light curve has also two peaks, one
around 0 lag and another secondary peak around 2090 days. Thus, the IR W1
light curve will also correlate with the optical light curve at the lag
$\tau$=109 days and at a latter time $\tau$ \+ 2090 days, resulting two peaks
in the CCF. The same is applicable for the target OB64 shown in Figure 3.
Thus, we calculated lags around the primary peaks in these cases. We obtained
te lags between the optical and WISE W1 and W2- band light curves for all of
those 729 (142+587) AGNs.

Figure 4:  The probability $p(r_{max})$ represents the likelihood of
uncorrelated light curves producing a cross-correlation coefficient equal to
or exceeding $r_{max}$. This probability is presented as a function of the
maximum cross-correlation coefficient $r_{max}$ in the left and right panels
for W1 and W2- bands, respectively. The vertical dashed line represents
$r_{max}$=0.6, whereas $p(r_{max})$=0.2 is shown by the horizontal dashed
line. The reliable time lags are defined with $r_{max}$ $\geq$ 0.6, and
$p(r_{max})$ $<$ 0.2.

Figure 5: Examples of rejected targets based on ICCFs between the optical and
W1- band light curves. The two targets ID7907462 (left) and ID954476 (right)
from Sample 2 are the examples with negative lag and a flat CCF, respectively.
The JAVELIN lag measurements from Yang et al. (2020) for these two targets are
shown by the dashed red lines. The observed-frame time lags obtained from ICCF
with their uncertainties are also mentioned at each panel. Figure 6:
Comparison of the time lags in the W1- band obtained through ICCF and the
JAVELIN method in the observed-frame of the targets. The dashed blue line is
the 1:1 relation. The brown triangles and circles represent the accepted
targets from Sample 1 and Sample 2, respectively, for use in the RM analysis.
The black points are the rejected targets. For Sample 2, the JAVELIN lags are
obtained from Yang et al. (2020). Most of the rejected targets deviate from
the 1:1 relation.

###  4.2 Quality assessment

In this section, we assess the dependability of the lag measurements obtained
through cross-correlation analysis. The accuracy of the recovered time lags
between the optical light curves, derived from various surveys, and sparsely
sampled IR light curves is influenced by different sources of uncertainty. The
presence of substantial data gaps in the IR light curves, coupled with a
limited number of data points in the optical-IR light curve pairs, can lead to
false positive peaks near the boundaries of the lag search window when using
the ICCF method. These factors collectively impact the reliability of the
obtained time lags. Consequently, it is crucial to exercise caution when
selecting targets with significant and reliable lags.

The maximum cross-correlation coefficient ($r_{max}$) is a key factor to
assess the reliability in lag measurements. This coefficient signifies the
strength of the correlation between two light curves. Previous DRM studies
have often employed specific thresholds for $r_{max}$, such as $r_{max}$ $>$
0.5 (Yang et al., 2020) or $r_{max}$ $>$ 0.4 (Chen et al., 2023), as a means
to determine the reliability of lag measurements.

However, it is crucial to acknowledge that relying solely on $r_{max}$ is
insufficient for guaranteeing lag reliability. There are situations where
$r_{max}$ might be relatively high even when the light curves possess only a
limited number of data points but exhibit strong linear variation. Even small
variability amplitude compared to the flux uncertainties in the light curves
may result high $r_{max}$ as long as certain features in the light curves
align well. Moreover, it is crucial to emphasize that the calculations for
$r_{max}$ do not incorporate the uncertainties in flux.

To address these issues and evaluate the reliability of lag measurements, a
new method called $\tt PyI^{2}CCF$ 666https://github.com/legolason/PyIICCF has
been introduced (Guo & Barth, 2021; U et al., 2022, Guo et al., in prep). This
method employs a null-hypothesis approach to assess whether the observed
cross-correlation could be equally derived from two uncorrelated red-noise
light curves. To put this method into practice, we utilized the publicly
available $\tt PyI^{2}CCF$ tool to investigate the reliability of our lag
measurements. Our methodology involved generating $10^{3}$ mock light curve
pairs for each object, maintaining the same sampling and cadence as the
observed light curves while utilizing DRW models. Subsequently, we calculated
the lag significance indicator, referred to as $p(r_{max})_{\tau>0}$, by
tallying the occurrences of positive lags ($\tau\,>$ 0) where the $r_{max}$
value exceeded the observed $r_{max}$ across all simulations. It is important
to mention that we retained the same searching window that was initially used
during the lag calculations. Finally, we apply a set of criteria to exclude
measurements deemed unreliable, such as:

  1. 1.

We define the lag measurements as reliable if $r_{max}$ $\geq$ 0.6 and
$p(r_{max})_{\tau>0}$ $<$ 0.2. These thresholds are adopted following U et al.
(2022) and Woo et al. (2024) (See Figure 4). Using these criteria, we identify
a total of 85 targets in the W1- band and 79 in the W2- band within Sample 1
as having reliable lags. Similarly, for Sample 2, out of 587 targets, we
select 399 in the W1- band and 372 in the W2- band.

  2. 2.

However, we discard 4 and 3 targets from Sample 1 in the W1 and W2- bands,
respectively, due to negative lags within a 1$\sigma$ uncertainty range, as
negative lags are physically unacceptable (See Figure 5, left). Similarly,
from Sample 2, 16 and 14 targets are rejected with negative lags within
1$\sigma$ uncertainty in W1 and W2- bands, respectively. Additionally, we
eliminate targets displaying a flat CCF spanning over thousands of days, as
this condition results in an average lag computed over the entire search range
(refer to Figure 5, right). Consequently, 18 targets from Sample 2 in both W1
and W2- bands are excluded from the analysis.

After applying the criteria, we narrowed down from 142 targets to 81 targets
in the W1- band and 76 targets in the W2- band in Sample 1. In the case of
Sample 2, we obtained 365 AGNs with W1- band light curves and 340 with W2-
band light curves, out of the initial sample of 587 targets. In summary, our
final sample consists of 446 (416) AGNs with a reliable and significant dust
lag measurements with W1- (W2) band light curves.

In Figure 6, we present a comparison of the W1- band lags in the observed-
frame between the ICCF and JAVELIN methods. The final sample is represented by
brown data points, while the rejected AGNs are shown in black. The majority of
the rejected targets deviate from the expected 1:1 relation. In contrast, the
lags obtained for the final sample demonstrate a good consistency, with a
scatter of 0.17 dex. We present a few examples of the selected ICCF and
JAVELIN analyses for Sample 1 and Sample 2 in Figures 18 and 19, respectively.

We checked that shifting the optical light curve by the median lag value
($\tau_{cent}$) obtained from ICCF resulted in a good alignment with the IR
light curves, as demonstrated in Figures 2, 16, and 17. Moreover, the error
estimation in ICCF analysis is more conservative (Guo et al., 2022), and the
ICCF measurements incorporate a cross-correlation reliability test.
Consequently, in all subsequent discussions, our primary focus will be on the
lags derived from the ICCF analysis.

= -75mm

Table 4: Details of the final sample and the results of their time lag
measurements for W1 and W2- bands are provided. The complete table is
available online in a machine-readable format, separately for Sample 1 and
Sample 2.

ID no | ID | $z$ | log($L_{bol}$) | log(Edd) | $r_{max,W1}$ | pval_W1 | $\tau_{W1}$ | $\tau_{W1,ADC}$ |  $\tau_{W1,ADC}$ (JAV) | $r_{max,W2}$ | pval_W2 | $\tau_{W2}$ | $\tau_{W2,ADC}$ |  $\tau_{W2,ADC}$ (JAV)
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
|  |  | $\mathrm{erg\,s^{-1}}$ |  |  |  | days | days | days |  |  | days | days | days
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10) | (11) | (12) | (13) | (14) | (15)
|  |  |  |  |  | Sample 1 |  |  |  |  |  |  |  |
OB1 | $J014017.07-00503.0$ | 0.334 | $46.14\pm 0.01$ | $-0.68\pm 0.14$ | 0.85 | 0.0811 | $617^{+105}_{-99}$ | $773^{+189}_{-178}$ | $1059^{+97}_{-188}$ | 0.82 | 0.1421 | $619^{+117}_{-113}$ | $795^{+209}_{-208}$ | $1076^{+104}_{-211}$
OB2 | $J031027.83-004950.8$ | 0.080 | $45.05\pm 0.01$ | $-0.83\pm 0.14$ | 0.74 | 0.0089 | $73^{+51}_{-33}$ | $133^{+91}_{-78}$ | $163^{+103}_{-91}$ | 0.73 | 0.0093 | $98^{+74}_{-39}$ | $165^{+87}_{-95}$ | $179^{+102}_{-63}$
OB3 | $J032213.90+005513.5$ | 0.185 | $45.69\pm 0.01$ | $-0.23\pm 0.14$ | 0.71 | 0.1689 | $375^{+110}_{-150}$ | $400^{+116}_{-145}$ | $616^{+29}_{-27}$ | 0.70 | 0.1953 | $422^{+130}_{-164}$ | $434^{+130}_{-172}$ | $619^{+23}_{-13}$
… | … | … | … | … | … | … | .. | … | … | … | … | … | … | …
… | … | … | … | … | … | … | .. | … | … | … | … | … | … | …
|  |  |  |  |  | Sample 2 |  |  |  |  |  |  |  |
5871 | $J000535.23+000021.6$ | 0.736 | $45.28\pm 0.01$ | $-0.82\pm 0.57$ | 0.78 | 0.1716 | $594^{+204}_{-267}$ | $707^{+291}_{-353}$ | $985^{+179}_{-188}$ | … | … | … | … | …
2401 | $J000654.10-001533.5$ | 1.725 | $46.99\pm 0.01$ | … | … | … | … | … | … | 0.83 | 0.0760 | $2571^{+287}_{-233}$ | $2500^{+464}_{-511}$ | $1510^{+1072}_{-653}$
9836 | $J001735.69-011325.0$ | 0.805 | $45.89\pm 0.01$ | $-0.27\pm 0.18$ | 0.87 | 0.0829 | $437^{+162}_{-178}$ | $585^{+208}_{-259}$ | $795^{+200}_{-243}$ | 0.75 | 0.1909 | $782^{+290}_{-453}$ | $914^{+354}_{-424}$ | $971^{+252}_{-231}$
… | … | … | … | … | … | … | .. | … | … | … | … | … | … | …
… | … | … | … | … | … | … | .. | … | … | … | … | … | … | …

Note: Columns are (1) object ID number, (2) SDSS object ID, (3) redshift, (4)
AGN bolometric luminosity, (5) Eddington ratio obtained from this work, (6)
peak correlation coefficient from ICCF between optical and W1- band light
curve, (7) $p(r_{max})_{\tau>0}$ value obtained between optical and W1- band
light curve, (8) observed-frame W1- band time lag from ICCF, (9) observed-
frame AD contamination corrected W1- band time lag from ICCF, (10) observed-
frame AD contamination corrected W1- band time lag from JAVELIN, (11) peak
correlation coefficient from ICCF between optical and W2- band light curve,
(12) $p(r_{max})_{\tau>0}$ value obtained between optical and W2- band light
curve, (13) observed-frame W2- band time lag from ICCF, (14) observed-frame AD
contamination corrected W2- band time lag from ICCF, and (15) observed-frame
AD contamination corrected W2- band time lag from JAVELIN.

Figure 7:  The distribution of the lag-ratio ($\mathrm{\tau_{ratio}}$ =
$\tau_{W2}/\tau_{W1}$) between the dust lags in W2 and W1- bands with and
without AD contamination correction shown by red and black lines,
respectively. The median values with their standard deviations are mentioned
in the figure.

###  4.3 Correcting for redshift effects on dust lags

The AGNs under investigation in this study span a wide range of redshifts as
illustrated in Figure 1. Because of the wavelength-dependent nature of dust
lag, objects at higher redshifts will exhibit shorter rest-frame time lags in
specific wavelength bands, than the same luminosity AGNs at lower $z$ (Sitko
et al., 1993; Oknyanskij & Horne, 2001; Yoshii et al., 2014; Minezaki et al.,
2019; Sobrino Figaredo et al., 2020; Chen et al., 2023). This redshift effect
must be corrected, especially when dealing with a sample over a large redshift
range.

We used the correction factor in the form of (1+$z$)γ, as introduced by Yoshii
et al. (2014) and Minezaki et al. (2019), to account for wavelength-dependent
dust lags in order to remove the redshift effect. First, we determined the lag
ratio, $\mathrm{\tau_{ratio}}$ = $\tau_{W2}/\tau_{W1}$, between dust lags in
W2 and W1- bands for our selected final sample and computed $\gamma$ using the
formula, $\gamma$ =
log($\tau_{W2}/\tau_{W1}$)/log($\lambda_{W2}/\lambda_{W1}$). Figure 7
represents the distribution of $\mathrm{\tau_{ratio}}$, both with and without
AD contamination correction (see Section 4.4 for details of AD contamination
correction) obtained from ICCF analysis. Our analysis reveals median values of
$\mathrm{\tau_{ratio}}$ to be $1.19\pm 0.43$ ($1.16\pm 0.55$) and $1.22\pm
0.50$ ($1.26\pm 0.85$) for the AD contamination corrected and uncorrected
cases, respectively, obtained from ICCF (JAVELIN) analysis. These values
closely align with the median values of $\mathrm{\tau_{ratio}}$ = 1.15
reported by Lyu et al. (2019) and 1.26 reported by Chen et al. (2023) within
the $1\sigma$ range of the distribution.

By adopting a median $\mathrm{\tau_{ratio}}$ of 1.20 for both the AD
contamination corrected and uncorrected cases, we find that $\gamma\,\sim
0.62$. Previously Minezaki et al. (2019) empirically determined $\gamma$ as
1.18 based on the $K$\- and $H$\- band lag measurements in six nearby Seyfert
galaxies. In the case of W1 and W2-band lag measurements, Chen et al. (2023)
reported a $\gamma$ $\sim$ 0.76 using a sample of 78 AGNs. Thus, we conclude
that our study on the redshift dependency is broadly consistent with the
previous studies. Combined with the cosmological time dilation correction
factor of (1+$z$)-1, the total correction factor for the observed-frame time
lag becomes approximately (1+$z$)γ-1 $\sim$ (1+$z$)-0.38. We uniformly applied
this total correction factor to the measured time lag in the observed-frame,
thereby obtaining the rest-frame time lags for the total sample.

Figure 8: Comparisons of torus sizes with and without AD contamination
correction are presented in W1 (left) and W2- (right) bands, color-coded by
Eddington ratio. The AGNs which do not have Eddington ratio measurements are
shown by the open grey circles. The best fits to all the data points are shown
by the dashed brown lines, while the dashed black lines represent the 1:1
relation. The intrinsic scatter from the best fit is also mentioned at each
panel. In the bottom panels, we represent the logarithm of the size ratio as a
function of Eddington ratio.

###  4.4 Correction for the accretion disk (AD) contamination in the IR W1 and
W2- band light curves

The observed IR fluxes in W1 and W2- bands contain contribution from the torus
as well as from the AD (Tomita et al., 2006; Kishimoto et al., 2008; Lira et
al., 2011), the latter of which makes the delay time between the optical and
IR light curves shorter than the actual lag. Thus, the AD contamination needs
to be removed from the observed IR fluxes to get the actual time lag. The AD
emission in the IR can be estimated by considering a power-law spectrum of the
accretion disk (Koshida et al., 2014) as follows

$F_{IR}^{AD}(t)=F_{OP}(t)\left(\frac{\nu_{IR}}{\nu_{OP}}\right)^{a}$ (2)

where, $F_{IR}^{AD}(t)$ and $F_{OP}(t)$ are the AD component of the IR flux
and the optical flux in V- band ($g$\- band for Sample 2) at time t,
respectively. $\nu_{OP}$ and $\nu_{IR}$ are the effective frequencies of
optical in V- band ($g$\- band for Sample 2) and IR (W1, W2) bands,
respectively, and $a$ is the power-law index. We assumed a constant power-law
index of $a=1/3$ as suggested by the standard disk model (Shakura & Sunyaev,
1973) to calculate $F_{IR}^{AD}(t)$. Note that the polarized light spectra in
the NIR strongly suggest intrinsic AD spectra with $a=1/3$ (Kishimoto et al.,
2008). In contrast, Koshida et al. (2014) utilized not only $a=1/3$ but also
$a=0.0$, taking into account the potential maximal contribution of AD
contamination in the IR bands, while Minezaki et al. (2019) used $a=0.1$
following Yoshii et al. (2014). If $a$ is smaller than $a=1/3$, the time lag
generally becomes larger, because the AD contamination in the IR-band flux is
larger (i.e., a larger correction). To empirically test the effect, we
performed the ADC with $a=0.0$ using a representative subsample of 35 AGNs,
which covers the entire range of luminosity and Eddington ratio of the total
sample. We obtained insignificant effect, with a 0.03 dex increase in both W1
and W2- band lag measurements (see also Koshida et al., 2014; Minezaki et al.,
2019).

We estimated $F_{IR}^{AD}(t)$ for each epoch of IR (W1 and W2) observations
using the average flux values of the quasi-simultaneous epochs (within $\pm$ 5
days) in the optical light curve. If for any epoch in IR light curve, there is
no quasi-simultaneous observations in the optical, we obtained the
corresponding optical flux from the DRW modeled optical light curve using
JAVELIN. It is worth noting that the influence of host galaxy flux
contamination in the optical band is more evident in lower redshift AGNs. This
may introduce a spurious redshift dependence in the dust lags. However, this
host-galaxy flux is assumed to remain constant throughout the entire
observational period. Finally, $F_{IR}^{AD}(t)$ was subtracted from the
observed IR (W1, W2) flux ($F_{IR,obs}$) to get the IR flux, which is coming
from the torus ($F_{IR,dust}$) as

$F_{IR,dust}(t)=F_{IR,obs}(t)-F_{IR}^{AD}(t)$ (3)

Following this analysis, our next objective is to investigate how the AD
contamination corrected and uncorrected dust lags are influenced by AGN
luminosity and the Eddington ratio for the selected AGNs in our study. The AGN
luminosity can be contaminated by stellar-light from the host galaxies. Liu et
al. (2019) modeled and subtracted the stellar continua from the spectra and
derived the continuum luminosity at rest-frame 5100 Å ($L_{5100{\AA}}$). We
used their estimated $L_{5100{\AA}}$ for Sample 1 ($z<$ 0.33) and obtained the
AGN bolometric luminosity ($L_{bol}$) by using a bolometric correction factor
of 9.26 (Richards et al., 2006). The AGN bolometric luminosities for Sample 2
were collected from Yang et al. (2020). The same bolometric correction factor
was used for Sample 2. In order to obtain the Eddington ratios for our own
sample, we retrieved FWHM of $H\beta$ values from Liu et al. (2019) and Shen
et al. (2011). Subsequently, the $H\beta$ BLR size ($R_{H\beta}$) were
estimated by utilizing the $R_{H\beta}$–$L_{5100{\AA}}$ relation from Woo et
al. (2024). Then, we employed the virial equation (Woo & Urry, 2002),
incorporating the FWHM of $H\beta$ and a virial factor of $f_{BLR}=1.12$ (Woo
et al., 2015), to derive the single epoch $M_{BH}$. Finally, we obtained the
Eddington ratio ($L_{bol}/L_{Edd}$) using Eddington luminosity
$L_{Edd}=1.26\times 10^{38}M_{BH}/M_{\odot}$. These Eddington ratios were
determined for approximately $\sim$72$\%$ of the AGNs in our final sample,
specifically those with available FWHM values for $H\beta$.

In Figure 8, we present a comparison between the redshift-corrected rest-frame
torus sizes, with and without AD contamination correction, denoted as
$R_{W1,ADC}$ ($R_{W2,ADC}$) and $R_{W1,non-ADC}$ ($R_{W2,non-ADC}$),
respectively, for the W1- (W2) band. A linear fit to the torus sizes obtained
with/without correction for AD contamination results

$\begin{split}&\parbox{1084.04999pt}{$\text{log}(R_{W1,\text{ADC}}/400\,\text{light-
day})=(1.12\pm 0.02)\text{log}(R_{W1,\text{non-ADC}}\\\ /400\,\text{light-
day})+(0.05\pm 0.01)$ }\end{split}$ (4)

and

$\begin{split}&\parbox{1084.04999pt}{$\text{log}(R_{W2,\text{ADC}}/400\,\text{light-
day})=(1.12\pm 0.02)\text{log}(R_{W2,\text{non-ADC}}\\\ /400\,\text{light-
day})+(0.02\pm 0.01)$}\end{split}$ (5)

While Lyu et al. (2019) performed a correction for the AD contamination in the
IR light curves, Yang et al. (2020) and Chen et al. (2023) did not take this
corrective approach into account. Consequently, in our study, we implemented
this correction for the W1 and W2- band light curves. The key findings from
this analysis are as follows: (1) correction of AD contamination leads to
larger torus sizes in AGNs, (2) because of the positive correlation between
torus size and AGN luminosity (see Section 5 for details), slopes in Equations
LABEL:adeq1 and LABEL:adeq2 greater than unity suggest that AD contamination
becomes more prominent as luminosity increases, (3) the size ratio
$\mathrm{D_{ratio}}$ ($\mathrm{R_{IR,ADC}/R_{IR,non-ADC}}$) between the size
with/without the AD contamination correction shows a negligible correlation
with the Eddington ratio: Spearman’s rank correlation coefficient
$\mathrm{\rho=-0.003}$ and a p-value of 0.955 ($\mathrm{\rho=-0.041}$,
p=0.492) based on W1 (W2)- band lags, implying no significant dependence on
Eddington ratio, and (4) intercept values for the W1- band are larger than
those for the W2- band, exhibiting a similar slope of $\sim$ 1.12. This
indicates that AD contamination is more significant in the W1- band, aligning
with expectation as AD contaminates more at shorter IR wavelengths. All of
these findings indicate that it is necessary to correct for AD contribution to
the observed IR fluxes.

We compare our obtained redshift and AD contamination corrected torus sizes
with that available in the literature in Section B, finding good agreement
with Lyu et al. (2019) and Chen et al. (2023) (See Figure 20). Table 4
describes the details of the finally selected targets and results of their
lag-analysis.

##  5 Results

###  5.1 The torus size and AGN luminosity correlation

In this section, we investigate the correlation between the measured torus
size based on W1 and W2- band flux variations with AGN luminosity. First we
will use bolometric luminosity of each AGN, and then we will also investigate
the relation using the torus luminosity measured at W2- band. Note that we use
the W1 and W2- band lag measurements after correcting for AD contamination and
redshift effect as described in Sections 4.4, and 4.3, respectively.

We use the following linear equation to parameterize the
$R_{dust,IR}$–$L_{bol}$ relationship

$\text{log}(R_{dust,IR}/1\,\text{light-
day})=\beta\,+\,\alpha\,\text{log}\left(\frac{L_{AGN}}{L_{0}}\right)$ (6)

where $\alpha$ and $\beta$ are the slope and intercept, respectively. $L_{0}$
is the reference point, which is set close to the median of AGN luminosity
$L_{AGN}$ used in the fitting. In fitting Equation 6, we employ the Bivariate
Correlated Errors and intrinsic Scatter (BCES; Akritas & Bershady, 1996). This
method is particularly advantageous as it accommodates measurement errors with
heterogeneity in variance, incorporates intrinsic scatter, and accounts for
correlations within the measurement errors.

Figure 9: Correlation between the torus size and AGN bolometric luminosity,
measured in W1 (left) and W2- (right) bands after correcting for the AD
contamination. The dashed blue line represents the best fit to the data points
and the shaded blue region is the range driven by the uncertainties on the
best fit. The intrinsic scatter is also shown at each panel in the figure. The
dotted-dashed black colored line shows the best-fits for $\alpha_{fix}=0.5$.

####  5.1.1 The torus size correlation with AGN bolometric luminosity

In Figure 9 we present the best-fit relation between torus size and AGN
bolometric luminosity ($L_{bol}$) for W1(left) and W2-(right) bands, where
$L_{0}$ is set to $10^{45.7}\,erg\,s^{-1}$. The torus size scales with the
bolometric luminosity as, $R_{dust,W1}\,\propto\,L_{bol}^{0.39}$, and
$R_{dust,W2}\,\propto\,L_{bol}^{0.33}$ with an intrinsic scatter of 0.19 and
0.18 dex, in W1- and W2- bands, respectively. Note that if we use the lag
measurements without correcting for AD contamination, we obtain similar slopes
of $\alpha$ =0.38 and 0.35 within the uncertainties for W1 and W2- bands,
respectively, but with a smaller intercept ($\beta$) values and slightly
larger intrinsic scatters. This result indicates that the presence of AD
contamination in the observed IR light curves can cause an underestimation of
the torus size. Therefore, correcting for AD contamination is essential to
accurately recover the true torus sizes in AGNs.

To test the potential difference of the results due to the fitting method
(Park et al., 2012), we employ an additional approach: Bayesian linear
regression using the LINMIX$\textunderscore$ERR technique (Kelly, 2007). This
method excels in handling intrinsic scatter within the relationship,
addressing measurement errors in both the independent and dependent variables,
and effectively managing correlations among the measurement errors. We find
that these different fitting methods yield consistent results with a similar
slope and intercept within the uncertainties. We also note that two different
lag analysis methods, i.e., ICCF and JAVELIN, provide consistent results.

Furthermore, depending on the scatter in the observed $R_{dust,IR}$–$L_{bol}$
relation and the varied uncertainties in luminosities and torus sizes, we
perform inverse linear regression fitting, i.e., AGN bolometric luminosity as
a function of torus sizes. The inverse linear fitting results indicate a
relationship: $L_{bol}\,\propto\,R_{W1}^{1.91\pm 0.08}$, and
$L_{bol}\,\propto\,R_{W2}^{2.04\pm 0.11}$ derived from the AD contamination
corrected W1 and W2- band lags, respectively, with a considerably larger
intrinsic scatter of $\sigma\sim 0.4$ dex. From these findings, we derive
dependencies of torus sizes on $L_{bol}$ as $R_{W1}\,\propto\,L_{bol}^{0.52}$
and $R_{W2}\,\propto\,L_{bol}^{0.49}$, closely aligning with the expected
value of 0.5 from the dust-sublimation equilibrium model. It is important to
note that the determination of the dust torus size occurs at a distance from
the accretion disk, where dust particles evaporate because of the central
ionizing luminosity. Hence, the AGN luminosity determines the torus size in
AGN. Consequently, in alignment with our scientific objectives, we assert that
our choice of the forward linear regression fitting between torus size as a
function of AGN luminosity is more appropriate and is consequently employed
for further analysis (refer to Isobe et al. (1990) and Feigelson & Babu (1992)
for the choice of fitting method).

The slopes ($\alpha$) we obtained in our analysis exhibit values ranging from
0.33 to 0.39. These slopes are notably shallower than the expected value of
0.5 predicted by the dust sublimation model. In contrast, Lyu et al. (2019)
and Li & Shen (2023) reported $R_{dust,W1}\,\propto\,L_{bol}^{0.47}$
($R_{dust,W2}\,\propto\,L_{bol}^{0.45}$), and
$R_{dust,W1}\,\propto\,L_{bol}^{0.51}$, respectively, in their findings, which
differ considerably from our results. However, it is noteworthy that Chen et
al. (2023) reported a value of $\alpha=0.37\pm 0.03$ for both W1- and W2-
bands based on the study of DRM in AGNs with $H\beta$ RM data. Our obtained
slopes align with their findings within $1\sigma$ uncertainty. Furthermore,
Minezaki et al. (2019) also reported a consistent slope of 0.42 from DRM based
on K- band dust lags. On the other hand, recent work by Gravity Collaboration
et al. (2023) mentioned a much shallower slope of approximately $\alpha=0.4$
from optical/NIR interferometry. Consequently, the observation of smaller
slope values in the relation between torus size and luminosity in DRM is
consistent with these findings and is not unexpected.

Our results indicate that the slope is likely shallower than the expected
value of 0.5 from the dust radiation equilibrium model. This deviation could
be attributed to the possibility that the optical luminosity at 5100 Å, which
is used to obtain AGN bolometric luminosity, may not accurately represent the
ionizing luminosity (Du & Wang, 2019).

Figure 10: Torus sizes in W1 and W2- bands after correcting for the AD
contamination plotted against the AGN IR luminosity in W2- band in the left
and right panels, respectively. The dashed blue line represents the best fit
to the data points and the shaded blue region is the range driven by the
uncertainties on the best fit. The intrinsic scatter is also shown at each
panel in the figure. The dotted-dashed black colored line shows the best-fits
for $\alpha_{fix}=0.5$.

####  5.1.2 The torus size correlation with AGN IR luminosity

To further investigate the potential discrepancy in the observed
$R_{dust,IR}$–$L_{bol}$ relationship compared to the expected relation, we
obtain a torus size – luminosity relationship ($R_{dust,IR}$–$L_{W2}$) based
on the luminosity of the infrared continuum emission in the W2- band, as it
serves as luminosity indicators for AGNs, less obscured by dust. Thus, these
relations can be applicable for obscured AGNs. We employ Equation 6 to
characterize the relationship between $R_{dust,IR}$ and $L_{W2}$, where
$L_{0}$ is set to $10^{45.0}\,erg\,s^{-1}$. $L_{W2}$ is luminosity measured
with W2- band, which was derived from the average flux values of the W2- band
light curves studied in this work, without correction for host contamination
or galactic extinction. In Figure 10, we present the $R_{dust,IR}$–$L_{W2}$
relationships for AD contamination corrected W1 and W2- band light curves in
the left and right panels, respectively. All the fitting results are
summarized in Table 5.

Interestingly, we observe similar slopes of $\sim$ 0.39 (W1) and 0.33 (W2) in
the $R_{dust,IR}$–$L_{W2}$ relations as those obtained from
$R_{dust,IR}$–$L_{bol}$ within the uncertainties, irrespective of whether we
apply corrections for AD contamination or not. Our obtained slope differs
considerably from the value of 0.5 reported in the torus size – IR luminosity
relation at 12 $\mu m$, as per the K- band DRM measurements involving only 10
AGNs by Minezaki et al. (2019). This discrepancy is likely due to the limited
quantity of data and luminosity range in their study.

The similarity in shallower slopes observed in the torus size – luminosity
relations using both optical and IR luminosities can be attributed to the
self-shadowing effect of the slim disk model (Kawaguchi & Mori, 2011; Wang et
al., 2014a; Du et al., 2018). Since dust particles within the torus reach
thermal equilibrium with the IR radiation re-emitted by dust exposed to
ionizing UV radiation, we would expect a slope of 0.5 in the
$R_{dust,IR}$–$L_{W2}$ relationship (Barvainis, 1987). However, the self-
shadowing effect of slim accretion disks leads to a contraction of the
ionization front of the BLR and the torus, resulting in a reduction in torus
size with increasing accretion rate. This effect yields the shallower slopes
observed in both the $R_{dust,IR}$–$L_{bol}$ and $R_{dust,IR}$–$L_{W2}$
relationships. For further details, please refer to Section 6.4.

= -28mm  Fitting method: BCES ICCF JAVELIN $L_{AGN}$ Case lag $\alpha$ $\beta$
$\sigma$ $\beta_{fix}$($\alpha_{fix}$=0.5) $\alpha$ $\beta$ $\sigma$ (1) (2)
(3) (4) (5) (6) (7) (8) (9) (10) $L_{bol}$ NON ADC optical–W1 0.38 $\pm$ 0.02
2.61 $\pm$ 0.01 0.20 2.66 $\pm$ 0.01 0.34 $\pm$ 0.02 2.62 $\pm$ 0.01 0.20
optical–W2 0.35 $\pm$ 0.02 2.72 $\pm$ 0.01 0.20 2.76 $\pm$ 0.01 0.30 $\pm$
0.02 2.76 $\pm$ 0.01 0.19 $L_{bol}$ ADC optical–W1 0.39 $\pm$ 0.02 2.66 $\pm$
0.01 0.19 2.69 $\pm$ 0.01 0.35 $\pm$ 0.02 2.69 $\pm$ 0.01 0.19 optical–W2 0.33
$\pm$ 0.02 2.76 $\pm$ 0.01 0.18 2.78 $\pm$ 0.01 0.30 $\pm$ 0.02 2.78 $\pm$
0.01 0.19 Fitting method: LINMIX ICCF JAVELIN $L_{AGN}$ Case lag $\alpha$
$\beta$ $\sigma$ $\beta_{fix}$($\alpha_{fix}$=0.5) $\alpha$ $\beta$ $\sigma$
$L_{bol}$ NON ADC optical–W1 0.37 $\pm$ 0.02 2.63 $\pm$ 0.01 0.20 2.68 $\pm$
0.01 0.35 $\pm$ 0.02 2.62 $\pm$ 0.01 0.20 optical–W2 0.34 $\pm$ 0.02 2.74
$\pm$ 0.01 0.20 2.78 $\pm$ 0.01 0.30 $\pm$ 0.02 2.76 $\pm$ 0.01 0.19 $L_{bol}$
ADC optical–W1 0.38 $\pm$ 0.02 2.69 $\pm$ 0.01 0.19 2.72 $\pm$ 0.01 0.35 $\pm$
0.02 2.70 $\pm$ 0.01 0.19 optical–W2 0.32 $\pm$ 0.02 2.78 $\pm$ 0.01 0.18 2.81
$\pm$ 0.01 0.30 $\pm$ 0.02 2.79 $\pm$ 0.01 0.19 Fitting method: BCES ICCF
JAVELIN $L_{AGN}$ Case lag $\alpha$ $\beta$ $\sigma$
$\beta_{fix}$($\alpha_{fix}$=0.5) $\alpha$ $\beta$ $\sigma$ $L_{W2}$ NON ADC
optical–W1 0.40 $\pm$ 0.02 2.67 $\pm$ 0.01 0.20 2.75 $\pm$ 0.01 0.36 $\pm$
0.02 2.67 $\pm$ 0.01 0.19 optical–W2 0.36 $\pm$ 0.02 2.77 $\pm$ 0.01 0.19 2.83
$\pm$ 0.01 0.30 $\pm$ 0.02 2.80 $\pm$ 0.01 0.20 $L_{W2}$ ADC optical–W1 0.39
$\pm$ 0.02 2.73 $\pm$ 0.01 0.19 2.78 $\pm$ 0.01 0.36 $\pm$ 0.02 2.75 $\pm$
0.01 0.19 optical–W2 0.33 $\pm$ 0.02 2.81 $\pm$ 0.01 0.18 2.86 $\pm$ 0.01 0.30
$\pm$ 0.02 2.83 $\pm$ 0.01 0.19 Fitting method: LINMIX ICCF JAVELIN $L_{AGN}$
Case lag $\alpha$ $\beta$ $\sigma$ $\beta_{fix}$($\alpha_{fix}$=0.5) $\alpha$
$\beta$ $\sigma$ $L_{W2}$ NON ADC optical–W1 0.40 $\pm$ 0.02 2.69 $\pm$ 0.01
0.20 2.79 $\pm$ 0.01 0.36 $\pm$ 0.02 2.68 $\pm$ 0.01 0.19 optical–W2 0.35
$\pm$ 0.02 2.79 $\pm$ 0.01 0.19 2.89 $\pm$ 0.01 0.30 $\pm$ 0.02 2.81 $\pm$
0.01 0.20 $L_{W2}$ ADC optical–W1 0.40 $\pm$ 0.02 2.75 $\pm$ 0.01 0.19 2.86
$\pm$ 0.01 0.36 $\pm$ 0.02 2.76 $\pm$ 0.01 0.19 optical–W2 0.33 $\pm$ 0.02
2.81 $\pm$ 0.01 0.18 2.93 $\pm$ 0.01 0.30 $\pm$ 0.02 2.83 $\pm$ 0.01 0.19

Table 5: Linear regression fitting results of the torus size – luminosity
relation.

Note: Columns are (1) luminosity used in fitting, (2) different cases, NON
ADC: without AD contamination correction, ADC: correction for AD contamination
applied, (3) name of the bands used to measure the lags, (4) slope of the
best-fits with $\alpha$ as free parameter using lags from ICCF, (5) intercept
of the best-fits with $\alpha$ as free parameter using lags from ICCF, (6)
intrinsic scatter in dex, (7) intercept of the best-fits with $\alpha$ fixed
to 0.5 using lags from ICCF, (8) slope of the best-fits with $\alpha$ as free
parameter using lags from JAVELIN, (9) intercept of the best-fits with
$\alpha$ as free parameter using lags from JAVELIN, and (10) intrinsic scatter
in dex.

Figure 11: Comparison between the torus sizes in W1 and W2- bands obtained
from ICCF in the rest-frame of the targets based on their AD contamination
corrected IR light curves. The best-fit to the data points is shown by the
dashed brown line, while the dashed black line represents the 1:1 relation.
Figure 12: Comparison between the BLR and torus sizes. The red circular points
represent targets with torus sizes in W1 and W2- bands from this work. The
dashed red lines are the linear fits between the BLR and the torus sizes in W1
and W2- bands, respectively. Measurements from Chen et al. (2023) are shown as
grey points. The dotted black lines represent the best-fits for the total
sample. The intrinsic scatter for the total sample is also mentioned at each
panel.

##  6 Discussion

###  6.1 Comparison between torus sizes in W1 and W2- bands

In Figure 11, we compare torus sizes in the W1 and W2- bands, revealing that
$R_{dust,W2}$ exceeds that in the W1- band. This result aligns with our
expectations, considering that lag generally increases with wavelength.
Moreover, the mean torus size ratio between W2 and W1- bands is found to be
$\sim$ 1.20. A linear fit to the torus sizes between W1 and W2 yields

$\begin{split}&\parbox{650.43pt}{$\text{log}(R_{dust,W2}/400\,\text{light-
day})=(0.90\pm 0.02)\,\text{log}(R_{dust,W1}/\\\ 400\,\text{light-
day})+(0.09\pm 0.01)$}\end{split}$ (7)

and after correcting for the AD contamination, we found

$\begin{split}&\parbox{650.43pt}{$\text{log}(R_{dust,W2}/400\,\text{light-
day})=(0.90\pm 0.02)\,\text{log}(R_{dust,W1}/\\\ 400\,\text{light-
day})+(0.08\pm 0.01)$}\end{split}$ (8)

Indeed, the fitting results indicate that $R_{dust,W2}$ is statistically
larger than that in the W1- band.

Table 6: Details of BLR–RM sample

= -20mm  Object name ID no $R_{H\beta}$ log$L_{5100{\AA}}$ $M_{BH}$ Eddington
ratio Reference (light-days) ($\mathrm{erg\,s^{-1}}$) ($\mathrm{\times
10^{7}M_{\odot}}$) (1) (2) (3) (4) (5) (6) (7) PG 0923+201 OB11
$108.20^{+6.60}_{-12.30}$ $45.31\pm 0.06$ $118.00^{+11.00}_{-16.00}$ $-0.90\pm
0.09$ Li et al. (2021) PG 1001+291 OB16 $37.30^{+6.90}_{-6.00}$ $45.59\pm
0.03$ $3.33^{+0.62}_{-0.54}$ $0.93\pm 0.11$ Li et al. (2021) PG 1001+054 OB17
$65.50^{+5.60}_{-3.90}$ $44.64\pm 0.03^{a}$ $10.73^{+1.15}_{-1.17}$ $-0.52\pm
0.07$ Bao et al. (2022) Mrk 142 OB24 $2.74^{+0.73}_{-0.83}$ $43.54\pm 0.04$
$0.17^{+0.06}_{-0.07}$ $0.18\pm 0.24$ Bentz et al. (2009) Mrk 142 OB24
$7.90^{+1.20}_{-1.10}$ $43.57\pm 0.05$ $0.17^{+0.06}_{-0.07}$ $0.21\pm 0.24$
Du et al. (2014) PG 1048+342 OB27 $36.80^{+2.40}_{-3.40}$ $44.55\pm 0.05^{a}$
$4.44^{+0.31}_{-0.32}$ $-0.23\pm 0.07$ Bao et al. (2022) PG 1115+407 OB31
$54.40^{+10.20}_{-6.70}$ $44.39\pm 0.06$ $7.76^{+2.23}_{-1.95}$ $-0.63\pm
0.18$ Hu et al. (2021) PG 1121+422 OB32 $115.80^{+24.20}_{-20.20}$ $44.83\pm
0.01$ $10.10^{+3.00}_{-2.00}$ $-0.31\pm 0.16$ Woo et al. (2024) Mrk 40 OB33
$3.99^{+0.49}_{-0.68}$ $42.48\pm 0.11$ $0.53^{+0.06}_{-0.07}$ $-1.38\pm 0.13$
Bentz et al. (2009) Mrk 1310 OB41 $3.66^{+0.59}_{-0.61}$ $42.23\pm 0.17$
$0.17^{+0.07}_{-0.07}$ $-1.13\pm 0.30$ Bentz et al. (2009) PG 1202+281 OB43
$38.50^{+9.10}_{-8.50}$ $44.43\pm 0.01$ $12.00^{+3.00}_{-3.00}$ $-0.78\pm
0.15$ Woo et al. (2024) PG 1307+085 OB48 $105.60^{+36.00}_{-46.60}$ $44.79\pm
0.02$ $29.35^{+10.75}_{-13.53}$ $-0.81\pm 0.26$ Kaspi et al. (2000) PG
1404+226 OB55 $21.40^{+2.80}_{-6.50}$ $44.11\pm 0.03$ $0.68^{+0.14}_{-0.23}$
$0.14\pm 0.17$ Hu et al. (2021) PG 1415+451 OB56 $30.00^{+5.40}_{-4.60}$
$44.23\pm 0.04$ $1.75^{+0.36}_{-0.32}$ $-0.15\pm 0.13$ Hu et al. (2021) PG
1427+480 OB59 $33.30^{+20.60}_{-19.40}$ $44.79\pm 0.01$ $5.00^{+3.00}_{-3.00}$
$-0.04\pm 0.37$ Woo et al. (2024) PG 1448+273 OB62 $30.10^{+10.70}_{-5.90}$
$44.24\pm 0.03$ $1.01^{+0.38}_{-0.23}$ $0.10\pm 0.19$ Hu et al. (2021) Mrk
1392 OB65 $26.70^{+3.50}_{-3.90}$ $43.20\pm 0.16$ $6.30^{+0.80}_{-0.90}$
$-0.17\pm 0.18$ U et al. (2022) Mrk 290 OB70 $8.72^{+1.21}_{-1.02}$ $43.11\pm
0.06$ $1.97^{+0.30}_{-0.26}$ $-1.32\pm 0.11$ Denney et al. (2010) Mrk 493 OB75
$11.60^{+1.20}_{-2.60}$ $43.11\pm 0.08$ $0.15^{+0.04}_{-0.15}$ $-0.20\pm 0.46$
Wang et al. (2014b) PG 2349-014 OB83 $55.90^{+11.00}_{-10.80}$ $44.66\pm 0.02$
$63.00^{+14.00}_{-14.00}$ $-1.27\pm 0.14$ Woo et al. (2024)

Note: Columns are (1) object name, (2) object ID number used in this work, (3)
$H\beta$ BLR size (rest-frame) obtained from BLR–RM, (4) host-corrected AGN
optical luminosity at 5100 Å collected from the literature with BLR–RM
measurements, (5) $M_{BH}$ obtained from RM, (6) Eddington ratio in log scale
derived from reverberation-mapped $M_{BH}$, and (7) references.
$L_{5100{\AA}}$ marked with $a$ is not corrected for host-starlight
contamination.

###  6.2 Comparison between torus and BLR sizes

We have identified a total of 19 AGNs within our sample, all of which possess
prior measurements of BLR size and $M_{BH}$ based on $H\beta$ line obtained
through BLR–RM analysis available in the existing literature (See Table 6).

Figure 12 illustrates a comparative analysis between the torus sizes observed
in the W1 and W2- bands and their corresponding $H\beta$ BLR sizes
($R_{H\beta}$) as determined by BLR–RM in the left and right panel,
respectively for those 19 AGNs as shown by the red points. Notably, our
examination reveals a discernible correlation between BLR sizes and our
obtained torus sizes with a Spearman’s rank correlation coefficient of
$\mathrm{\rho=0.73}$ with a p-value of $2.56\times 10^{-4}$ for the W1- band,
and $\mathrm{\rho=0.71}$ with a p-value of $1.04\times 10^{-3}$ for the W2-
band. We employ the following linear equation to parameterize the relationship
between $R_{H\beta}$ and $R_{dust,IR}$

$\text{log}(R_{H\beta}/1\,\text{light-
day})=\epsilon\,\text{log}(R_{dust,IR}/1\,\text{light-day})+K$ (9)

where $\epsilon$ denotes the slope and $K$ is the intercept. We find that the
BLR size scales with the torus size as,
$R_{H\beta}\,\propto\,R_{dust,W1}^{1.42\pm 0.32}$ with $K=2.18\pm 0.88$,
$\sigma=0.37$ dex, and $R_{H\beta}\,\propto\,R_{dust,W2}^{1.42\pm 0.31}$ with
$K=2.29\pm 0.87$, $\sigma=0.39$ dex. The median size ratios between torus and
BLR are approximately $\sim$ 9.5 and 15.3 for W1 and W2- bands, respectively.
In addition, Koshida et al. (2014) reported a significant four to five-fold
disparity in the mean reverberation radii between the innermost dust torus
based on K- band lag and the BLR. A median size ratio between the torus and
BLR, $R_{dust,K}/R_{H\beta}$, was determined to be approximately $\sim$6.5
based on 10 objects (Gandhi et al., 2015). Du et al. (2015) reported
$R_{dust,K}/R_{H\beta}$ $\sim$4 for 10 low-accreting AGNs, and
$R_{dust,K}/R_{H\beta}$ $\sim$7 for high-accreting AGNs with a smaller sample
size (4 objects). Kokubo & Minezaki (2020) also reported
$R_{dust,K}/R_{H\beta}$ $\sim$4 based on measurements from 15 Seyfert
galaxies. Moreover, Mandal et al. (2018) found about 11.5 times larger torus
size (based on K-band) than BLR size in H0507+164. Recently, Chen et al.
(2023) found $R_{dust,K}/R_{H\beta}$ $\sim$6.2, $R_{dust,W1}/R_{H\beta}$
$\sim$9.2, and $R_{dust,W2}/R_{H\beta}$ $\sim$11.2 based on a significantly
larger sample of 78 AGNs with $H\beta$-RM measurements. Our findings indicate
that the ratio of $R_{dust,W1}/R_{H\beta}$ is $\sim$ $0.98\pm 0.36$ dex,
suggesting a torus size roughly 9.5 times larger than that of the BLR, which
closely aligns with the results of Chen et al. (2023). In contrast, we observe
a slightly higher ratio of $R_{dust,W2}/R_{H\beta}$ at around $1.18\pm 0.36$
dex, indicating a torus size approximately 15.1 times larger than the BLR size
based on W2-band lags, derived from a sample of 19 AGNs. The discrepancy
observed in W2-band data may be attributed to the smaller number of AGNs with
$H\beta$-RM data in our sample compared to their sample. Nonetheless, taking
into account the rms scatter of approximately 0.36 dex, our measurements are
consistent with previous findings reported in the literature. It is worth
highlighting a consistent trend observed across these 19 AGNs: the torus sizes
consistently surpass the corresponding BLR sizes.

Since, the above $R_{H\beta}$ – $R_{dust,IR}$ relations are based on limited
number of AGNs, to expand the sample size we incorporate measurements from
Chen et al. (2023), who reported torus sizes for 78 AGNs with $R_{H\beta}$
measurements. Their sample does not cover the S82 region, and only 28 among
the 76 targets from their sample overlap with our initial Sample 1 before
applying the lag-quality assessment. After implementing our lag-quality
assessment (see Section 4.2), only 12 targets remain, that overlap with their
sample. Note that, Chen et al. (2023) did not correct for AD contamination in
their measured lags. A direct comparison of the torus sizes based on W1- (W2)
band lags for those 12 targets between our measurements and that from Chen et
al. (2023) is discussed in Section B, which shows an offset of 0.16 (0.14) dex
between the two sets of measurements. Additionally, we find that, among the
remaining 66 (78-12) targets, 60 have torus size measurements with
$r_{max}\,\geq\,0.6$. Thus, we include these additional 60 targets from Chen
et al. (2023) in the $R_{H\beta}$ – $R_{dust,IR}$ plane represented by the
grey points in Figure 12 and obtain the best fit relations using Equation 9.
We use same redshift correction factor of $(1+z)^{-0.38}$ to get the rest-
frame torus size, while disregarding the extra uncertainties caused by the AD
contamination correction which was not performed for the additional 60
targets. For the total BLR–RM sample (19+60), $R_{H\beta}$ varies with
$R_{dust,IR}$ as, $R_{H\beta}\,\propto\,R_{dust,W1}^{1.20\pm 0.09}$ with
$K=1.44\pm 0.23$, $\sigma=0.34$ dex, and
$R_{H\beta}\,\propto\,R_{dust,W2}^{1.22\pm 0.08}$ with $K=1.65\pm 0.21$,
$\sigma=0.33$ dex.

While Chen et al. (2023) reported a linear relationship between $R_{H\beta}$
and $R_{dust,IR}$, our findings reveal a super-linear relationship
characterized by a slope exceeding unity. This aligns with the findings of
Gravity Collaboration et al. (2023), who also identified a slope greater than
unity of $\sim$ 1.12 and 1.10 in the $R_{H\beta}$ – $R_{dust,IR}$ relationship
using K-band torus sizes derived from RM and optical/NIR interferometry,
respectively. According to the photoionization model, the relationship between
$R_{H\beta}$ and $R_{dust,IR}$ can be expressed as (Chen et al., 2023):

$\frac{R_{H\beta}}{R_{dust,IR}}=\left(\frac{L_{ion}}{L_{OUV}}\right)^{1/2}\left(\frac{\sigma_{s}T_{IR}^{4}}{cn_{e}kT_{c}}\right)^{1/2}$
(10)

where, $L_{ion}$ and $L_{OUV}$ represent the ionizing and optical-UV
luminosity, respectively. $\sigma_{s}$, $n_{e}$, $k$, $T_{c}$, and $T_{IR}$
denote the Stefan-Boltzmann constant, electron number density, Boltzmann
constant, temperature of the ionized clouds, and temperature emitting IR
photons, respectively. The deviation from a linear relationship between
$R_{H\beta}$ and $R_{dust,IR}$ may occur for two primary reasons. First,
within AGNs, BLR gas clouds could exhibit a range of $n_{e}$ (Baldwin et al.,
1995). Second, the optical luminosity measured at 5100 Å, utilized to derive
$L_{bol}$, might not precisely represent the ionizing luminosity, implying
that $L_{bol}$ is not directly proportional to $L_{ion}$ (refer to the
discussion by Czerny et al., 2019; Fonseca Alvarez et al., 2020).
Consequently, the observed super-linear relationship between $R_{H\beta}$ and
$R_{dust,IR}$ is not unexpected.

Moreover, when considering the relationship between torus size in W1- (W2)
band and AGN bolometric luminosity, we find that the BLR size scales with
$L_{bol}$ as $R_{H\beta}\,\propto\,L_{bol}^{0.47}$
($R_{H\beta}\,\propto\,L_{bol}^{0.40}$). Thus, relying on torus sizes in W1
and W2- bands, the slope of $R_{H\beta}$–$L_{bol}$ relation ranges between
0.47 and 0.40 with a mean value of $\sim 0.43$, which is smaller than the
expected value of 0.5 from photoionization model (Davidson, 1972). However, it
closely aligns with the observed slope of 0.402 obtained from BLR–RM, as
reported by Woo et al. (2024) within the uncertainty.

Consequently, these findings imply that the BLR may either constitute an inner
extension of the torus or originate from infalling material associated with
the torus.

Figure 13: Left: $R_{dust}$–$L_{bol}$ relations based on W1, W2 and K-band
lags, which are corrected for the AD contamination. Theoretical estimates are
calculated with Equations 11 and 12, respectively, for silicate with a
sublimation temperature $T_{sub}$ = 1500 K (black dashed line) and for
graphite grains with $T_{sub}$ = 1800 K (black dotted-dashed line). The
H$\beta$ BRL size – luminosity relation adopted from Woo et al. (2024) is also
presented. Right: Torus size measured at K, W1, and W2- band for three
different values of $L_{bol}$ = $10^{44}$, $10^{46}$ and
$10^{48}\,erg\,s^{-1}$, and the best-fit relation for each luminosity case.
The error bars represent the intrinsic scatter in the corresponding torus size
– luminosity relation.

###  6.3 Dependency of torus size on wavelength

The availability of nearly simultaneous monitoring observations across various
IR- bands, coupled with optical data, enables us to determine torus sizes
across distinct IR wavelengths and, in turn, examine whether there exists any
wavelength-dependent trend in torus sizes. However, in our sample, only two
targets have previous DRM measurements based on K-band data in the literature:
ID2149795 (J213818.96$+$011222.3, with W1 and W2) and ID392158
(J232640.01$-$003041.4, with W1). Thus, we cannot establish a one-to-one
relationship between the lags in the K- band and those in the other IR W1 and
W2- bands.

Instead, we present the torus size – AGN bolometric luminosity relations based
on our work using W1 and W2-band, along with the work by Minezaki et al.
(2019) using K-band in Figure 13 (left). All these $R_{dust}$–$L_{bol}$
relations have been adjusted for the AD contamination in their respective IR
(W1, W2 and K- band) light curves. The $R_{dust,K}$–$L_{V}$ relation is based
on K-band lag and V-band luminosity of a sample of 36 quasars with relatively
moderate to low luminosities
($10^{42.7}\,erg\,s^{-1}<L_{bol}<10^{46.4}\,erg\,s^{-1}$), the majority of
which have bolometric luminosities below $10^{46}\,erg\,s^{-1}$. We converted
the V- band luminosity of the $R_{dust,K}$–$L_{V}$ relation to bolometric
luminosity using a correction factor of 10 (Yang et al., 2020). If we consider
the same $R_{dust,K}$–$L_{bol}$ relation to hold for high luminosity AGNs,
Figure 13 (left) illustrates wavelength-dependent torus sizes: $R_{dust,W2}$
$>$ $R_{dust,W1}$ $>$ $R_{dust,K}$ for a given luminosity. We notice that as
luminosity increases, the torus sizes in different IR bands gradually
converge, suggesting that higher luminosity objects tend to exhibit a more
compact torus structure, a notion in alignment with the findings of Kishimoto
et al. (2011b) from interferometric observations.

The innermost radius of the torus composed of silicate grains is defined by
(Barvainis, 1987; Nenkova et al., 2008b)

$R_{dust,Si}\simeq
1.3\left(\frac{L_{bol}}{10^{46}erg\,s^{-1}}\right)^{1/2}\left(\frac{1500K}{T_{sub}}\right)^{2.6}pc,$
(11)

where $T_{sub}$ represents the dust sublimation temperature. Similarly, the
sublimation radius for pure graphite grains is given by (Mor & Netzer, 2012)

$R_{dust,C}\simeq
0.5\left(\frac{L_{bol}}{10^{46}erg\,s^{-1}}\right)^{1/2}\left(\frac{1800K}{T_{sub}}\right)^{2.8}pc.$
(12)

The hot dust spectrum obtained by Mor & Netzer (2012) suggests that the
existence of hot graphite dust makes a substantial contribution to the
observed luminosity at 3.6 and 4.5 $\mu m$. On the other hand, using $T_{sub}$
= 1700 K and a dust grain size of $\sim$ 0.1 $\mu m$, Koshida et al. (2014)
found that the sublimation radius is close to the $K$-band dust reverberation
radius. For comparision we over-plot the theoretical prediction of the
$R_{dust}$–$L_{bol}$ relations as defined by the Equations 11 and 12 for two
distinct sublimation temperatures $T_{sub}$ = 1500 K and 1800 K, respectively.

Using Equations 11 and 12, we estimate the sublimation radius of a low-
luminosity AGN with $L_{bol}=10^{44}\,erg\,s^{-1}$. For silicate dust with
$T_{sub}$ = 1500 K, $R_{dust,Si}$ is 0.13 pc, while for pure graphite with
$T_{sub}$ = 1800 K, $R_{dust,C}$ is 0.05 pc. Whereas, for a high-luminosity
AGN with $L_{bol}=10^{48}\,erg\,s^{-1}$, the sublimation radius is
significantly larger as $\sim$ 13.0 pc and 5.0 pc, respectively, for silicate
with $T_{sub}$ = 1500 K and pure graphite grains with $T_{sub}$ = 1800 K.

We also empirically determine the torus size from the observed
$R_{dust}$–$L_{bol}$ relationships (Equation 6). The derived radii are
$R_{dust,W1}=0.08$ pc ($R_{dust,W2}=0.13$ pc), and $R_{dust,W1}=3.0$ pc
($R_{dust,W2}=2.8$ pc), respectively, for AGNs with $L_{bol}=10^{44}$ and
$10^{48}\,erg\,s^{-1}$. By comparing the theoretical estimates and the
empirically determined torus sizes, we find that for low-luminosity AGNs
($\sim 10^{43.4}-10^{45}\,erg\,s^{-1}$), the primary source of the 3.4 and 4.6
$\mu m$ emission lies within the region bound by the graphite and silicate
sublimation radii, respectively. These results align with the previous
findings with the Spitzer Space Telescope data. For example, Vazquez et al.
(2015) derived the torus size of NGC 6418 ($L_{bol}\geq
10^{43.3}\,erg\,s^{-1}$) as $R_{dust,3.6\mu m}=0.031$ pc with 3.6 $\mu m$
data, and $R_{dust,4.5\mu m}=0.040$ pc with 4.5 $\mu m$ data.

In contrast, for very high luminosity AGNs ($L_{bol}\geq
10^{47}\,erg\,s^{-1}$), W1 and W2-band emissions originate from a similar
region, which is roughly consistent with the sublimation radius of graphite
grains with $T_{sub}=1800$ K (See left panel of Figure 13). This alignment is
well-supported by the computed emission spectrum of graphite dust in the model
by Mor & Netzer (2012), which peaks in the 2–3 $\mu m$ range and exhibits
substantial emission in the 3.4$-$4.6 $\mu m$ range as well. Note that the
observed $R_{dust}$–$L_{bol}$ relations show shallower slopes than expected
from theoretical relations, resulting in a deviation of the observed torus
sizes from theoretical predictions as luminosity increases. In addition, we
present the H$\beta$ BLR size – luminosity relation from Woo et al. (2024),
demonstrating that the torus lies outside the BLR in AGNs.

There are only a few cases with all three (i.e., K, W1, and W2-band) lag
measurements. For example, a wavelength dependent torus size is observed in
ID2149795 (J213818.96$+$011222.3, $L_{bol}=10^{45.69}\,erg\,s^{-1}$) from our
sample as $R_{dust,K}$ = $0.16\pm 0.01$ pc (Minezaki et al., 2019),
$R_{dust,W1}$ = $0.20^{+0.10}_{-0.15}$ pc, and $R_{dust,W2}$ =
$0.35^{+0.10}_{-0.16}$ pc, indicating a stratified torus structure. In
contrast, Z229–15($L_{bol}\sim 10^{44.11}\,erg\,s^{-1}$) shows somewhat
different trend. Guise et al. (2022) reported $R_{dust,3.6\mu m}\sim 0.009$
pc, and $R_{dust,4.5\mu m}\sim 0.023$ pc, while Mandal et al. (2021) reported
K-band lag as $R_{dust,K}\sim 0.017$ pc. Similarly, ID392158
(J232640.01$-$003041.4) with $L_{bol}=10^{45.64}\,erg\,s^{-1}$ from our sample
shows comparable torus sizes in K and W1-bands with $R_{dust,K}=0.19\pm 0.06$
pc (Minezaki et al., 2019) and $R_{dust,W1}=0.18^{+0.11}_{-0.13}$ pc. However,
considering the relatively large uncertainty of the measured lags and the
scatter in the $R_{dust}$–$L_{bol}$ relationships, we expect such exceptions.

By combining the best-fit parameters obtained by Minezaki et al. (2019) using
K- band lags and our measurements in the W1 and W2- bands, we determine torus
sizes from the best-fit torus size – AGN bolometric luminosity relations at K,
W1 and W2- bands at luminosities $L_{bol}$ = $10^{44}$, $10^{46}$, and
$10^{48}\,erg\,s^{-1}$ spanning the entire dynamical range used in this study
(see right panel of Figure 13). At $L_{bol}$ = $10^{46}\,erg\,s^{-1}$, nearly
the median luminosity of our sample, we obtain the ratios of the dust torus
sizes: $R_{dust,K}:R_{dust,W1}:R_{dust,W2}$ = 1.0:1.5:1.8. We find that the
torus size scales with the wavelength as $R_{dust}\,\propto\,\lambda^{0.80}$
at $L_{bol}$ = $10^{46}\,erg\,s^{-1}$, indicating a stratified torus
structure. Note that the dependence of $R_{dust}$ on wavelength becomes weaker
as luminosity increases as indicated by the shallower slopes observed in the
$R_{dust}$–$L_{bol}$ relations.

Extensive efforts, comprising both observational studies (Tomita et al., 2006;
Yoshii et al., 2014; Minezaki et al., 2019) and theoretical investigations
(Barvainis, 1987; Kawaguchi & Mori, 2011), have been dedicated to gain
insights into the structure of the dust torus in AGNs. The presence of a
radial temperature gradient in the dust torus was reported by Tomita et al.
(2006). Additionally, a clumpy dust torus model (Nenkova et al., 2008a;
Kawaguchi & Mori, 2011) proposes a vertical stratification of rotating dust
cloud distributions within galactic nuclei (García-Burillo et al., 2019).
Furthermore, Almeyda et al. (2017, 2020) through simulations of the
reverberation response of a clumpy torus, found that the torus response
exhibits a strong wavelength dependency. This phenomenon arises from the
temperature gradient within the cloud surfaces of the torus and the highly
anisotropic cloud emission at shorter wavelengths. Additionally, our findings
support the idea that different IR emissions originate from distinct regions
of the torus, thus lending support to the notion of a stratified torus
structure.

Figure 14: The deviations of the dust-torus sizes from the corresponding best-
fitting relationships ($\Delta R_{W1}$, and $\Delta R_{W2}$) as a function of
Eddington ratio, $\mathcal{R}_{[OIII]}$, $\mathcal{R}_{FeII}$, and redshift
$z$. The black points represent DRM sample corrected for AD contamination in
the IR W1, W2 light curves. The Spearman’s rank correlation coefficients and
the corresponding probabilities are also mentioned in each panel. The blue
dashed lines represent $\Delta R_{IR}=0$. Figure 15: Same as Figure 14, but
for a fixed slope $\alpha_{fix}=0.5$.

###  6.4 Scatter in the $R_{dust,IR}$–$L_{bol}$

In this section, we explore the potential correlation between the scatter
around the $R_{dust,IR}-L_{bol}$ relations and various other AGN parameters,
including the Eddington ratio, the flux ratios [O iii]/$H\beta$
($\mathcal{R}_{[OIII]}$), Fe ii/$H\beta$ ($\mathcal{R}_{FeII}$), and redshift
$z$. We define the deviations of the measured torus sizes from the
corresponding best-fitting relationships as $\Delta R_{IR}$ =
$\mathrm{log}R_{IR}$ \- $\mathrm{log}R_{best-fit}$.

We compare $\Delta R_{IR}$ with the Eddington ratio in Figure 14 (top left)
considering the slope in the $R_{dust,IR}$–$L_{bol}$ relation as free
parameter. It is important to note that the Eddington ratio in this analysis
is derived from a single-epoch $M_{BH}$ (see Section 4.4) using the BLR size –
luminosity relation (Woo et al., 2024), which is based on AGNs in the
luminosity range $10^{42}<L_{5100{\AA}}<10^{46}\,erg\,s^{-1}$ with the
majority of the AGNs having luminosities below $10^{44.5}\,erg\,s^{-1}$. This
limitation imposes constraints on the precision of $M_{BH}$ measurements at
higher redshift, where high-luminosity AGNs are more commonly observed. As a
result, the canonical single-epoch $M_{BH}$ can lead to substantial over-
estimations and, consequently, underestimations in the Eddington ratio,
especially for AGNs with high accretion rates frequently found in the high-
luminosity regime (also see section C).

Nonetheless, despite these limitations, we identify a mild negative
correlation between $\Delta R_{IR}$ and the Eddington ratio, with Spearman’s
rank correlation coefficients of $\mathrm{\rho=-0.13}$ and a probability of p
= $1.87\times 10^{-2}$ in the W1- band, and $\mathrm{\rho=-0.13}$ with a
p-value of $2.61\times 10^{-2}$ in the W2- band.

We also compare $\Delta R_{IR}$ with the Eddington ratio in Figure 15
considering a fixed slope of $\alpha_{fix}=0.5$ in the $R_{dust,IR}-L_{bol}$
relation as predicted by dust radiation equilibrium theory. This time, we find
a moderate negative correlation between them, with $\mathrm{\rho=-0.23}$, p =
$2.36\times 10^{-5}$ and $\mathrm{\rho=-0.29}$, p = $3.66\times 10^{-7}$ in W1
and W2- bands, respectively.

However, it is important to note that the Eddington ratio scales with
$L/R_{H\beta}$ $\sim$ $L/R_{dust,IR}$ due to the linear correlation present
between $R_{H\beta}$ and $R_{dust,IR}$ (Chen et al., 2023). Further, the
change in $\Delta R_{IR}$ varies as $R_{dust,IR}/L^{0.5}$. As a result, the
negative correlation observed between $\Delta R_{IR}$ and the Eddington ratio
may be attributed to the self-correlation between the Eddington ratios and the
deviation from the best-fit relation (Fonseca Alvarez et al., 2020).

Instead of Eddington ratio, we also utilized $\mathcal{R}_{[OIII]}$ and
$\mathcal{R}_{FeII}$, as indicators of Eddington ratio (Boroson & Green, 1992;
Sulentic et al., 2000; Ludwig et al., 2009; Marziani et al., 2003; Sun & Shen,
2015; Du & Wang, 2019). We collected $H\beta$, [O iii], and Fe ii luminosities
from Liu et al. (2019) for Sample 1. For Sample 2, the $H\beta$, and [O iii]
luminosities were retrieved from Shen et al. (2011), whereas Fe ii
luminosities were obtained from Rakshit et al. (2020). Hence, we present
$\Delta R_{IR}$ as a function of $\mathcal{R}_{[OIII]}$ and
$\mathcal{R}_{FeII}$ in Figure 14, and 15 for varying $\alpha$ and a fixed
$\alpha_{fix}=0.5$, respectively. Notably, when $\alpha$ is treated as a free
parameter, $\Delta R_{IR}$ exhibits relatively weak correlations with
$\mathcal{R}_{[OIII]}$ and $\mathcal{R}_{FeII}$—positive and negative,
respectively. However, these correlations become more pronounced when $\alpha$
is held constant at 0.5, as seen in Figures 14 and 15 for comparison.

In addition, we examine the dependence of $\Delta R_{IR}$ with redshift $z$ to
investigate any residual redshift dependence. We find no significant
dependence when $\alpha$ is treated as a free parameter. This suggests that
the redshift dependence of dust lags is negligible after the $z$ correction
applied to our sample as detailed in Section 4.3. However, a noticeable
correlation between $\Delta R_{IR}$ and $z$ emerges when we use
$\alpha_{fix}=0.5$, indicating a reduction in torus sizes for AGNs at higher
redshifts. However, this is mainly due to the luminosity effect since on
average higher luminosity AGNs are selected at higher $z$.

Our findings are consistent with the previous studies by Chen et al. (2023).
Notably, Yang et al. (2020) reported no dependency of $\Delta R_{W1}$ on
Eddington ratio. Note that they used Eddington ratio measured from single
epoch $M_{BH}$. In contrast, Chen et al. (2023) found that $\Delta R_{W1}$ and
$\Delta R_{W2}$ systematically decrease with increasing accretion rate,
suggesting the shortening of time lags for both BLRs and the dusty torus. They
explained such behaviors by the self-shadowing effect of the slim disk model.
Similarly, Fonseca Alvarez et al. (2020) found that AGNs have a shorter BLR
size with lower $\mathcal{R}_{[OIII]}$. According to them, the deviation from
the radius–luminosity relation depends on the distribution of the UV/optical
SED and the relative contribution of ionizing radiation (for the case of
H$\beta$ BLR, see Woo et al., 2024). All of these findings align with our own
results, implying that lag shortening of both BLRs and torus is caused by the
same physical mechanism, i.e., the accretion rate. Additionally, the weakening
of the overall correlation between $\Delta R_{IR}$ and the Eddington ratio,
$\mathcal{R}_{[OIII]}$, and $\mathcal{R}_{FeII}$ as the slope becomes flatter
suggests that the flattening of the $R_{dust,IR}$–$L_{bol}$ relation from the
dust radiation equilibrium model is influenced by the self-shadowing effect of
the accretion disk as accretion rate increases.

However, the obtained $R_{dust,IR}$–$L_{AGN}$ relations may suffer from
several selection biases, such as, the IR lags may be underestimated for high
luminosity sources because of the limited optical and IR baseline. Yang et al.
(2020) mentioned that IR lags below 7 yr in the observed-frame can be detected
within the base-line, whereas the lag detection probability will decrease
toward longer lags in the observed-frame.

Utilizing high cadence light curve data in optical with strong features
present in both optical and IR light curves, our analysis has allowed us to
detect lags that are shorter than the cadence of the IR light curves ($\sim$
180 days) for a few targets. This lag detection was achieved through the
application of both ICCF and JAVELIN techniques. It is worth noting that
JAVELIN can detect lags shorter than the observing cadence (Li et al., 2019).
However, it is important to consider that the IR lags might be prone to
overestimation when the observed-frame lags are significantly shorter than the
cadence of the IR light curves, for instance, falling below 100 days.

To ensure the reliability of our lag measurements, we employed a comprehensive
validation process in addition to assessing the quality of measured lags based
on $r_{max}$ and $p(r_{max})_{\tau>0}$ values, as described in Section 4.2.
This validation procedure involved a systematic comparison of each lag
measurement by shifting the optical light curves according to the calculated
lag value and overlaying them onto the observed IR light curve. The results of
this comparison demonstrated a consistent agreement between the shifted
optical light curves and the observed IR light curves. As an illustrative
example, please refer to Figure 2 for the target OB76 and Figure 16 for OB45
and OB60, all of which exhibit observed-frame W1- band lags of less than 100
days. These figures provide clear visual evidence of the robustness of our lag
measurements in such cases. The combined effect of these selection biases on
the $R_{dust,IR}$–$L_{bol}$ relation can be demonstrated properly by detail
simulation, which is beyond the scope of this paper.

##  7 Summary

We present DRM results for a total of $\sim$ 446 (416) AGNs with W1- (W2) band
light curves, covering a large dynamic range in luminosities from
$10^{43.4}\,\mathrm{to}\,10^{47.6}\,erg\,s^{-1}$ with redshift $<2$ using
16–20 yr optical light curve data from different ground based optical surveys
and IR data in W1 and W2- bands from WISE during 2010–2020. We explored the
$R_{dust,IR}$–$L_{bol}$ relations for different combinations of data. Our main
findings are summarized below.

  * •

We find that the torus size scales with the AGN bolometric luminosity as
$R_{dust,W1}\propto L_{bol}^{0.38}$ and $R_{dust,W2}\propto L_{bol}^{0.35}$ in
the W1 and W2- bands, respectively. Our best-fit slopes are much shallower
than that expected from the thermal equilibrium model of $R_{dust}\propto
L_{bol}^{0.5}$. However, our results are consistent with that obtained from
the optical/IR interferometric observations, which also predict shallower
slope than expected from the radiation equilibrium model.

  * •

We remove the AD contamination from the observed IR (W1 and W2- bands) fluxes
and obtained the $R_{dust,IR}$–$L_{bol}$ relation. The correction of the AD
contamination in IR fluxes results slightly different slope with
$R_{dust,W1}\propto L_{bol}^{0.39}$ and $R_{dust,W2}\propto L_{bol}^{0.33}$.
As a result of this correction for AD contamination in the IR light curves, we
observe larger torus sizes compared to cases where AD contamination is not
accounted for in the analysis.

  * •

Utilizing IR luminosity in the W2- band as the indicator for AGN luminosity,
we identified analogous slopes in the $R_{dust,IR}$–$L_{W2}$ relationships as
those obtained for $R_{dust,IR}$–$L_{bol}$ associations. This observation
implies that both unobscured and obscured AGNs exhibit a consistent trend
within the torus size – luminosity plane.

  * •

The torus sizes in W1 and W2- bands are found to be strongly correlated with
$R_{dust,W2}$ larger than $R_{dust,W1}$. We also compare our obtained IR
$R_{dust,IR}$–$L_{bol}$ relations in W1 and W2- bands with the
$R_{dust,K}$–$L_{bol}$ relation based on K- band lags available in the
literature. We find a wavelength dependent torus size that increases with
wavelengths as $R_{dust,K}:R_{dust,W1}:R_{dust,W2}$ = 1.0:1.5:1.8
($R_{dust}\,\propto\,\lambda^{0.80}$) at $L_{bol}$ = $10^{46}\,erg\,s^{-1}$.
This suggests a stratified dust torus. Additionally, a noteworthy observation
is that objects with higher luminosities tend to display a more compact torus
structure.

  * •

The sizes of torus are consistently determined to be larger than that of the
BLR, indicating that the dust torus exists outside of the BLR. Our analysis
reveals a super-linear relationship between BLR size and torus size, as
$R_{H\beta}\,\propto\,R_{dust,W1}^{1.20}$ and
$R_{H\beta}\,\propto\,R_{dust,W2}^{1.22}$, demonstrating a consistent slope
observed in the $R_{H\beta}$–$L_{5100{\AA}}$ relationship available in the
literature.

  * •

We find that the scatter observed in the $R_{dust,IR}$–$L_{bol}$ relation
depends on Eddington ratio, while a mild correlation is evident with
$\mathcal{R}_{[OIII]}$ and $\mathcal{R}_{FeII}$. These findings collectively
suggest that as the accretion rate increases, the torus size tends to shorten
because of the self-shadowing effect of the slim disk, and the
$R_{dust,IR}$–$L_{bol}$ relation flattens compared to the predictions from the
dust radiation equilibrium model. Nevertheless, we cannot definitively
establish this correlation, except for the one arising due to the inherent
relationship between AGN luminosity and Eddington ratio.

The uncertainty in the measured lags and the scatter in the obtained torus
size – luminosity relations can be improved using dedicated AGN monitoring
survey, simultaneously, in optical and IR wavelengths in future with proper
cadence. In this case the Vera C. Rubin Observatory Legacy Survey of Space and
Time (LSST; Ivezić et al., 2019) and the Nancy Grace Roman Space Telescope
(NGRST), which was previously known as WFIRST (Green et al., 2012) can play an
important role. The LSST+NGRST data set can be used to measure the size of the
dust torus in local AGNs with redshift $z$ $<$ 0.5 since NGRST can observe up
to $\sim$ 2 $\mu m$ (Yang et al., 2020). In order to employ dust lag as a
standard candle for cosmology and constrain cosmological parameters, Hönig et
al. (2017) started a large DRM program, namely ’VEILS’ (VISTA Extragalactic
Infrared Legacy Survey) that will observe about 1350 Seyfert 1 galaxies in the
redshift range of 0.1 $<$ $z$ $<$ 1.2. The Spectro-Photometer for the History
of the Universe, Epoch of Reionization and Ices Explorer (SPHEREx; Kim et al.,
2021) mission aims to provide optical and NIR multi-epoch spectroscopic data,
which can be utilized in RM studies for bright AGNs. These RM studies can
generate a unique data set that, when combined with complementary optical
observations, can aid in comprehending the physical properties of the dust
torus as well as the central structures of AGNs.

We thank the anonymous referee for valuable comments and suggestions. This
publication uses data products from the WISE, a joint project between the
University of California, Los Angeles, and the Jet Propulsion
Laboratory/California Institute of Technology, and is funded by the National
Aeronautics and Space Administration. The publication also leverages data
products from NEOWISE, a Jet Propulsion Laboratory/California Institute of
Technology project, funded by the Planetary Science Division of the National
Aeronautics and Space Administration. The CRTS survey is conducted by the U.S.
National Science Foundation under grants AST-0909182 and AST-1313422. The
ASAS-SN is supported by the Gordon and Betty Moore Foundation through grant
GBMF5490 to the Ohio State University and NSF grant AST-1515927. The
Development of ASAS-SN was supported by the NSF grant AST- 0908816, the Mt.
Cuba Astronomical Foundation, the Center for Cosmology and AstroParticle
Physics at the Ohio State University, the Chinese Academy of Science South
America Center for Astronomy (CASSACA), the Villum Foundation, and George
Skestos. ZTF is funded by the National Science Foundation through grant number
AST-1440341 and a partnership with Caltech, IPAC, the Weizmann Institute for
Science, the Oskar Klein Center at Stockholm University, the University of
Maryland, the University of Washington, Deutsches Elektronen-Synchrotron and
Humboldt University, Los Alamos National Laboratories, the TANGO Consortium of
Taiwan, the University of Wisconsin at Mil- waukee, and Lawrence Berkeley
National Laboratories. COO, IPAC, and UW manage operations. This paper also
used data observed with the Samuel Oschin Telescope at the Palomar Observatory
as part of the Palomar Transient Factory project, a scientific collaboration
between the California Institute of Technology, Columbia University, Las
Cumbres Observatory, the Lawrence Berkeley National Laboratory, the National
Energy Research Scientific Computing Center, the University of Oxford, and the
Weizmann Institute of Science. This work has been supported by the Basic
Science Research Program through the National Research Foundation of the
Korean Government (grant No. NRF-2021R1A2C3008486) and Samsung Science and
Technology Foundation under Project Number SSTF –BA1501–05. S.W. acknowledges
the support from the National Research Foundation of Korea (NRF) grant funded
by the Korean government (MEST) (No. 2019R1A6A1A10073437). Software: JAVELIN
(Zu et al., 2011), PyCALI (Li et al., 2014)
Facilities: WISE, NEOWISE, CRTS, ASAS-SN, PTF, ZTF

##  Appendix A Examples of light curves and their lag measurements

We present examples of the light curves from a selection of 12 randomly chosen
targets in Sample 1 and Sample 2, showcasing their wide range of bolometric
luminosities ($10^{44.3}-10^{47.4}\,erg\,s^{-1}$). These examples, displayed
in Figures 16 and 17, include the light curves in optical and IR W1 and W2-
bands for Sample 1 and Sample 2, respectively. They serve as representative
illustrations of our total sample. Furthermore, the lag measurements analyses
between the optical and W1- band light curves for these targets, utilizing
ICCF and JAVELIN, are depicted in Figure 18 and 19.

Figure 16: Examples of light curves from Sample 1 in optical and IR W1- and
W2- bands for the targets (a) OB01, (b) OB14, (c) OB32, (d) OB45, (e) OB60,
and (f) OB68. The IR W1- and W2- band light curves represented in the plot are
not corrected for the AD contamination. Same notations are used as Figure 2.

Figure 17: Examples of light curves from Sample 2 in optical and IR W1- and
W2- bands for the targets (a) ID47417, (b) ID1124174, (c) ID1189298, (d)
ID1524735, (e) ID2023374, and (f) ID4065339. The IR W1- and W2- band light
curves represented in the plot are not corrected for the AD contamination.
Same notations are used as Figure 2.

Figure 18: The ICCFs and ${\tt JAVELIN}$ analysis between optical and W1- band
light curves of the example targets (a) OB01, (b) OB14, (c) OB32, (d) OB45,
(e) OB60, and (f) OB68 from Sample 1. Same notations are used as Figure 3.

Figure 19: The ICCFs and ${\tt JAVELIN}$ analysis between optical and W1- band
light curves of the example targets (a) ID47417, (b) ID1124174, (c) ID1189298,
(d) ID1524735, (e) ID2023374, and (f) ID4065339 from Sample 2. Same notations
are used as Figure 3.

##  Appendix B Comparison of our lag measurements with earlier work

In the context of our study, it is important to note that Lyu et al. (2019)
previously reported lag measurements for 67 PG quasars. We identified 18
targets from our selected targets that overlap with the sample reported by Lyu
et al. (2019). To assess the consistency between our new measurements and
those obtained by Lyu et al. (2019) using their independent $\chi^{2}-$fit
method, we compare the torus sizes in W1- band in Figure 20 (top left). With
the exception of two targets, OB17 and OB32, our obtained torus sizes exhibit
consistency (within 2$\sigma$) with those reported by Lyu et al. (2019),
displaying an offset of 0.07 dex.

In Figure 16, we present the light curves of one of the most deviant targets,
OB32, along with the corresponding results from the ICCF analysis shown in
Figure 18. Notably, the optical light curve of OB32, shifted by our recovered
time lag, aligns well with the W1-band light curve. For determining AGN
bolometric luminosity, Lyu et al. (2019) employed optical-to-IR SED fitting
and adopted an uncertainty of 0.3 dex in the luminosity. Figure 20 (top right)
demonstrates the comparable bolometric luminosities we obtained for these 18
targets, in line with those derived from optical-to-IR SED fitting by Lyu et
al. (2019).

In addition, from our analysis, we identified a subset of 12 targets that are
common to both our sample and the sample reported by Chen et al. (2023). They
used MICA (Li et al., 2016) in their time-lag measurements. Despite the robust
agreement observed between the measurements obtained through our approach and
those reported by Chen et al. (2023) (See the middle panels of Figure 20), our
derived torus sizes consistently appear larger with an offset of 0.16 (0.14)
dex for the W1- (W2) bands. This discrepancy arises from the fact that Chen et
al. (2023) did not perform an AD contamination correction in their analysis,
whereas our obtained torus sizes are corrected for AD contamination in the W1
and W2- band light curves.

Finally, we compare our measured torus sizes which are corrected for the AD
contamination with those reported by Yang et al. (2020), who used JAVELIN
without employing AD contamination correction for Sample 2. Despite most of
our measurements lying above the 1:1 line due to the AD contamination
correction, the two sets of measurements exhibit good consistency, showing an
offset of 0.06 dex and a scatter of 0.15 dex.

Note that to ensure consistency in our comparison, we obtained the
corresponding rest-frame time lags by multiplying the estimated observed-frame
lags both from the literature and our own analysis, by the redshift correction
factor (1+$z$)-0.38.

Figure 20: Top panel: (Left) Comparison of the torus sizes in W1- band between
our measurements and that obtained by Lyu et al. (2019). These two
measurements are found to be consistent with each other, except for two
targets, namely OB17 and OB32, shown by individually marking them in the
figure. (Right) Comparison of our obtained bolometric luminosities with that
from Lyu et al. (2019) who performed optical-to-IR SED fitting to derive their
bolometric luminosities. Middle panel: a comparison of the torus sizes in both
the W1 and W2- bands of our measurements (AD contamination corrected) with the
findings of Chen et al. (2023) (AD contamination not corrected). Bottom panel:
a comparison of the torus sizes in the W1-band between our measurements,
corrected for AD contamination, and those from Yang et al. (2020) without AD
contamination correction. The black dashed lines represent the 1:1 relation.
The offset values, along with their rms scatters as $1\sigma$ uncertainties
are also shown at each panel in the figure.

##  Appendix C Comparison between Eddington ratios derived from single epoch
and reverberation mapped $M_{BH}$

We calculated the Eddington ratio for our sample by employing the single epoch
$M_{BH}$ estimation method outlined in Section 4.4. In Figure 21, we present a
comparative analysis of the Eddington ratios derived from the single epoch
$M_{BH}$ estimates in our study (represented by red circles) and those
obtained by Liu et al. (2019) (depicted by blue circles) with those obtained
by using $M_{BH}$ measurements from BLR–RM. Since the black hole mass derived
from RM is considered to be more accurate, Figure 21 clearly illustrates that
the single epoch Eddington ratios tend to be underestimated as the accretion
rate increases. This observation provides evidence for the need to account for
the accurate measurement of $M_{BH}$ from BLR–RM to derive Eddington ratio.

Figure 21: Comparison between the Eddington ratios obtained from single epoch
$M_{BH}$ and $M_{BH}$ from BLR–RM. The red and blue circles represent the
objects with $\mathrm{Eddington\,ratio_{single\,epoch}}$ calculated in this
work (see Section 4.4) and retrieved from Liu et al. (2019), respectively,
with the linear-fits shown by the dotted lines. The black dashed line is the
1:1 relation.

## References

  * Abbott et al. (2018) Abbott, T. M. C., Abdalla, F. B., Allam, S., et al. 2018, ApJS, 239, 18, doi: 10.3847/1538-4365/aae9f0
  * Akritas & Bershady (1996) Akritas, M. G., & Bershady, M. A. 1996, ApJ, 470, 706, doi: 10.1086/177901
  * Almeyda et al. (2020) Almeyda, T., Robinson, A., Richmond, M., Nikutta, R., & McDonough, B. 2020, ApJ, 891, 26, doi: 10.3847/1538-4357/ab6aa1
  * Almeyda et al. (2017) Almeyda, T., Robinson, A., Richmond, M., Vazquez, B., & Nikutta, R. 2017, ApJ, 843, 3, doi: 10.3847/1538-4357/aa7687
  * Baldwin et al. (1995) Baldwin, J., Ferland, G., Korista, K., & Verner, D. 1995, ApJ, 455, L119, doi: 10.1086/309827
  * Bao et al. (2022) Bao, D.-W., Brotherton, M. S., Du, P., et al. 2022, ApJS, 262, 14, doi: 10.3847/1538-4365/ac7beb
  * Barvainis (1987) Barvainis, R. 1987, ApJ, 320, 537, doi: 10.1086/165571
  * Bellm et al. (2019) Bellm, E. C., Kulkarni, S. R., Graham, M. J., et al. 2019, PASP, 131, 018002, doi: 10.1088/1538-3873/aaecbe
  * Bentz et al. (2009) Bentz, M. C., Walsh, J. L., Barth, A. J., et al. 2009, ApJ, 705, 199, doi: 10.1088/0004-637X/705/1/199
  * Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393, doi: 10.1051/aas:1996164
  * Bessell et al. (1998) Bessell, M. S., Castelli, F., & Plez, B. 1998, A&A, 333, 231
  * Blandford & McKee (1982) Blandford, R. D., & McKee, C. F. 1982, ApJ, 255, 419, doi: 10.1086/159843
  * Boroson & Green (1992) Boroson, T. A., & Green, R. F. 1992, ApJS, 80, 109, doi: 10.1086/191661
  * Burtscher et al. (2013) Burtscher, L., Meisenheimer, K., Tristram, K. R. W., et al. 2013, A&A, 558, A149, doi: 10.1051/0004-6361/201321890
  * Chambers et al. (2016) Chambers, K. C., Magnier, E. A., Metcalfe, N., et al. 2016, arXiv e-prints, arXiv:1612.05560, doi: 10.48550/arXiv.1612.05560
  * Chen et al. (2023) Chen, Y.-J., Liu, J.-R., Zhai, S., et al. 2023, MNRAS, doi: 10.1093/mnras/stad1136
  * Czerny et al. (2019) Czerny, B., Wang, J.-M., Du, P., et al. 2019, ApJ, 870, 84, doi: 10.3847/1538-4357/aaf396
  * Davidson (1972) Davidson, K. 1972, ApJ, 171, 213, doi: 10.1086/151274
  * Denney et al. (2010) Denney, K. D., Peterson, B. M., Pogge, R. W., et al. 2010, ApJ, 721, 715, doi: 10.1088/0004-637X/721/1/715
  * Drake et al. (2009) Drake, A. J., Djorgovski, S. G., Mahabal, A., et al. 2009, ApJ, 696, 870, doi: 10.1088/0004-637X/696/1/870
  * Du & Wang (2019) Du, P., & Wang, J.-M. 2019, ApJ, 886, 42, doi: 10.3847/1538-4357/ab4908
  * 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
  * Du et al. (2015) —. 2015, ApJ, 806, 22, doi: 10.1088/0004-637X/806/1/22
  * Du et al. (2018) Du, P., Zhang, Z.-X., Wang, K., et al. 2018, ApJ, 856, 6, doi: 10.3847/1538-4357/aaae6b
  * Elvis et al. (1994) Elvis, M., Wilkes, B. J., McDowell, J. C., et al. 1994, ApJS, 95, 1, doi: 10.1086/192093
  * Feigelson & Babu (1992) Feigelson, E. D., & Babu, G. J. 1992, ApJ, 397, 55, doi: 10.1086/171766
  * Fonseca Alvarez et al. (2020) Fonseca Alvarez, G., Trump, J. R., Homayouni, Y., et al. 2020, ApJ, 899, 73, doi: 10.3847/1538-4357/aba001
  * Fukugita et al. (1996) Fukugita, M., Ichikawa, T., Gunn, J. E., et al. 1996, AJ, 111, 1748, doi: 10.1086/117915
  * Gandhi et al. (2015) Gandhi, P., Hönig, S. F., & Kishimoto, M. 2015, ApJ, 812, 113, doi: 10.1088/0004-637X/812/2/113
  * García-Burillo et al. (2019) García-Burillo, S., Combes, F., Ramos Almeida, C., et al. 2019, A&A, 632, A61, doi: 10.1051/0004-6361/201936606
  * Gaskell & Peterson (1987) Gaskell, C. M., & Peterson, B. M. 1987, ApJS, 65, 1, doi: 10.1086/191216
  * Gaskell & Sparke (1986) Gaskell, C. M., & Sparke, L. S. 1986, ApJ, 305, 175, doi: 10.1086/164238
  * Gravity Collaboration et al. (2020) Gravity Collaboration, Dexter, J., Shangguan, J., et al. 2020, A&A, 635, A92, doi: 10.1051/0004-6361/201936767
  * Gravity Collaboration et al. (2023) Gravity Collaboration, Amorim, A., Bourdarot, G., et al. 2023, A&A, 669, A14, doi: 10.1051/0004-6361/202244655
  * Green et al. (2012) Green, J., Schechter, P., Baltay, C., et al. 2012, arXiv e-prints, arXiv:1208.4012, doi: 10.48550/arXiv.1208.4012
  * Grier et al. (2017) Grier, C. J., Trump, J. R., Shen, Y., et al. 2017, ApJ, 851, 21, doi: 10.3847/1538-4357/aa98dc
  * Grier et al. (2019) Grier, C. J., Shen, Y., Horne, K., et al. 2019, ApJ, 887, 38, doi: 10.3847/1538-4357/ab4ea5
  * Guise et al. (2022) Guise, E., Hönig, S. F., Gorjian, V., et al. 2022, MNRAS, 516, 4898, doi: 10.1093/mnras/stac2529
  * Guo & Barth (2021) Guo, H., & Barth, A. J. 2021, in American Astronomical Society Meeting Abstracts, Vol. 53, American Astronomical Society Meeting Abstracts, 226.08
  * Guo et al. (2022) Guo, H., Barth, A. J., & Wang, S. 2022, ApJ, 940, 20, doi: 10.3847/1538-4357/ac96ec
  * Homayouni et al. (2020) Homayouni, Y., Trump, J. R., Grier, C. J., et al. 2020, ApJ, 901, 55, doi: 10.3847/1538-4357/ababa9
  * Hönig et al. (2017) Hönig, S. F., Watson, D., Kishimoto, M., et al. 2017, MNRAS, 464, 1693, doi: 10.1093/mnras/stw2484
  * Hovatta et al. (2014) Hovatta, T., Pavlidou, V., King, O. G., et al. 2014, MNRAS, 439, 690, doi: 10.1093/mnras/stt2494
  * Hu et al. (2021) Hu, C., Li, S.-S., Yang, S., et al. 2021, ApJS, 253, 20, doi: 10.3847/1538-4365/abd774
  * Isobe et al. (1990) Isobe, T., Feigelson, E. D., Akritas, M. G., & Babu, G. J. 1990, ApJ, 364, 104, doi: 10.1086/169390
  * Ivezić et al. (2019) Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al. 2019, ApJ, 873, 111, doi: 10.3847/1538-4357/ab042c
  * Jarrett et al. (2011) Jarrett, T. H., Cohen, M., Masci, F., et al. 2011, ApJ, 735, 112, doi: 10.1088/0004-637X/735/2/112
  * Jayasinghe et al. (2019) Jayasinghe, T., Stanek, K. Z., Kochanek, C. S., et al. 2019, MNRAS, 485, 961, doi: 10.1093/mnras/stz444
  * Jha et al. (2022) Jha, V. K., Joshi, R., Chand, H., et al. 2022, MNRAS, 511, 3005, doi: 10.1093/mnras/stac109
  * Jiang et al. (2017) Jiang, Y.-F., Green, P. J., Greene, J. E., et al. 2017, ApJ, 836, 186, doi: 10.3847/1538-4357/aa5b91
  * Kaspi et al. (2000) Kaspi, S., Smith, P. S., Netzer, H., et al. 2000, ApJ, 533, 631, doi: 10.1086/308704
  * Kawaguchi & Mori (2011) Kawaguchi, T., & Mori, M. 2011, ApJ, 737, 105, doi: 10.1088/0004-637X/737/2/105
  * Kelly (2007) Kelly, B. C. 2007, ApJ, 665, 1489, doi: 10.1086/519947
  * Kelly et al. (2009) Kelly, B. C., Bechtold, J., & Siemiginowska, A. 2009, ApJ, 698, 895, doi: 10.1088/0004-637X/698/1/895
  * Kim et al. (2021) Kim, M., Jeong, W.-S., Yang, Y., et al. 2021, Journal of Korean Astronomical Society, 54, 37, doi: 10.5303/JKAS.2021.54.2.37
  * Kishimoto et al. (2008) Kishimoto, M., Antonucci, R., Blaes, O., et al. 2008, Nature, 454, 492, doi: 10.1038/nature07114
  * 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. (2009) —. 2009, A&A, 507, L57, doi: 10.1051/0004-6361/200913512
  * Kishimoto et al. (2011b) —. 2011b, A&A, 536, A78, doi: 10.1051/0004-6361/201117367
  * Kishimoto et al. (2007) Kishimoto, M., Hönig, S. F., Beckert, T., & Weigelt, G. 2007, A&A, 476, 713, doi: 10.1051/0004-6361:20077911
  * Kokubo & Minezaki (2020) Kokubo, M., & Minezaki, T. 2020, MNRAS, 491, 4615, doi: 10.1093/mnras/stz3397
  * Koshida et al. (2009) Koshida, S., Yoshii, Y., Kobayashi, Y., et al. 2009, ApJ, 700, L109, doi: 10.1088/0004-637X/700/2/L109
  * Koshida et al. (2014) Koshida, S., Minezaki, T., Yoshii, Y., et al. 2014, ApJ, 788, 159, doi: 10.1088/0004-637X/788/2/159
  * Law et al. (2009) Law, N. M., Kulkarni, S. R., Dekany, R. G., et al. 2009, PASP, 121, 1395, doi: 10.1086/648598
  * Li & Shen (2023) Li, J., & Shen, Y. 2023, arXiv e-prints, arXiv:2302.12437, doi: 10.48550/arXiv.2302.12437
  * Li et al. (2019) Li, I-Hsiu, J., Shen, Y., Brandt, W. N., et al. 2019, ApJ, 884, 119, doi: 10.3847/1538-4357/ab41fb
  * Li et al. (2021) Li, S.-S., Yang, S., Yang, Z.-X., et al. 2021, ApJ, 920, 9, doi: 10.3847/1538-4357/ac116e
  * Li et al. (2016) Li, Y.-R., Wang, J.-M., & Bai, J.-M. 2016, ApJ, 831, 206, doi: 10.3847/0004-637X/831/2/206
  * Li et al. (2014) Li, Y.-R., Wang, J.-M., Hu, C., Du, P., & Bai, J.-M. 2014, ApJ, 786, L6, doi: 10.1088/2041-8205/786/1/L6
  * Lira et al. (2011) Lira, P., Arévalo, P., Uttley, P., McHardy, I., & Breedt, E. 2011, MNRAS, 415, 1290, doi: 10.1111/j.1365-2966.2011.18774.x
  * Liu et al. (2019) Liu, H.-Y., Liu, W.-J., Dong, X.-B., et al. 2019, ApJS, 243, 21, doi: 10.3847/1538-4365/ab298b
  * Ludwig et al. (2009) Ludwig, R. R., Wills, B., Greene, J. E., & Robinson, E. L. 2009, The Astrophysical Journal, 706, 995, doi: 10.1088/0004-637X/706/2/995
  * Lyu et al. (2019) Lyu, J., Rieke, G. H., & Smith, P. S. 2019, ApJ, 886, 33, doi: 10.3847/1538-4357/ab481d
  * MacLeod et al. (2010) MacLeod, C. L., Ivezić, Ž., Kochanek, C. S., et al. 2010, ApJ, 721, 1014, doi: 10.1088/0004-637X/721/2/1014
  * MacLeod et al. (2012) MacLeod, C. L., Ivezić, Ž., Sesar, B., et al. 2012, ApJ, 753, 106, doi: 10.1088/0004-637X/753/2/106
  * Mainzer et al. (2014) Mainzer, A., Bauer, J., Cutri, R. M., et al. 2014, ApJ, 792, 30, doi: 10.1088/0004-637X/792/1/30
  * Mandal et al. (2018) Mandal, A. K., Rakshit, S., Kurian, K. S., et al. 2018, MNRAS, 475, 5330, doi: 10.1093/mnras/sty200
  * Mandal et al. (2021) Mandal, A. K., Rakshit, S., Stalin, C. S., et al. 2021, MNRAS, 501, 3905, doi: 10.1093/mnras/staa3828
  * Marziani et al. (2003) Marziani, P., Zamanov, R. K., Sulentic, J. W., & Calvani, M. 2003, MNRAS, 345, 1133, doi: 10.1046/j.1365-2966.2003.07033.x
  * Masci et al. (2019) Masci, F. J., Laher, R. R., Rusholme, B., et al. 2019, PASP, 131, 018003, doi: 10.1088/1538-3873/aae8ac
  * Minezaki et al. (2019) Minezaki, T., Yoshii, Y., Kobayashi, Y., et al. 2019, ApJ, 886, 150, doi: 10.3847/1538-4357/ab4f7b
  * Mor & Netzer (2012) Mor, R., & Netzer, H. 2012, MNRAS, 420, 526, doi: 10.1111/j.1365-2966.2011.20060.x
  * Mudd et al. (2018) Mudd, D., Martini, P., Zu, Y., et al. 2018, ApJ, 862, 123, doi: 10.3847/1538-4357/aac9bb
  * Nenkova et al. (2008a) Nenkova, M., Sirocky, M. M., Ivezić, Ž., & Elitzur, M. 2008a, ApJ, 685, 147, doi: 10.1086/590482
  * Nenkova et al. (2008b) Nenkova, M., Sirocky, M. M., Nikutta, R., Ivezić, Ž., & Elitzur, M. 2008b, ApJ, 685, 160, doi: 10.1086/590483
  * Oknyanskij & Horne (2001) Oknyanskij, V. L., & Horne, K. 2001, in Astronomical Society of the Pacific Conference Series, Vol. 224, Probing the Physics of Active Galactic Nuclei, ed. B. M. Peterson, R. W. Pogge, & R. S. Polidan, 149
  * Park et al. (2012) Park, D., Kelly, B. C., Woo, J.-H., & Treu, T. 2012, ApJS, 203, 6, doi: 10.1088/0067-0049/203/1/6
  * Peterson (1993) Peterson, B. M. 1993, PASP, 105, 247, doi: 10.1086/133140
  * Peterson et al. (1998) Peterson, B. M., Wanders, I., Horne, K., et al. 1998, PASP, 110, 660, doi: 10.1086/316177
  * Peterson et al. (2004) Peterson, B. M., Ferrarese, L., Gilbert, K. M., et al. 2004, ApJ, 613, 682, doi: 10.1086/423269
  * Rakshit et al. (2020) Rakshit, S., Stalin, C. S., & Kotilainen, J. 2020, ApJS, 249, 17, doi: 10.3847/1538-4365/ab99c5
  * Rees (1984) Rees, M. J. 1984, ARA&A, 22, 471, doi: 10.1146/annurev.aa.22.090184.002351
  * Richards et al. (2006) Richards, G. T., Lacy, M., Storrie-Lombardi, L. J., et al. 2006, ApJS, 166, 470, doi: 10.1086/506525
  * Salpeter (1964) Salpeter, E. E. 1964, ApJ, 140, 796, doi: 10.1086/147973
  * Sanders et al. (1989) Sanders, D. B., Phinney, E. S., Neugebauer, G., Soifer, B. T., & Matthews, K. 1989, ApJ, 347, 29, doi: 10.1086/168094
  * Seyfert (1943) Seyfert, C. K. 1943, ApJ, 97, 28, doi: 10.1086/144488
  * Shakura & Sunyaev (1973) Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 500, 33
  * Shappee et al. (2014) Shappee, B., Prieto, J., Stanek, K. Z., et al. 2014, in American Astronomical Society Meeting Abstracts, Vol. 223, American Astronomical Society Meeting Abstracts #223, 236.03
  * Shen et al. (2011) Shen, Y., Richards, G. T., Strauss, M. A., et al. 2011, ApJS, 194, 45, doi: 10.1088/0067-0049/194/2/45
  * Sitko et al. (1993) Sitko, M. L., Sitko, A. K., Siemiginowska, A., & Szczerba, R. 1993, ApJ, 409, 139, doi: 10.1086/172649
  * Sobrino Figaredo et al. (2020) Sobrino Figaredo, C., Haas, M., Ramolla, M., et al. 2020, AJ, 159, 259, doi: 10.3847/1538-3881/ab89b1
  * Suganuma et al. (2006) Suganuma, M., Yoshii, Y., Kobayashi, Y., et al. 2006, ApJ, 639, 46, doi: 10.1086/499326
  * Sulentic et al. (2000) Sulentic, J. W., Zwitter, T., Marziani, P., & Dultzin-Hacyan, D. 2000, ApJ, 536, L5, doi: 10.1086/312717
  * Sun & Shen (2015) Sun, J., & Shen, Y. 2015, The Astrophysical Journal Letters, 804, L15, doi: 10.1088/2041-8205/804/1/L15
  * Tomita et al. (2006) Tomita, H., Yoshii, Y., Kobayashi, Y., et al. 2006, ApJ, 652, L13, doi: 10.1086/509878
  * U et al. (2022) U, V., Barth, A. J., Vogler, H. A., et al. 2022, ApJ, 925, 52, doi: 10.3847/1538-4357/ac3d26
  * Urry & Padovani (1995) Urry, C. M., & Padovani, P. 1995, PASP, 107, 803, doi: 10.1086/133630
  * Vazquez et al. (2015) Vazquez, B., Galianni, P., Richmond, M., et al. 2015, ApJ, 801, 127, doi: 10.1088/0004-637X/801/2/127
  * Wandel et al. (1999) Wandel, A., Peterson, B. M., & Malkan, M. A. 1999, ApJ, 526, 579, doi: 10.1086/308017
  * Wang et al. (2014a) Wang, J.-M., Qiu, J., Du, P., & Ho, L. C. 2014a, ApJ, 797, 65, doi: 10.1088/0004-637X/797/1/65
  * Wang et al. (2014b) Wang, J.-M., Du, P., Hu, C., et al. 2014b, ApJ, 793, 108, doi: 10.1088/0004-637X/793/2/108
  * Woo & Urry (2002) Woo, J.-H., & Urry, C. M. 2002, ApJ, 579, 530, doi: 10.1086/342878
  * Woo et al. (2015) Woo, J.-H., Yoon, Y., Park, S., Park, D., & Kim, S. C. 2015, ApJ, 801, 38, doi: 10.1088/0004-637X/801/1/38
  * Woo et al. (2024) Woo, J.-H., Wang, S., Rakshit, S., et al. 2024, ApJ, 962, 67, doi: 10.3847/1538-4357/ad132f
  * Wright et al. (2010) Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868, doi: 10.1088/0004-6256/140/6/1868
  * Yang et al. (2020) Yang, Q., Shen, Y., Liu, X., et al. 2020, ApJ, 900, 58, doi: 10.3847/1538-4357/aba59b
  * Yoshii et al. (2014) Yoshii, Y., Kobayashi, Y., Minezaki, T., Koshida, S., & Peterson, B. A. 2014, ApJ, 784, L11, doi: 10.1088/2041-8205/784/1/L11
  * Zu et al. (2011) Zu, Y., Kochanek, C. S., & Peterson, B. M. 2011, ApJ, 735, 80, doi: 10.1088/0004-637X/735/2/80