The Seoul National University AGN Monitoring Project. IV. Hα Reverberation Mapping of Six AGNs and the Hα Size–Luminosity Relation

The broad-line region (BLR) size–luminosity relation has paramount importance for estimating the mass of black holes in active galactic nuclei (AGNs). Traditionally, the size of the Hβ BLR is often estimated from the optical continuum luminosity at 5100 Å, while the size of the Hα BLR and its correlation with the luminosity is much less constrained. As a part of the Seoul National University AGN Monitoring Project, which provides 6 yr photometric and spectroscopic monitoring data, we present our measurements of the Hα lags of high-luminosity AGNs. Combined with the measurements for 42 AGNs from the literature, we derive the size–luminosity relations of the Hα BLR against the broad Hα and 5100 Å continuum luminosities. We find the slope of the relations to be 0.61 ± 0.04 and 0.59 ± 0.04, respectively, which are consistent with the Hβ size–luminosity relation. Moreover, we find a linear relation between the 5100 Å continuum luminosity and the broad Hα luminosity across 7 orders of magnitude. Using these results, we propose a new virial mass estimator based on the Hα broad emission line, finding that the previous mass estimates based on scaling relations in the literature are overestimated by up to 0.7 dex at masses lower than 107 M ⊙.


INTRODUCTION
Mass is the most important physical property of a black hole that we can measure.While the mass of a Corresponding author: Jong-Hak Woo woo@astro.snu.ac.kr black hole can be determined by measuring its gravitational radius, this method thus far has been applied to only two black holes (Event Horizon Telescope Collaboration et al. 2019, 2022).For most extragalactic black holes, the mass is instead measured by observing the kinematics of orbiting bodies near the black hole.While the mass of several black holes has been measured via the kinematics of surrounding gas (e.g., Scharwächter et al. 2013;den Brok et al. 2015;Gravity Collaboration et al. 2018;Kabasares et al. 2022) or stars (e.g., van der Marel 1994;Nguyen et al. 2018Nguyen et al. , 2019)), this technique requires an exceptional angular resolution that can resolve the sphere of influence of the black hole (Peebles 1972).
The mass of active galactic nuclei (AGNs) can be measured without good spatial resolution via reverberation mapping (Blandford & McKee 1982).By measuring the time delay (τ ) of the broad emission line flux against the continuum, combined with its line width (∆V ) measured from the spectrum, the mass of the black hole (M • ) can be determined using where c is the vacuum speed of light and G is the gravitational constant.The virial factor, f , is a dimensionless scale factor reflecting the geometry and kinematics of the BLR.To date, the masses of more than 100 AGNs have been measured by applying this technique to broad Hβ lines (e.g., Bentz & Katz 2015).This method can be extended to a far larger number of black holes using the size-luminosity relation of the Hβemitting zone of the BLR by estimating the time lag from the 5100 Å luminosity of the AGN accretion disk (e.g., Kaspi et al. 2000;Bentz et al. 2013).It provides a short-cut to estimating the broad line time lag without going through a reverberation mapping campaign, enabling a way to estimate the mass of the black hole with a single spectroscopic observation; hence, it is called the single-epoch method.
Compared to Hβ, the Hα-emitting zone of BLR has been relatively unexplored, currently with the Hα lag measurements of only ∼50 AGNs (e.g., Kaspi et al. 2000;Bentz et al. 2010;Grier et al. 2017).This is because observing Hα poses more challenges than observing Hβ.For instance, Hα is in the redder part of the optical spectrum, making it vulnerable to the airglow lines as well as Fraunhofer A and B band telluric absorption lines, given appropriate redshifts.Moreover, there is no strong narrow line in the vicinity of Hα that could be used for flux calibration, whereas the calibration of Hβ emission lines can utilize the invariant and strong fluxes of [O III] narrow lines.
Nevertheless, the benefit of using broad Hα lines for the single-epoch method outweighs its difficulty.First, measuring Hα flux is more reliable than measuring Hβ.The Hα line is stronger than the Hβ line by at least a factor of 3, and this factor increases to 4-6 for broad emission lines (Netzer 1990).Some AGNs even exhibit a relatively weak, if present at all, broad Hβ emission (Osterbrock 1981).
Furthermore, broad line fluxes of Hα, as well as Hβ, can be measured with less degeneracy than the continuum luminosity, making it an ideal proxy for the sizeluminosity relation.The observed continuum luminosity in the AGN spectrum is contaminated by the starlight from its host galaxy or synchrotron radiation from the jet in the case of radio-loud AGNs, which must be removed to use the relation.The removal of host stellar emission can be achieved either by modeling the image of the AGN to determine the host galaxy fraction to the AGN spectrum with high-resolution images (e.g., Bentz et al. 2013) or by decomposing the continuum spectrum as a sum of the stellar and AGN components in the high S/N spectra (e.g., Park et al. 2012).The removal of jet contamination would require multi-wavelength observations (e.g., Paltani et al. 1998;Soldi et al. 2008).The broad Hα/Hβ lines, on the other hand, are purely from the BLR of the AGN and can be separated from narrow lines.
There is, however, one difficulty in using Hα for singleepoch mass estimation: a size-luminosity relation involving Hα luminosity has not yet been reported.As a workaround, Greene & Ho (2005) demonstrated the empirical relation between the broad Hα line luminosity and the 5100 Å continuum luminosity and proposed to use it in conjunction with the Hβ size-luminosity relation to construct a Hα-based single-epoch mass estimator.To date, it has been applied to a number of AGNs that are too faint to measure the AGN 5100 Å luminosity and/or Hβ line width correctly.In particular, the masses of low-luminosity AGNs and active intermediatemass black holes (IMBHs) were measured using this recipe (e.g., Reines et al. 2013;Shin et al. 2022), which suffers from substantial uncertainty due to the scatter of the scaling relations.Therefore, a relation between the size of the Hα BLR and the broad Hα luminosity will provide more robust estimations.
The Seoul National University AGN Monitoring Project (SAMP; Woo et al. 2019b;Rakshit et al. 2019) is a reverberation mapping campaign aimed at the Hβ time lags of dozens of high-luminosity AGNs to expand the size-luminosity relation toward a higher luminosity regime.In this paper, we present the SAMP results on Hα time lag measurements and demonstrate the new empirical relation between the Hα BLR size and the broad Hα luminosity.In section 2, we describe the data acquisition and reduction.In section 3, we perform spectral decomposition and Hα flux measurements.The time lag measurements are provided in section 4. Section 5 presents the size-luminosity relation of the Hα broad line.We discuss the implications of this sizeluminosity relation in section 6. Section 7 gives a brief summary of the paper.Throughout this paper, we adopt a flat ΛCDM cosmology with H 0 = 72 km s −1 Mpc −1 and Ω m = 0.3.

OBSERVATIONS AND DATA REDUCTION
The initial sample of SAMP observations consisted of 100 AGNs selected from the literature, which was described in detail by Woo et al. (2019b).To briefly summarize, 85 AGNs in local universe (z < 0.5) with V < 17 were selected from the Million Quasars catalog (MILLI-QUAS, Flesch 2015Flesch , 2021)), whose observed-frame lags were expected to be 40 < (1 + z)τ Hβ < 250 days based on the R-L relation of Bentz et al. (2013).The other 15 AGNs were selected from the Palomar-Green catalog (Boroson & Green 1992).During the first few years, we were able to identify AGNs with very low variability.Note that since the expected lag of the sample is relatively large, we were able to predict whether the line flux would vary at each epoch based on the photometric light curves.By selecting the most variable sources, we narrowed down the sample to 32 objects for continuous monitoring for six years, by excluding objects with weak variability.In this paper, we specifically focus on 13 objects of which Hα lines were observable with our spectral configurations.

Photometry
We carried out our photometric monitoring observations using several telescopes, including MDM 1.3 m and 2.4 m telescopes, the Lemmonsan Optical Astronomy Observatory (LOAO) 1 m telescope, the Lick observatory 1 m nickel telescope, the Las Cumbres Observatories Global Telescope (LCOGT) network, and the Deokheung Optical Astronomy Observatory (DOAO) 1 m telescope.The acquisition of the photometric images, reduction processes, and photometry are described in our previous paper (Woo et al. in prep.).Here, we used fully reduced and intercalibrated B and V band light curves.Typically, the B-band light curves spanned ∼2000 days with a median cadence of 4 days, resulting in ∼250 epochs, except for Mrk 1501, which was observed for 115 epochs over 1740 days with a median cadence of 6 days.The V-band light curves were acquired with a median cadence of one week.

Spectroscopy
Spectroscopic observations of Hα lines were carried out using the Shane 3 m telescope, located at the Lick observatory on Mt.Hamilton, California, USA.Note that while we used the Lick 3 m and MDM 2.4 m telescopes for the SAMP, we only use the Lick 3 m data, which covers the Hα line.The details of the spectroscopic observations of the SAMP were described by Rakshit et al. (2019).
We used the Kast double spectrograph, 1 which employs dichroic beamsplitters to acquire the red side and blue side spectra simultaneously.We used the red-side spectra with a 600 lines/mm grating.At the beginning of the campaign, the wavelength coverage of our spectra was 4450-7280 Å with 2.33 Å /pixel sampling, and the spectra obtained during this period are hereafter denoted as early configuration spectra.In September 2016, the detector was replaced with a 2K×4K CCD, covering 4750-8120 Å with 1.27 Å /pixel sampling, and the wavelength coverage was slightly adjusted to 5050-8424 Å in March 2019.Spectra obtained after September 2016 are hereafter denoted as late configuration spectra, which constitutes 80% of epochs.We used a 4 slit width to minimize the slit loss.The instrumental resolving power is R = 650, which was measured from unblended airglow lines near 7500 Å .This corresponds to the FWHM velocity of 460 km s −1 .Note that the actual resolution of AGN spectra would be better than what is measured from the night sky emissions since the slit width is wider than the seeing FWHM (1. 5-4 ).
Each night, we obtained the bias, arc, and flat frames at the beginning and end of the night.Note that the arc lamp images were taken using the 0. 5 width slit to improve the accuracy of the wavelength solution.We also observed at least one of the spectrophotometric standard stars listed by Oke (1990), and any spectra taken on nights without spectrophotometric stars were discarded from the Hα analysis.
The red-side spectra of Lick/Kast were preprocessed primarily using PypeIt v1.4 (Prochaska et al. 2020a,b).This pipeline was chosen to minimize the human intervention in the fitting of the wavelength solution and the sensitivity function.The latter is particularly susceptible to human factors due to the highly variable telluric OH absorption band near the red-side edge of the spectra.We created pixel flats and traced the slit using dome flat-field images.The wavelength solutions were derived using Ne and Ar lines in the arc frames in full template mode, and barycentric corrections were applied to each object frame.We used the optimal extraction algorithm (Horne 1986), implemented in PypeIt, to obtain photon-count spectra because optimally extracted spectra yield higher S/N than those produced with standard aperture extraction.The optimal extraction algorithm is generally not recommended for extended objects such as AGNs with resolved NLRs (e.g., Barth et al. 2015).However, we confirm that all objects in our campaign did not show extended narrow lines, even under the best seeing conditions (typically ≤1.5).Furthermore, the resulting optimally extracted spectra showed no differences compared to the aperture-extracted spectra, except for having a higher S/N.

Flux calibration
We visually inspected all spectra of spectrophotometric standard stars.After masking strong Balmer absorption lines, sensitivity functions were constructed from each reliable standard star spectrum by jointly fitting it with polynomial functions of wavelength and the model of telluric absorption lines using a script provided by PypeIt in IR mode.This yielded 1-6 different sensitivity functions per successful night.We derived the median value of the individual sensitivity function, and the function that showed a median sensitivity of the given night was chosen to be the representative sensitivity function of that night.
Spectra of AGNs and standard stars were then calibrated using the representative sensitivity for that night.We discarded any spectra taken on the nights that failed to produce at least one reliable sensitivity.Then, the atmospheric extinction was corrected based on the airmass difference between the sensitivity function and the object frame.
To further calibrate the flux in each spectrum, we compared the synthetic V-band flux obtained from the spectrum with the photometric light curves.First, we constructed the Javelin model (Zu et al. 2016) of the photometric V-band light curve and interpolated it onto each spectral epoch.The B-band light curve, which has a much shorter cadence, was jointly modeled with the Vband to improve the quality of the interpolation.Then, each spectrum, after being multiplied by the V-band filter transmission curve, was integrated to synthesize the V-band flux.Finally, we scaled the spectra so that the synthetic V-band flux was the same as the photometric flux.Note that the late configuration spectra did not cover the entirety of the V-band bandwidth.For these epochs, we calculated the portion of V-band flux that was included in each spectrum based on early configuration spectra, which was ∼80% on average, assuming that the spectral shape of AGNs did not change significantly over the campaign period.

Telluric correction
The telluric absorption features were corrected using a set of atmospheric absorption line models of the Lick observatory provided by PypeIt.This consists of 28,413 high-spectral resolution model spectra of telluric absorption lines.We first performed principal component analysis (PCA) on the models and chose the largest 25 components.Then, we fitted the calibrated spectra of spectrophotometric standards with PCA components using pPXF (Cappellari 2017) to obtain the weight for each component to reconstruct a high-resolution telluric model for each standard spectrum.We averaged the telluric models for each night to create the nightly telluric model spectrum.Then, for each object spectrum, we masked the narrow lines and fitted its red part (λ obs > 6500 Å ) with the nightly model and polynomial functions using pPXF.After shifting and broadening the model by the best-fit parameters provided by pPXF, we divided the spectrum by the model to generate a telluric corrected spectrum of our targets.An example of telluric correction is shown in Figure 1.After correcting for the telluric absorption lines, the corrected flux level is consistent with the nearby continuum.However, we note strong residuals for a couple of spectra for each object after removing the telluric lines, which is presumably due to the decreased signal-to-noise ratio.Nevertheless, these residual features did not affect our spectral analysis, as the residual features are far narrower than the AGN emission lines.

Shift correction
Despite the wavelength calibration and the barycentric corrections, the wavelength solutions deviate between the different spectra by several angstroms because of instrumental flexure and/or pointing accuracy within the slit.To compensate for this, we shifted each spectrum so that the peak of the Hα line falls exactly on the theoretical wavelength of Hα.Note that the Hα line is easier to use for this purpose compared to much weaker narrow emission lines.We first calculated the derivatives of the individual spectra by applying a Savitzky-Golay filter using SciPy (Virtanen et al. 2020) with a window width of ∼ 1260 km s −1 , and the peak of Hα was calculated by finding the root of the derivatives.The window width was chosen based on experiments so that the secondary peaks due to [N II] lines were smoothed out.Finally, we shifted and resampled the spectra.
Among the full sample of 32 AGNs that were monitored, 13 objects showed Hα in the observed spectra, depending on the redshift.Four of them suffered from strong telluric absorption because the Fraunhofer Aband fell on the very center of the Hα line.For these objects, the described telluric correction was unreliable, so we did not analyze it further.Additionally, there were fewer than 20 spectra available for 3 objects, which is unsuitable for time-series analysis.We discarded these objects as well.We present the analysis of the remain- ing 6 objects, whose properties are summarized in Table 1.

SPECTRAL ANALYSIS
Many of the early reverberation mapping studies measured the broad Hα line flux by directly integrating the spectrum within a fixed range without fitting the line profile with a model.In this case, the continuum below the emission line was fitted with a straight line at the two ends of the Hα line profile (e.g., Kaspi et al. 2000;Bentz et al. 2010).While this procedure is straightforward, we found it insufficient for our objects at higher redshift.First, the Hα emission line of our objects is located close to the edge of the detector, making it challenging to directly determine the continuum level.This is further complicated by the presence of telluric OH absorption lines near the edge of the detector.Although the absorption was averaged out through correction, some epochs exhibited strong residuals due to velocity mismatch.Finally, some of our AGNs displayed moderate Fe II lines in the blue-side continuum of Hα.While their fluxes are relatively small compared to the continuum or Hα emission, they are still strong enough to influence the slope of the continuum fit, thereby reducing the accuracy of the Hα flux.This issue is similar to what was pointed out by Barth et al. (2015) regarding the construction of Hβ light curves.It is therefore preferable to model the line profiles and measure the broad Hα flux at each epoch.
To do this, we first constructed the mean spectrum for each object.We modeled its continuum as a power law along with the Fe II lines based on the model by Boroson & Green (1992) using suitable windows, i.e., 4175-4250 Å , 4500-4725 Å , 5090-5780 Å , 6000-6280 Å , and 6800-7650 Å , in the rest frame, if covered by the spectrograph.The portions of the spectrum that showed either (1) strong telluric residual or (2) strong narrow lines were masked before the continuum fitting, leaving at least one window on the blue side of Hβ, one between Hα and Hβ, and one in the red side of Hα.We did not include the stellar host continuum in the model since our 6 objects do not exhibit strong stellar absorption features.We determined the best-fit model based on the maximum-likelihood method using the zeus MCMC sampler (Karamanis & Beutler 2020;Karamanis et al. 2021).The best-fit models of the continuum and Fe II lines were subtracted from the mean spectrum, leaving the line spectrum only.
To constrain the narrow line profile, we first modeled the [S II] λλ6717, 6731 doublet.First, we fitted the wing of the broad Hα as a cubic polynomial in the windows of 6685-6708 Å and 6760-6785 Å .After subtracting the model of the broad Hα wing, each of the [S II] lines was fitted with a single Gaussian profile using the maximum likelihood estimators.The acquired models of the [S II] were then subtracted from the observed spectrum.Then, the narrow Hα and [N II] lines were modeled as 1-2 Gaussian profiles, with shared shifts and widths among different lines, along with a sum of 2-4 independent Gaussian profiles for broad Hα.Upon fitting, we masked [O I]λλ6300, 6364 lines and, in the case of VIII Zw 218, the telluric line residuals as well.We imposed a prior such that the narrow line shifts and widths follow normal distributions centered at the measurement from [S II] with a standard deviation of ∼ 100 km s −1 .Furthermore, we restricted the parameters to have the following bounds: 1. the velocity shift of any component is within ±1600 km s −1 , 2. the narrow line dispersion is smaller than 800 km s −1 , and 3. the line dispersion of any Gaussian in broad line model is smaller than Note-Columns are (1) object name, (2) the SDSS identifier, (3) right ascension (J2000), ( 4) declination (J2000), ( 5) redshift (from NASA/IPAC Extragalactic Database (NED)), ( 6) galactic extinction in V-band by Schlafly & Finkbeiner ( 2011), ( 7) number of epochs in photometric light curve, ( 8) median cadence of photometric light curve, ( 9) number of epochs in spectroscopic light curve, (10) median cadence of spectroscopic light curve, and ( 11) SAMPID (refer to Woo et al. 2019b).
30, 000 km s −1 , and is further restricted based on visual inspection of the spectrum.We adopted the maximum a posteriori estimator from MCMC samplings as our mean spectrum model.We also measured the FWHMs of the model line profiles and found their uncertainties based on Monte Carlo randomization for 1000 iterations.
At each iteration, we added Gaussian random noise to the model parameters based on their measurement uncertainties and constructed a randomized profile.We measured the FWHMs of the randomized profiles and adopted their standard deviation as the uncertainty of the FWHM.The decomposed Hα and other narrow lines are shown in Figure 2, and the measurements after subtracting the instrumental resolving power, R = 650, are summarized in Table 2.Note that the line widths listed here should only be used as a reference to the fit quality since the resolving power measured from the airglow lines can be underestimated, as noted in § 2.2.
To measure the broad Hα flux from each epoch, we modeled the spectra from individual nights using the mean spectrum model.We first modeled and subtracted the continuum and Fe II lines using the same wavelength windows.We masked the [O I] lines and telluric residuals as we did upon fitting the mean spectra.Then, assuming the narrow lines did not change over the period of observations, we subtracted the narrow line model obtained from the mean spectrum.After subtracting the continuum, Fe II, and narrow lines, the residual spectrum contained the broad Hα line only.We constructed the prior for the multiple Gaussian models for broad Hα as follows.For the parameters that deter-mined the shape of the Hα line (i.e., the flux ratio and 1st/2nd-moment differences between any pair of Gaussian components), we imposed Gaussian priors with the mean and the standard deviation from the mean spectrum model.We did not favor any specific value for the total flux of the Hα line and adopted a flat prior.For each spectrum, we calculated the maximum a posteriori estimators from MCMC samplings.Finally, the flux of the best-fit model was taken as the Hα flux of each epoch.Note that the uncertainty of the narrow emission line fitting does not affect the Hα lag measurements as we subtracted a constant flux of narrow emission lines in each epoch, which is expected to be non-varying during the campaign.

TIME LAG MEASUREMENTS
We measured the time lags between the continuum and Hα line light curves using the interpolated crosscorrelation function (ICCF; White & Peterson 1994) with a modified averaging scheme.We first converted the Hα line fluxes into magnitudes.After we calculated the continuum-interpolated ICCF and the lineinterpolated ICCF, we took the Fisher transformation (Fisher 1921) on both one-sided ICCFs, averaged them, and then took the inverse transformation to obtain the z-transformed average ICCF.The ICCF centroid was calculated over the largest continuous interval containing the peak of the ICCF, where the ICCF values within the interval were larger than 80% of the peak value.We performed 10,000 realizations of the flux randomization/random-subset selection (FR/RSS; Peterson et al. 1998Peterson et al. , 2004)).For each realization, we resampled each light curve and added random Gaussian  noise to each flux value according to their measurement uncertainty.Duplicate points were averaged, and their uncertainties were divided by √ n to compensate for the duplication.The median of the distribution of the centroid and the central 68% confidence interval was taken as the time lag measurement and the associated uncertainties.To check the consistency with the ICCF, we also calculated the time lag using other commonly-used methods, the z-transformed discrete correlation function (zDCF; Alexander 1997) and the Javelin model (Zu et al. 2011) Finally, we rated the quality of our lag measurements based on the lag differences among different methods as follows: • Rating A if its ICCF lag, zDCF lag, and Javelin lag agree within 1-σ and the maximum difference between them is within two months (60 days).
• Rating B if its lag is measured (1-σ above zero lag) with all 3 methods, while the maximum difference between them is larger than two months (60 days).
• Rating C if its lag is not constrained or detected.
The measured lags are summarized in   lag.However, the CCFs at time lags between 100 days and 200 days were relatively unexplored due to the large seasonal gaps.This is reflected in the large difference of ICCFs using different interpolation methods, as well as the lack of points in the zDCF.The primary peak of the Javelin lag distribution falls in this seasonal gap, which questions the accuracy of the Javelin lag for this specific case.We assess the lag of this object as rating B.
For J0101+422, we measured 117 +17 −18 days of the time lag.This is consistent with its zDCF and Javelin lags.We assess the lag of this object as rating A.
For PG 0947+396, we measured 71 +16 −34 days of the time lag.While this is consistent with the zDCF and Javelin lags, they both showed bimodality in their lag, where zDCF preferred the smaller mode and Javelin preferred the larger mode.On the other hand, the ICCF lag captured the average between the two lags.We assess the lag of this object as rating B.
For J1217+333, we assess that our lag measurements are not reliable.Javelin is unconstrained, and the zDCF lag is consistent with 0 within its 68% confidence interval.While we measured 199 +137 −132 days from the ICCF, this is likely to be the average of the window size we used to find the lag.Supporting this, the ICCF and zDCF both show multiple modes in the given window.Moreover, it is far larger than what is measured using Hβ (Woo et al. in prep.).We assess the lag of this object as rating C and exclude it from further analysis.VIII Zw 218 showed a single, clean peak in CCFs, and we measured 140 +25 −25 days of the time lag.This value is consistent with zDCF and Javelin lags.We assess the lag of this object as rating A.
We had to detrend the light curves of PG 1440+356 before measuring the time lag since the Hα light curve showed a monotonic increase in the Hα light curve.Without detrending, the CCF values at any lag were higher than 0.4, rendering the lag and its uncertainty measurements unreliable.After detrending, we measured 79 +68 −29 days of time lag.This value is consistent with zDCF and Javelin lags.We note that both ICCF and zDCF are skewed toward higher lag values, and this is reflected in the measurements.We assess the lag of this object as rating A.
We find that the uncertainties of our ICCF-based lag measurements range from 10% to 70%, a considerable portion of which stems from the number of epochs in the light curve.The lack of distinct variability features in the light curve also increases the uncertainty.However, we note that the relative uncertainties of our lag measurements are comparable to those reported in the literature.For example, 95% of the lag measurements we collected in Table 4 have uncertainty between 8%-70%.

Determining the Hα BLR size-luminosity relation
To compare the BLR size of Hα with Hα luminosity, we compiled a sample of Hα lags and luminosities from this work and the literature (Kaspi et al. 2000;Bentz et al. 2010;Barth et al. 2011;Grier et al. 2017;Sergeev et al. 2017;Cho et al. 2020;Feng et al. 2021;Li et al. 2022).While we tried to include as many AGNs as possible, the following objects were excluded.First, Grier et al. (2017) provided quality ratings for the time lags, and 4 objects (out of 18 with Hα lags measured) with ratings of 1 or 2 were excluded.Additionally, we note that Hα time lags of Ark 564 (Shapovalova et al. 2012) and NGC 7469 (Shapovalova et al. 2017) were measured.However, the Hα lag of Ark 564 was consistent with 0 delay within the error bar, while the CCF between Hα and the continuum light curves of NGC 7469 exhibited multiple peaks with r max values smaller than 0.5.Thus, we did not include these two objects in our analysis.We adopted the time lag values and their 1-σ uncertainties from each paper.For objects presented in this paper, we adopted the ICCF time lag values and their 1-σ errors.The lags presented in the observed frame were divided by 1 + z to convert them into the rest-frame lags.We obtained the broad Hα luminosities from the fluxes by multiplying 4πd L 2 , where d L is the luminosity distance.We ignored the systematic differences in methods of measuring flux between the papers, which are discussed in § 5.2.We also collected the AGN continuum luminosities at 5100 Å for the sample.We applied galactic extinction correction based on the galactic extinction map by Schlafly & Finkbeiner (2011) and the extinction curve by Cardelli et al. (1989).The corrected luminosities, as well as their time lags, are listed in Table 4.  Schlafly & Finkbeiner (2011) and the extinction curve by Cardelli et al. (1989).Time lag values are presented in the rest frame.Uncertainties shown here are the 68% confidence intervals.
We model a generic relation between two variables x and y as where σ int denotes the intrinsic scatter of the relation.We should note that NGC 4395 has an extremely low luminosity, while PG 1226+023 has an extremely high luminosity compared to other AGNs in our sample.We tested whether removing these outliers from the sample would alter the fit.As presented in Table 5, we found that excluding or including NGC 4395 and PG 1226+023 yields similar best fits, with the largest difference between parameters within the 1-σ boundary.Thus, we present the best fit for the size-luminosity relation of the Hα BLR to be where σ int = 0.31 ± 0.03.The best fit, as well as the objects used to derive the fit, are plotted in Figure 9.

Systematic differences in luminosities
We note that the narrow line fluxes were handled differently in the literature when measuring the Hα luminosity.For example, in their work, Cho et al. (2021), Feng et al. (2021), andLi et al. (2022) modeled the broad Hα line and narrow line components separately with multiple Gaussians.Similarly, Grier et al. (2017) modeled the broad line after the narrow lines were modeled and removed using high-pass filtering, which was described by Shen et al. (2016).In contrast, Kaspi et al. (2000), Bentz et al. (2010), Barth et al. (2011), andSergeev et al. (2017) integrated the flux over certain wavelength windows without removing the narrow lines of Hα, [N II], and even [S II] depending on the object.In principle, direct integration leads to inaccurate measurements by including narrow line fluxes while excluding the wing of broad Hα.These two effects can work in the opposite way such that they could cancel each other out at a certain level.We assessed that the effect of the Hα wing seems negligible for the objects we included because all aforementioned references chose a sufficiently large interval for direct integration.Moreover, the narrow line fluxes seem negligible compared to their strong Hα luminosities for high luminosity objects that Kaspi et al. (2000) or Sergeev et al. (2017) presented.How-  4) size of the subset, ( 5), ( 6), and ( 7) the best fit parameters to the model described by Eq. 2. The best-fit parameters and their uncertainties are given in the median and the standard deviation.
ever, visually inspecting the spectra presented by Bentz et al. (2010) reveals that the Hα fluxes they measured do include a substantial portion of narrow line fluxes, resulting in a larger scatter.
The measurements of the continuum luminosity suffer from a similar issue as well.Specifically, Kaspi et al. (2000) and Sergeev et al. (2017) did not remove the host galaxy contamination when measuring the continuum flux at 5100 Å.However, the host contamination would be negligible for these high-luminosity AGNs.For instance, Jalan et al. (2023) proposed an empirical correction for host contamination based on the total luminosity at 5100 Å.According to their Equation 2, on average, 30% of the continuum luminosities of 15 AGNs in the aforementioned references can be attributed to their host galaxy starlight.This corresponds to the bias of 0.15 dex, which is far smaller than the intrinsic scatter of our fit.However, we do not use their empirical correction because of its large scatter.
Despite the aforementioned issues, we do not observe any systematic deviation from the size-luminosity relations of the AGNs with inaccurate luminosity measurements.Thus, we conclude that neither the narrow lines nor the host contamination induced bias in determining the size-luminosity relation, albeit the intrinsic scatter could be overestimated.
For consistency check, we investigate whether the Hα BLR size-luminosity relation changes by excluding our new measurements.We find that the best-fit slope remains the same and the intrinsic scatter increases slightly (i.e., 0.29 ± 0.03 dex), presumably due to the small size of our new sample.However, the RMS scatter of the SAMP AGNs is smaller (0.15 dex) than that of the total sample (i.e., 0.27 dex) as the measurement errors of the SAMP AGNs are much smaller than that of the literature sample.Note that previous Hα lag and error measurements were obtained in various ways.Thus, a uniform analysis of the entire sample is required to better constrain the size-luminosity relation and its scatter.

Comparison between the Hα and continuum luminosities
The derived size-luminosity relations in § 5 have slopes very similar to each other, suggesting that the relation between the Hα luminosity and the continuum luminosity at 5100 Å is close to a linear one.Also, theoretically, we expect the Hα emission line flux to scale linearly with the optical continuum luminosity, assuming that all AGNs share the same spectral slope of the powerlaw continuum (Yee 1980).We directly compare L Hα and λL λ 5100 Å in Figure 10 by fitting the relation using the generic relation model given by Eq. 2. The best-fit parameters are presented in Table 5.The slope β = 0.97 ± 0.03 for this relation is virtually the same as unity, and the fit can be simplified as λL λ 5100 Å = (19±1)L Hα with an intrinsic scatter of 0.18 dex.The fit parameter did not change even if we excluded NGC 4395 and PG 1226+023 from the fit, suggesting that a simple linear relation is valid across 7 orders of magnitude of luminosities (10 39 < λL λ 5100 Å /[erg s −1 ] < 10 46 , or 10 38 < L Hα /[erg s −1 ] < 10 45 ).The linearity of this relation contradicts the results by Greene & Ho (2005), who found a supra-linear relation between the continuum and the line luminosities, with a slope of 1.157 ± 0.005.Note that their sample consists of the AGNs from the Sloan Digital Sky Survey (SDSS), without reverberation mapping results.The major difference compared to our study is that they did not subtract the host galaxy contamination in the continuum luminosity or remove the narrow line components for measuring the broad Hα luminosity.In the case of the narrow line contamination, Greene & Ho (2005) discussed that the narrow component contribution is only ∼7% of the Hα flux, and the slope of the relation did not change even if they used the luminosity of the broad Hα line only.On the other hand, the host contamination can systematically change the slope.While the continuum subtraction was not performed for some AGNs in Table 4 in our study, these AGNs are mostly high-luminosity objects (λL λ 5100 Å > 10 44 erg s −1 ), for which host galaxy contribution would be negligible.However, Greene & Ho (2005) did not subtract the host galaxy contribution for their sample over the entire luminosity range, including low-luminosity AGNs where the host galaxy contribution is expected to be systematically larger.While they did exclude galaxies with high host contamination based on modeling the equivalent width of the Ca II K-line, the model was constructed from the library of early-type galaxies only, which are not expected to be the host galaxies of low-luminosity AGNs.Presumably these are why we obtained a different slope compared to that of Greene & Ho (2005).A more recent study by Rakshit et al. (2020), which decomposed SDSS spectra to acquire host-free AGN continuum and broad Hα line luminosities, also demonstrated a similar supra-linear slope of 1.126 ± 0.004.However, they decomposed the host spectra using the eigenspectra constructed from late-type galaxies, which can lead to a similar template mismatch.On the other hand, studies that did not use SDSS spectra deduced slopes close to unity.For example, Shen & Liu (2012) obtained a slope of 1.010 ± 0.042 from 60 objects in z ∼1.5-2.2.Similarly, Jun et al. (2015) demonstrated a slope of 1.044 ± 0.008 with AGNs in z ∼0-6.2 and continuum luminosity between 10 42 < λL λ 5100 Å /[erg s −1 ] < 10 47 .We discuss that the slope β is unity since (1) all contradicting studies carried out using SDSS spectra had issues with template mismatch, (2) other studies support our results, and (3) the slope of unity predicts the Hα luminosity of NGC 4395 well from its continuum luminosity, while the supra-linear slope would underestimate its Hα luminosity.

Broad line region stratification
We now compare the time lag of Hα with that of Hβ.For our 5 objects, we adopted the Hβ time lags in our previous paper (Woo et al. in prep.).Most of the other objects also have Hβ lags measured as well, which are summarized in Table 4.The comparison between two lags is presented in Figure 11.The best-fit parameters for these two lags are listed in Table 5, which can be summarized as τ Hα = 1.68τHβ with an intrinsic scatter of σ int = 0.23 ± 0.03 dex, where the ratio does not depend on the lag.This shows stratified BLRs across a wide range of time lags, which has been predicted and observed before (e.g., Bentz et al. 2010 and references therein), and can be attributed to the optical depth of Hα being larger than that of Hβ.

Slope of the size-luminosity relation
The slopes of the Hα size-luminosity relation we found in § 5, 0.61 ± 0.04 and 0.59 ± 0.04 are consistent with those of the Hβ size-continuum luminosity relation by Bentz et al. (2013, 0.55±0.03with Clean+ExtCorr sample).This is expected since (1) the Hα lags are proportional to the Hβ time lags, as discussed in § 6.2, and (2) the Hα luminosity is proportional to the continuum luminosity, as discussed in § 6.1.
While the Hα size -Hα luminosity relation has not previously been reported, some studies have compared the Hβ BLR size with Hβ luminosity.Wu et al. (2004) obtained the slope of 0.684 ± 0.106 based on available Hβ lags, including a large number of objects from Kaspi et al. (2000).A similar slope, 0.687±0.063,was obtained by Kaspi et al. (2005).On the other hand, more recent studies support the Hβ size-Hβ luminosity relation having a slope close to that of Bentz et al. (2013).For instance, Greene et al. (2010) obtained a much smaller slope of 0.53 ± 0.04.Du et al. (2015) also obtained a smaller slope of 0.51 ± 0.03 with a much larger sample.
Our result is broadly consistent with these studies, although further studies are needed to determine the slope of the size-luminosity relation using the broad Hβ luminosity.

Single-epoch mass comparison and its implications for IMBH studies
The mass of the black hole can be estimated by Eq. 1.We propose a new single-epoch mass estimator using the luminosity and velocity of broad Hα, where we adopted f = 4.47 for σ and f = 1.12 for FWHM from Woo et al. (2015).
In comparison, Greene & Ho (2005) proposed a method to estimate the size of Hα as follows.The continuum luminosity is first estimated from Hα broad line luminosities using Eq. 1 of their paper.Then, the Hβ size -continuum luminosity relation is used to estimate the Hβ BLR size.A commonly adopted relation is Eq. 2 from Bentz et al. (2013)  which will hereafter be denoted as GH+B.Using this along with Eq. 3 by Greene & Ho (2005), one can also construct mass estimates similar to Eqs. 5 and 6.An example of the analog relation can be found in Reines et al. (2013).
However, the new mass estimator predicts a substantially different mass from the one the GH+B estimator predicts.In Figure 9, the GH+B relation is plotted as a dashed line, along with the size-luminosity relation we obtained as a solid line.While the difference between two relations is larger than 0.5 dex if L Hα ≤ 10 38 erg s −1 , the Hα size of NGC 4395 is even smaller than what our relation predicts.This difference can be understood as the direct result of the supra-linear relation between L Hα and λL λ 5100 Å by Greene & Ho (2005), as discussed in § 6.1.
To quantify the difference, we compared our mass estimate with that based on the GH+B relation in Figure 12.While it is negligible for AGNs with masses above 10 9 M , the black hole mass can be overestimated by more than a factor of two with the GH+B relation for low-luminosity AGNs with L Hα ≤ 10 40 erg s −1 , which are strong candidates for IMBHs.
IMBHs are important objects because of their relevance to the formation of SMBHs since different scenarios of SMBH formation predict different IMBH mass functions and mass scaling relations even in the local universe (M • -M * , M • -σ * , etc.; e.g., Miller et al. 2015;Greene et al. 2020).In particular, according to the direct collapse scenario, IMBHs at the centers of galaxies are predicted to have minimum masses of about 10 4 -10 5 M (e.g., Ferrara et al. 2014).
If our new mass estimator is accurate, the mass measurements of active IMBHs in the literature may have been biased toward larger values.For example, using the GH+B relation, Reines et al. (2013) estimated the masses of AGNs in dwarf galaxies to be in the range of 10 5 M -10 6 M .However, considering that the Hα luminosities of these AGNs range from 10 38 erg s −1 to 10 40.5 erg s −1 , it is possible that the reported masses are overestimated by an average factor of 2-3, with some potentially having a mass closer to 30, 000 M .The slope of the Hα-based mass estimator is of critical importance.However, our combined sample does not include any object with 10 39 erg s −1 < L Hα < 10 41 erg s −1 , except for one object (NGC 4395) with a luminosity smaller than this range.It is necessary to include lowluminosity AGNs, which are more representative of active IMBHs, to properly constrain the slope of the sizeluminosity relation and the masses of IMBHs.
Conducting a reverberation mapping campaign for AGNs with Hα luminosities below 10 41 erg s −1 poses substantial challenges.First, these low luminosity AGNs are expected to exhibit very short Hα lags, requiring intra-night monitoring campaigns with a cadence of several hours or even minutes and continuous observations in a several-day time baseline.Such campaigns require the coordination of multiple telescopes.Second, observing these low luminosity AGNs requires higher sensitivity, which implies the use of larger aperture tele-scopes and/or longer exposure times.Despite these difficulties, monitoring campaigns for these low luminosity AGNs will be crucial to better constrain the Hα S-L relation and to calibrate the Hα-based mass estimators, particularly for IMBHs.

CONCLUSIONS
We present the Hα reverberation mapping results from the Seoul National University AGN Monitoring Project.While the SAMP mainly aims at performing Hβ reverberation mapping for more than 30 highluminosity AGNs (Woo et al. in preparation), we additionally obtained time series of Hα spectra for 6 objects and performed the reverberation mapping analysis.By combining our new measurements with the Hα lag measurements of other AGNs in the literature, we investigated the size-luminosity relation of the broad Hα line.Our main results are summarized as follows.
• We produced Hα light curves based on the spectral modeling of Hα emission line and measured the time lags against B band continua for 6 new objects, 5 of which we consider to be reliable.
• We collected a sample of AGNs with the lag and flux measurements of broad Hα of 47 AGNs, consisting of our 5 new objects and 42 from the literature.We calculated the Hα luminosities after correcting for Galactic extinctions.
• The size-luminosity relation we obtained based on the reverberation mapping results deviates from what is proposed based on single-epoch spectra by Greene & Ho (2005).We demonstrate that for AGNs in the IMBH mass regime, the black hole mass could be substantially overestimated by a factor of 3 on average if the relation by Greene & Ho (2005) is used.
• We propose two mass estimators based on Hα broad lines assuming f = 4.47 for σ and f = 1.12 for FWHM (Woo et al. 2015) SRG/2021/001334.We thank the staff of the observatories where data were collected for their assistance.

Figure 1 .
Figure 1.An example of telluric correction.The black line represents the flux-calibrated spectrum of PG 0947+396 without the telluric correction.Red represents the telluric model in the optical depth unit (arbitrarily scaled), with the velocity and velocity dispersion adjusted to fit the AGN spectrum.Blue is the AGN spectrum after correcting for telluric absorption.

Figure 2 .
Figure 2. Mean spectra of 6 objects around the Hα line and their best-fit models.Black represents the mean spectrum, green represents the power-law continuum, and brown represents the Fe II model from Boroson & Green (1992).Each blue line represents a narrow emission line.The thick red line shows the broad Hα model, whereas thin red lines show the individual Gaussian components of the model.

Figure 3 .Figure 4 .Figure 5 .Figure 6 .Figure 7 .Figure 8 .
Figure 3.Light curves and the time lag measurements of Mrk 1501.Upper left: B-band light curve, with the Javelin model and its uncertainty shown as a solid line and shaded region.Lower left: Broad Hα light curve, with the Javelin model shown similarly.Upper right: Blue lines indicate the ICCF, where the dashed line shows the ICCF upon interpolating the continuum, and the dash-dotted line shows the ICCF upon interpolating the line flux.The solid line represents the z-transformed average between two interpolations, as described in §4.Red dots indicate the zDCF, where the horizontal error bar indicates the bin size of the zDCF, and the vertical error bars indicate the standard deviation within the bin.Lower right: The time lag measurements.ICCF, zDCF, and Javelin are indicated by blue, red, and green colors, respectively.Vertical solid lines represent the median of each measurement, while dashed lines mark the 16th and 84th percentiles as 1-σ uncertainties.

log 10 LFigure 9 .
Figure 9. Left: The Hα time lag against the Hα line luminosity.The best-fit relation is denoted with a thick solid line, with the shaded region representing the 1-σ intrinsic scatter.The dashed line represents the estimated Hβ lag from the Hα luminosity by combining Equation 1 from Greene & Ho (2005) and Equation 2 from Bentz et al. (2013) with Clean+ExtCorr fit parameters (See §6.4).Right: The Hα time lag against the continuum luminosity at 5100 Å.The best-fit relation is denoted with a thick solid line, with the shaded area representing the 1-σ intrinsic scatter.The dashed line represents the estimated Hβ lag from the continuum luminosity using Equation 2 from Bentz et al. (2013) with Clean+ExtCorr fit parameters.

Figure 10 .
Figure10.The comparison between the continuum luminosity, λL λ 5100 Å , and the broad Hα luminosity, LHα.The best-fit relation is denoted with a thick solid line, with the shaded area representing the 1-σ intrinsic scatter.The dashed line along the diagonal represents Equation 1 fromGreene & Ho (2005).

Figure 11 .
Figure 11.The comparison between the Hα and Hβ lags.The best-fit relation is denoted with a thick solid line, with the shaded area representing the 1-σ intrinsic scatter.The dashed along the diagonal represents the τHα = τ Hβ line.

Figure 12 .
Figure 12.Demonstration of M • overestimation depending on different parameters.The abscissa represents the singleepoch mass of the black hole estimated from Hα luminosities using the relation we obtained, while the ordinate represents the excess of the mass when estimated using Greene & Ho (2005) and Bentz et al. (2013) compared to our estimate.The red lines show the contours with the same Hα line widths, while the blue lines show the contours with the same Hα luminosity.
We thank the anonymous referee for constructive comments that improved the manuscript.This work has been supported by the Basic Science Research Program through the National Research Foundation of Korean Government (2021R1A2C3008486) and the Samsung Science & Technology Foundation under Project Number SSTF-BA1501-05.The work of H.C. was supported by an NRF grant funded by the Korean Government (NRF-2018H1A2A1061365-Global Ph.D. Fellowship Program).S.W. acknowledges the support from the National Research Foundation of Korea (NRF) grant funded by the Korean government (MEST) (No. 2019R1A6A1A10073437). V.N.B. gratefully acknowledges assistance from the National Science Foundation (NSF) Research at Undergraduate Institutions (RUI) grant AST-1909297.Note that the findings and conclusions do not necessarily represent the views of the NSF.Research at UCLA was supported by the National Science Foundation through grant NSF-AST 1907208.Research at UC Irvine was supported by NSF grant AST-1907290.V.U.acknowledges funding support from NASA Astrophysics Data Analysis Program (ADAP) grant 80NSSC20K0450.Support for program HST-AR-17063.005was provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Associations of Universities for Research in Astronomy, Incorporated, under NASA contract NAS526555.S.R. acknowledges the partial support of SRG-SERB, DST, New Delhi through grant no.

Table 2 .
Line Width Measurements

Table 4 .
The Time Lags and Luminosities

Table 5 .
The Best-fit Parameters