Time Delay of Mg ii Emission Response for the Luminous Quasar HE 0435-4312: toward Application of the High-accretor Radius–Luminosity Relation in Cosmology

The photometric data were combined from a few dedicated monitoring projects, described in more detail in Zajaček et al. (2020). The source has been monitored in the V band as part of the OGLE-IV survey (Udalski et al. 2015) done with the 1.3 m Warsaw telescope at the Las Campanas Observatory, Chile. The exposure times were typically around 240 s, and the photometric errors were small, of the order of 0.005 magnitude (see Table 3). In the later epochs the quasar was observed, again in the V band, with the 40 cm Bochum Monitoring Telescope (BMT). ¹¹ These data showed a systematic offset of 0.2 magnitudes toward larger magnitudes with respect to the overlapping OGLE data, which we corrected for by the shift of all of the BMT points. SALT spectroscopic observations were also supplemented, whenever possible with SALTICAM imaging in the g band. We have analyzed all of these data; however, two data sets (2013 August 20 and 2019 January 27) showed a significant discrepancy with the other measurements. The first of the two sets of observations was done during full moon, and with the moon–target separation relatively large; spectroscopic observations were not affected, but the photometric observations were. During the second set of observations the night was dark but seeing was about 25 during the photometry, dropping down to 2'' during the spectroscopic exposures. We were unable to correct the data for these effects, and we did not include these data points in further considerations. Because of the collection of the data in the g band, we allowed for the grayshift of all the SALTICAM data, and the shift was determined using epochs when they coincided with the more precise OGLE set collected in the V band. Finally, we supplemented our photometry with the light curve from the Catalina Sky Survey, ¹² which is not of a very high quality (with uncertainties of 0.02–0.03 mag) but nicely covers the early period from 2005 until 2013. We have binned these data to reduce the scatter. Table 3 contains only the data points that were used in time-delay measurements. All the data points are presented in the upper panel of Figure 1.

Zoom In Zoom Out Reset image size

Recently, the quasar was monitored in the V band with a median sampling of 14 days using Lesedi, a 1 m telescope at the South African Astronomical Observatory (SAAO), with the Sutherland High Speed Optical Camera (SHOC). SHOC has a 5.7 × 5.7 arcmin² field of view (FoV). Each observation consists of nine dithered 60 s exposures. They are corrected for bias and flatfields (using dusk or dawn sky flats). Astrometry is obtained using the SCAMP tool. ¹³ The resulting median-combined image has a 7.5 × 7.5 arcmin² FoV centered on the quasar. The light curves were created using five calibration stars located on the same image as the quasar. The preliminary results are consistent with the last photometric point from SALT/SALTICAM.

Spectroscopic observations of HE 0435-4312 were performed with the 11 m telescope SALT, with the RSS (Burgh et al. 2003; Kobulnicky et al. 2003; Smith et al. 2006) in a long-slit mode and a slit width of 2''. We used the RSS PG1300 grating with a grating tilt angle of 2675. Order blocking was done with the blue PC04600 filter. Two exposures were always made, each of about 820 s. The same setup has been used throughout the whole monitoring period, from 2012 December 23 until 2020 August 20. Observations were always performed in the service mode.

The raw data reduction was performed by the SALT observatory staff, and the final reduction was performed by us with the use of the IRAF package. All the details were given in Sredzińska et al. (2017), where the results from the first three years of this campaign were presented. We followed exactly the same procedure for all 25 observations to minimize the possibility of unwanted systematic differences.

In order to get the flux calibration, we performed a weighted spline interpolation of the first order (with inverse measurement errors as weights) between the epochs of photometric and spectroscopic observations, thus an apparent V magnitude was assigned to each spectrum. Taking as a reference the composite spectrum and the continuum with a slope of α_λ = −1.56 from Vanden Berk et al. (2001), we just normalized each spectrum to the V magnitudes (5500 Å).

The reduced and calibrated spectra were fitted in the 2700–2900 Å range in the rest frame. The basic model components were as in Sredzińska et al. (2017). The data were represented by (i) continuum in the form of a power law with arbitrary slope and normalization, (ii) the Fe II pseudocontinuum, and (iii) two kinematic components representing the Mg ii line; each of the kinematic components was represented by two doublet components. The line is unresolved, and the doublet ratio could not be well constrained, so it was fixed at 1:1 (see Sredzińska et al. 2017 for the discussion). For the kinetic shapes we used Lorentzian profiles since they provided somewhat better representation of the data in χ² terms than Gaussian ones. All the parameters were fitted simultaneously. In order to determine the redshift and the most appropriate Fe II template, we studied in detail observation 23, which also covered the region around 3000 Å in the rest frame (see Appendix B). Thus for the adopted redshift z = 1.2231, the Fe II template was very slightly modified in comparison with Sredzińska et al. (2017), and the pseudocontinuum was kept at FWHM of 3100 km s⁻¹. Thus there were eight free parameters in the model: power-law normalization and slope, normalization of the Fe II pseudocontinuum, the width and normalizations of the two kinematic components of Mg ii, and the position of the second component, with the other one set at the quasar rest frame, together with Fe II.

We determined the mean and the rms spectra for HE 0435-4312 as we did previously for quasar HE 0413-4031 (Zajaček et al. 2020). Due to the particularly low signal-to-noise ratio (∼7.5) shown by spectrum no. 19, which is clearly an outlier in the light curve (Figure 1), we do not consider it for the estimation of the mean and rms spectra. We followed the methodology for constructing the mean and the rms spectra as outlined in Peterson et al. (2004). In particular, we formed the mean spectrum (without spectrum no. 19) using

$$\begin{eqnarray}&&\overline{F(\lambda )}=\displaystyle \frac{1}{N}\displaystyle \sum _{i=1}^{N}{F}_{i}(\lambda ),\end{eqnarray} \tag{ 1 }$$

where F_i (λ) is the ith spectrum out of a total of N spectra in the measured database. Subsequently, the rms spectrum, initially taking into account all spectral components (Mg ii line + continuum + Fe II pseudocontinuum), was constructed using

$$\begin{eqnarray}&&S(\lambda )={\left\{\displaystyle \frac{1}{N-1}\displaystyle \sum _{i=1}^{N}{\left[{F}_{i}(\lambda )-\overline{F(\lambda )}\right]}^{2}\right\}}^{1/2}.\end{eqnarray} \tag{ 2 }$$

The constructed mean and rms spectra are shown in Figure 2 (black solid lines) in the upper and upper middle panels, respectively. The quasar is not strongly variable, so the normalization of the rms is very low, and the spectrum is noisy, with visible effects of occasional imperfect sky subtraction. However, the overall quality of the rms spectrum is still suitable for the analysis. In both the mean and the rms spectra, the Mg ii line modeling required two kinematic components, since line asymmetry is clearly visible. The result is shown in Figure 2 with spectral components depicted with different colored lines described in the figure caption. For the fitting, we finally used the redshift as determined by Sredzińska et al. (2017), but with a slightly modified Fe II template based on the d11 template of Bruhweiler & Verner (2008). We also compared and analyzed other Fe II templates based on the updated CLOUDY model (Ferland et al. 2017) as well as the six-transition model by Kovačević-Dojčinović & Popović (2015) and Popović et al. (2019). For a detailed discussion of different Fe II templates—a total of eight setups with different redshifts as well as Lorentzian or Gaussian line component profiles—see Appendix B.

Zoom In Zoom Out Reset image size

The overall shapes of the mean and the rms spectra are similar (see Figure 2). The FWHM of the Mg ii line in the mean spectrum is ${3695}_{-21}^{+21}\,\mathrm{km}\,{{\rm{s}}}^{-1}$; the line in the rms spectrum might be slightly broader, ${3886}_{-341}^{+143}\,\mathrm{km}\,{{\rm{s}}}^{-1}$, but is consistent within the error margins. A much larger difference is seen in the line dispersion, which is much smaller in the mean spectrum than in the rms spectrum (${2707}_{-6}^{+10}$ km s⁻¹ versus ${3623}_{-412}^{+76}$ km s⁻¹). The FWHM and σ are larger in the rms spectrum, most likely due to its noisy nature, although there is an indication of a trend of both Mg ii FWHM and σ being larger in the rms spectrum than in the mean spectrum based on the analysis of the SDSS-RM sample (Wang et al. 2019). However, the ratio of FWHM to σ for both the mean and the rms spectra is far from the value 2.35 expected for a Gaussian shape. The source can be classified as A-type according to the classification of Sulentic et al. (2000). This is consistent with the Eddington ratio 0.58 determined by Sredzińska et al. (2017) for this object. There is also an interesting change in the line position between the mean and the rms spectra, as determined from the first moment of the distribution: 2805 Å versus 2792 Å.

Following the spectral studies by Barth et al. (2015) and Wang et al. (2019), we show in Figure 2 (lower middle and lower panels) the mean and the rms spectra constructed taking into account all the spectral components (black lines; Mg ii line + continuum + Fe II pseudocontinuum), the mean and the rms spectra that have the continuum power law subtracted from each spectrum (red lines), and finally the spectra with Fe II pseudocontinuum subtracted from individual spectra, which represents the true, line-only mean and the rms spectrum (blue lines). For the rms spectra, we do not detect any significant difference in the line width, which is expected to be smaller for the total-flux rms than for the line-only rms spectrum (Barth et al. 2015; Wang et al. 2019). This can be attributed to the overall noisy nature of our rms spectra, which are constructed from only 24 individual spectra. According to Barth et al. (2015) and Wang et al. (2019), the difference in the line width becomes smaller for a sufficiently long duration of the campaign. However, in our case, the spectroscopic monitoring was only ∼4.25 times longer than the emission-line time lag in the observer's frame (see Section 4). For such a short duration, the distribution of the ratios of line width between the total-flux and the line-only rms spectra is broad—between 0.5 and 1.0, with the peak close to 1.0; see Appendix C in Barth et al. (2015).

The fits to individual spectroscopic data sets were done in the same way as for the mean spectrum. In all 25 data sets, two kinematic components of the Mg ii line were needed to represent the line shape well. The normalization of Fe II pseudocontinuum was allowed to vary for each individual spectrum, while the Fe II width was fixed to FWHM = 3100 km s⁻¹. The Mg ii component kinematically related to the Fe II is slightly narrower, having an average FWHM of 2128 ± 28 km s⁻¹, while the second shifted component is somewhat broader, with FWHM of 2262 ± 90 km s⁻¹. However, we stress here that the broadband modeling by Sredzińska et al. (2017) did not support any identification of these components with separate regions, so the two components are just a mathematical representation of the slightly asymmetric line shape.

The parameters for observation 19 were considerably different, but as we already mentioned in Section 3.1, this observation was of a particularly low quality despite the fact that actually three exposures were made this time. Cirrus clouds were present during the whole night, and even more clouds started arriving during the second exposure of the quasar, so the third exposure was done. Nevertheless, all the parameters were determined with errors a few times larger than for the remaining observations.

The average value of the equivalent width of the Mg ii line is 23.6 ± 0.5 Å if observation 19 is included. If observation 19 is not taken into account, the mean value increases to 24.1 ± 0.6 Å, and this is similar to the value determined from the mean spectrum: 23.9 Å. Such values are perfectly in agreement with the properties of the bright quasar sample studied by Forster et al. (2001), where the mean EW of the Mg ii broad component was found to be ${27.4}_{-6.3}^{+8.5}$ Å. We do not see any traces of the narrow component in either the mean spectrum or the individual data sets.

The average value of the EW of Fe II, 18.2 ± 0.7 Å, is also similar to the value in the mean spectrum, 18.9 Å. The relative error is larger than for the Mg ii line since Fe II pseudocontinuum is more strongly coupled to the power-law continuum during the fitting procedure.

The dispersion in the measurements between observations partially represents the statistical errors, but partially reflect the intrinsic evolution of the source in time. This evolution is studied in the next section.

The quasar HE 0435-4312 is not strongly variable in terms of the continuum emission. The whole photometric light curve, including the CATALINA data, covers 14 yr, and the fractional variability amplitude for the continuum is 8.9% in flux. During the period covered by SALT data (7 yr) it is ∼4.9%. Fortunately, the Mg ii line flux shows significant variability, comparable to the continuum, at the level of 5.4% during this period. We do not observe a suppressed Mg ii variability in this source, unlike that seen in much larger samples (Goad et al. 1999; Woo 2008; Zhu et al. 2017). Yang et al. (2020) also detect the response of Mg ii emission to the continuum for 16 extreme-variability quasars, but with a smaller variability amplitude, ${\rm{\Delta }}\mathrm{log}L(\mathrm{Mg}\,{\rm\small{II}})=(0.39\pm 0.07){\rm{\Delta }}\mathrm{log}L(3000\,\mathring{\rm A} )$. A low Mg ii variability was modeled to be the result of a relatively large Eddington ratio (Guo et al. 2020) but HE 0435-4312 is also a source with a rather large Eddington ratio of 0.58 (Sredzińska et al. 2017). Overall, the fractional variability of the Mg ii line for our source is comparable to a variability of ∼10% on 100 day timescales as inferred for the SDSS-RM ensemble study (Sun et al. 2015).

The variations in quasar parameters are smooth overall. The quality of observation 6 was not very high, as discussed in Sredzińska et al. (2017), but it still allowed the Mg ii line parameters to be obtained properly. However, observation 19 created considerable problems. The registered number of photons was much lower (by a factor of a few) than in typical observations, even the background was rather low, and clouds were apparently present in the sky. We did the data fitting, and the derived values of the model parameters form clear outliers when compared with the trends (see Figure 1). Therefore, we did not use this observation in the remaining analysis and the determination of the time delay.

In Sredzińska et al. (2017) the change of the line position was discussed in much detail, since during the first three years the Mg ii line seems to move systematically toward longer wavelength with a surprising speed of 104 ± 14 km s⁻¹ yr⁻¹ with respect to the quasar rest frame. However, now this trend has seemingly stopped, and in recent observations it seems to have reversed. Such emission line behavior is frequently considered as a signature of a binary black hole (e.g., Popović 2012 for a review). However, to claim such a phenomenon would require extensive tests that are beyond the current paper, which is aimed at measurement of the time delay of the Mg ii line.

Combining all the seven methods listed in Table 1, the mean value in the observer's frame is ${\overline{\tau }}_{\mathrm{obs}}={658}_{-31}^{+29}$ days. We visually compare the continuum light curve and the original as well as the time-shifted light curve of the Mg ii line in Figure 3. For an easier comparison, the Mg ii line is shifted toward larger magnitudes by the difference in the mean values of the two light curves (1.74 mag). The correlation between the continuum and the shifted Mg ii light curve, although present, is not visually improved with respect to the zero time lag. This is not unexpected since our source does not exhibit such a large variability amplitude as sources with a low Eddington ratio. In addition, even when the line-emission light curve is shifted by the fiducial time delay, it can intrinsically exhibit periods when the line emission is decorrelated with respect to the continuum emission, which is referred to as the emission-line or the BLR holiday (studied in more detail for NGC 5548, Dehghanian et al. 2019). This justifies the need to use several robust statistical methods to assess the best time delay. We also tried to subtract a linear trend from both light curves, but since for both of them the slope is consistent with zero within fitting uncertainties, it did not yield an improvement. Even for a noticeable linear trend, as for instance for HE 0413-4031 (Zajaček et al. 2020), detrending actually led to a decrease in the correlation coefficient at the fiducial time delay. Hence, the subtraction of a linear trend should be tested on a larger set of sources to assess statistical relevance in terms of time-delay determination.

Zoom In Zoom Out Reset image size

Subsequently, we obtain the mean rest-frame value of ${\overline{\tau }}_{\mathrm{rest}}={\overline{\tau }}_{\mathrm{obs}}/(1+z)={296}_{-14}^{+13}$ days for the redshift of z = 1.2231. The light-travel distance of the Mg ii emission zone can then be estimated as ${R}_{\mathrm{Mg}{\rm\small{II}}}\sim c{\overline{\tau }}_{\mathrm{rest}}\,={7.7}_{-0.4}^{+0.3}\times {10}^{17}\,\mathrm{cm}={0.249}_{-0.012}^{+0.011}\,\mathrm{pc}$. The rest-frame time delay is comparable within uncertainties to the time delay of the previously analyzed highly accreting quasar HE 0413-4031 (z = 1.38, Zajaček et al. 2020).

Taking into account the rest-frame time delay of ${\overline{\tau }}_{\mathrm{rest}}={296}_{-14}^{+13}$ days and the Mg ii FWHM in the mean spectrum of ${{\rm{FWHM}}}_{{\rm{Mg}}{\rm\small{II}}}={3695}_{-21}^{+21}\,{\rm{km}}\,{{\rm{s}}}^{-1}$, we can estimate the central black hole mass under the assumption that Mg ii emission clouds are virialized. Using the anticorrelation between the virial factor and the line FWHM according to Mejía-Restrepo et al. (2018),

$$\begin{eqnarray}&&{f}_{\mathrm{MR}}={\left(\displaystyle \frac{{\mathrm{FWHM}}_{\mathrm{Mg}{\rm\small{II}}}}{3200\pm 800\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{-1.21\pm 0.29},\end{eqnarray} \tag{ 3 }$$

we obtain f_MR ∼ 0.84. The virial black hole mass can then be calculated using the Mg ii FWHM in the mean spectrum and the measured time delay as ${M}_{\mathrm{vir}}^{\mathrm{FWHM}}({f}_{\mathrm{MR}})=({6.6}_{-0.3}^{+0.3})\,\times {10}^{8}\,{M}_{\odot }$. Adopting the virial factor according to Woo et al. (2015), f_Woo = 1.12 (based on FWHM of the Hβ line), we obtain the virial black hole mass ${M}_{\mathrm{vir}}^{\mathrm{FWHM}}({f}_{\mathrm{Woo}})=({8.8}_{-0.4}^{+0.4})\,\times {10}^{8}\,{M}_{\odot }$. Here we adopted the FWHM from the mean spectrum since it is better constrained than the rms FWHM. The mean values, however, are consistent within uncertainties, which is in agreement with the general correlation of Mg ii line widths measured in the mean and the rms spectra using the SDSS-RM sample (Wang et al. 2019).

Using instead the Mg ii line dispersion in the mean spectrum, $\sigma ={2707}_{-6}^{+10}\,\mathrm{km}\,{{\rm{s}}}^{-1}$, and the associated virial factor f_σ = 4.47 (based on the Hβ line dispersion) according to Woo et al. (2015), the black hole mass is estimated as ${M}_{\mathrm{vir}}^{\sigma }({f}_{\sigma })=({1.89}_{-0.09}^{+0.08})\times {10}^{9}\,{M}_{\odot }$. This value is consistent with the value inferred from the broadband SED fitting using the model of a thin accretion disk according to Sredzińska et al. (2017), where they obtained M_SED = 2.2 × 10⁹ M_⊙. Hence for our source, using the line dispersion inferred from the mean spectrum (the rms spectrum is too noisy to reliably measure σ) appears to be beneficial for constraining the virial SMBH mass. Using the FWHM yields a virial mass below that inferred from broadband fitting.

Next, we estimate the Eddington ratio. Using our continuum light curve, the mean V-band magnitude is 17.15 ± 0.09 mag. With the redshift of z = 1.2231 and the mean Galactic foreground extinction in the V band of 0.045 mag according to NED, ¹⁴ we determine the luminosity at 3000 Å, $\mathrm{log}({L}_{3000}[\mathrm{erg}\,{{\rm{s}}}^{-1}])={46.359}_{-0.034}^{+0.038}$, for which we apply the conversions of Kozłowski (2015). To estimate the bolometric luminosity, the corresponding bolometric correction can be obtained via the simple power-law scaling by Netzer (2019), ${\kappa }_{\mathrm{bol}}=25\times {[{L}_{3000}/{10}^{42}\mathrm{erg}{{\rm{s}}}^{-1}]}^{-0.2}\sim 3.36$, which yields L_bol = 3.36L₃₀₀₀ ∼ 7.68 × 10⁴⁶ erg s⁻¹. The Eddington limit can be estimated for M_• ≃ 2 × 10⁹ M_⊙, which was obtained via the SED fitting as well as the virial mass using the line dispersion, L_Edd ≃ 2.5 × 10⁴⁷ × (M_•/(2 × 10⁹ M_⊙)) erg s⁻¹. Finally, the Eddington ratio is η = L_bol/L_Edd ∼ 0.31, which is about a factor of two smaller than the Eddington ratio of 0.58 obtained by Sredzińska et al. (2017) using the SED fitting. Still, the source is highly accreting, with η comparable to HE 0413-4031 (η = 0.4, Zajaček et al. 2020).

The highly accreting sources exhibit the largest scatter with respect to the standard RL relation, with a trend of shorter time delays by a factor of a few than expected based on their monochromatic luminosity. This was studied for the Hβ reverberation-mapped sources (Du et al. 2015, 2016; Martínez-Aldama et al. 2019), and confirmed for the higher-redshift Mg ii reverberation-mapped sources as well (Homayouni et al. 2020; Martínez-Aldama et al. 2020; Zajaček et al. 2020), which suggests a common mechanism for time-delay shortening driven by the accretion rate. The relation between the rest-frame time delay of the Hβ line and the linear combination of the monochromatic luminosity at 5100 Å and the relative strength of the Fe II line (flux ratio between Fe II and Hβ lines) yields a low scatter of only σ_rms = 0.196 dex (Du & Wang 2019), which suggests that the relative Fe II strength and hence the accretion rate can account for most of the scatter.

In addition, highly accreting Mg ii reverberation-mapped sources exhibit a generally lower scatter along the multidimensional RL relations (Martínez-Aldama et al. 2020). This motivates us to further study how high-luminosity and highly accreting sources such as HE 0435-4312 affect the scatter along the radius–luminosity relation on their own as well as in combination with independent observables, such as the Mg ii line FWHM, the Fe II strength R_{Fe II} , and the fractional variability F_var of the continuum, where the latter two parameters are correlated with the accretion rate.

The high luminosity of HE 0435-4312 is expected to be beneficial for constraining the Mg ii-based radius–luminosity relation. We add HE 0435-4312 to the original sample of 68 Mg ii reverberation-mapped sources studied in Martínez-Aldama et al. (2020). To characterize the accretion rate of our source and in order to compare it with other sources, we make use of the dimensionless accretion rate $\dot{{ \mathcal M }}$ expressed specifically for 3000 Å as (Wang et al. 2014; Martínez-Aldama et al. 2020)

$$\begin{eqnarray}&&\dot{{ \mathcal M }}=26.2{\left(\displaystyle \frac{{L}_{44}}{\cos \theta }\right)}^{3/2}{m}_{7}^{-2},\end{eqnarray} \tag{ 4 }$$

where L₄₄ is the luminosity at 3000 Å expressed in units of 10⁴⁴ erg s⁻¹ and m₇ is the central black hole mass expressed in units of 10⁷ M_⊙. The angle θ is the inclination angle with respect to the accretion disk and we set $\cos \theta =0.75$, which represents the mean inclination angle for type I AGNs according to studies of the covering factor of the dusty torus (see, e.g., Lawrence & Elvis 2010; Ichikawa et al. 2015).

Using Equation (4) and the estimate of black hole mass based on the line dispersion in the mean spectrum and the broadband SED fitting, we obtain $\dot{{ \mathcal M }}={4.0}_{-0.6}^{+0.7}$ or $\mathrm{log}\dot{{ \mathcal M }}\,={0.61}_{-0.06}^{+0.07}$ for HE 0435-4312, which puts it into the high-accretor category according to Martínez-Aldama et al. (2020), where all the Mg ii reverberation-mapped sources with $\mathrm{log}\dot{{ \mathcal M }}\gt 0.2167$ belong (median value of their Mg ii sample). Using the smaller SMBH mass based on Mg ii FWHM in the mean spectrum yields a larger $\dot{{ \mathcal M }}$ by a factor of a few, $\mathrm{log}\dot{{ \mathcal M }}={1.51}_{-0.06}^{+0.07}$ (based on f_MR) and $\mathrm{log}\dot{{ \mathcal M }}={1.26}_{-0.06}^{+0.07}$ (based on f_Woo). In the further analysis, we adopt the $\dot{{ \mathcal M }}$ value based on the line dispersion in the mean spectrum because of the consistency of ${M}_{\mathrm{vir}}^{\sigma }({f}_{\sigma })$ with the SMBH mass inferred from the broadband fitting of Sredzińska et al. (2017). In Figure 4 (left panel), we plot the RL relation for all 69 sources including HE 0435-4312 (green circle). In the right panel of Figure 4, we restrict the RL relation to only highly accreting sources, which results in a significantly reduced rms scatter of σ_rms = 0.1991 dex versus σ_rms = 0.2994 dex for the full sample. Adding HE 0435-4312 results in a small but detectable reduction of scatter for both the full sample (0.2994 versus 0.3014) and the highly accreting subsample (0.1991 versus 0.2012) with respect to the original Mg ii sample of 68 sources (Martínez-Aldama et al. 2020).

Zoom In Zoom Out Reset image size

In the extended RL relations that include other independent observables in the linear combination of logarithms, namely Mg ii line FWHM, the relative Fe II strength with respect to the Mg ii line R_{Fe II} , and the continuum fractional variability F_var, HE 0435-4312 lies within 1σ of the mean relation for both the full sample and the highly accreting subsample; see Figure 5 for the combination with FWHM, Figure 6 that includes R_{Fe II} , and Figure 7 that utilizes the continuum F_var. When one restricts the analysis to the highly accreting subsample, the rms scatter drops below 0.2 dex in all combinations, with the smallest scatter exhibited by the RL relation including R_{Fe II} (σ_rms = 0.1734 dex). In Table 2, we summarize the best-fit parameters as well as the rms scatter for all the RL relations for the highly accreting subsample including HE 0435-4312. Our high-luminosity source is beneficial for the reduction of rms scatter, except for the combination including R_{Fe II} , where the rms scatter marginally increases in comparison with the best-fit result of Martínez-Aldama et al. (2020). For all RL relations, adding the new quasar leads to an increase in the Pearson correlation coefficient. For the comparison of the rms scatter and the correlation coefficients with and without the inclusion of HE 0435-4312, see the last four columns of Table 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2. Parameters for The Standard RL Relation as Well as Multidimensional RL Relations Including Independent Observables, Namely FWHM, R_{Fe II} , and F_var, for the Highly Accreting Subsample (35 Sources Including HE 0435-4312)

${\rm{log}}{\tau }_{{\rm{obs}}}$	K1	K2	K3	σrms [dex]	r	σrms [dex] (without)	r (without)
${K}_{1}\mathrm{log}{L}_{44}+{K}_{2}$	0.422 ± 0.055	1.374 ± 0.082	⋯	0.1991	0.80	0.2012	0.78
${K}_{1}\mathrm{log}{L}_{44}+{K}_{2}{\mathrm{logFWHM}}_{3}+{K}_{3}$	0.43 ± 0.06	−0.13 ± 0.31	1.44 ± 0.17	0.1986	0.80	0.2007	0.78
${K}_{1}\mathrm{log}{L}_{44}+{K}_{2}\mathrm{log}{R}_{\mathrm{Fe}{\rm\small{II}}}+{K}_{3}$	0.45 ± 0.05	0.84 ± 0.29	1.28 ± 0.08	0.1734	0.85	0.1718	0.84
${K}_{1}\mathrm{log}{L}_{44}+{K}_{2}\mathrm{log}{F}_{\mathrm{var}}+{K}_{3}$	0.39 ± 0.06	−0.38 ± 0.18	0.99 ± 0.19	0.1861	0.83	0.1863	0.82

Note. The last four columns show rms scatter and Pearson's correlation coefficient with and without HE 0435-4312 for comparison.

Download table as:  ASCIITypeset image

The dimensionless accretion rate $\dot{{ \mathcal M }}$ of Mg ii sources defined by Equation (4) is intrinsically correlated with the rest-frame time delay via the black hole mass ($\dot{{ \mathcal M }}\propto {\tau }^{-2}$, while the Eddington ratio η ∝ τ⁻¹) as discussed by Martínez-Aldama et al. (2020). Therefore we do not use $\dot{{ \mathcal M }}$ in extended RL relations, because the correlation would be artificially enhanced. We merely use $\dot{{ \mathcal M }}$ for the division of the sample into the high and low accretors. For the extended RL relations, we prefer the independent observables Mg ii FWHM, UV R_{Fe II} , and F_var. Characterizing the accretion-rate intensity by $\dot{{ \mathcal M }}$ is justified by its correlation with the relative UV Fe II strength R_{Fe II} , which in turn is related to the accretion rate (Dong et al. 2011; Martínez-Aldama et al. 2020), which was analogously shown for optical Fe II strength (Shen & Ho 2014; Du & Wang 2019). For the whole sample of Mg ii sources, including HE 0435-4312, for which R_{Fe II} can be estimated (in total 66 sources), the Spearman correlation coefficient between $\dot{{ \mathcal M }}$ and R_{Fe II} is positive with ρ = 0.440 (with a p-value of 2.19 × 10⁻⁴). For the same sample, the relative Fe II strength is in turn positively correlated with the Eddington ratio η = L_bol/L_Edd with ρ = 0.447 (p = 1.71 × 10⁻⁴). This justifies the division into low and high accretors according to $\dot{{ \mathcal M }}$.

On the other hand, $\dot{{ \mathcal M }}$ is anticorrelated with respect to the continuum variability F_var with ρ = −0.480 (p = 2.97 × 10⁻⁵, 69 sources). Furthermore, the parameter $\dot{{ \mathcal M }}$ is anticorrelated with Mg ii FWHM, although the (anti)correlation is weaker than for R_{Fe II} and F_var with ρ = −0.386 (p = 1.05 × 10⁻³, 69 sources). The anticorrelation between Mg ii FWHM and $\dot{{ \mathcal M }}$, and hence also R_{Fe II} , is expected from optical studies of eigenvector 1 of the quasar main sequence (Boroson & Green 1992; Sulentic et al. 2000), which can be traced in the UV plane as well (Kozłowski et al. 2020). These (anti)correlations between $\dot{{ \mathcal M }}$ and other quantities across classical and extended RL relations are also apparent in the color-coded plots in Figures 4–7. The reduction in scatter for the high-$\dot{{ \mathcal M }}$ subsample could therefore arise due to a lower variability and the stabilisation of the luminosity and black hole mass ratio for the sources accreting close to the Eddington limit as discussed in detail in Martínez-Aldama et al. (2020); see also Ai et al. (2010) and Marziani & Sulentic (2014).

In order to clarify our idea of the sample division quantitatively, we divided the sample analyzed by Martínez-Aldama et al. (2020), with the studied source HE 0435-4312 included, considering the median values of R_{Fe II} and F_var ($\mathrm{log}{R}_{\mathrm{Fe}{\rm\small{II}}}=0.02$, $\mathrm{log}{F}_{\mathrm{var}}=-1.05$) instead of $\dot{{ \mathcal M }}$ ($\mathrm{log}\dot{{ \mathcal M }}\,=0.23$). Since the correlation between $\dot{{ \mathcal M }}$ and R_{Fe II} is positive, we expect that the high Fe II emitters show the highest $\dot{{ \mathcal M }}$ values, i.e., above the corresponding median value. On the other hand, the correlation between F_var and $\dot{{ \mathcal M }}$ is negative, hence highly accreting AGNs show a lower variability or F_var values with respect to the median value of each parameter. In the case of the $\dot{{ \mathcal M }}$–R_{Fe II} correlation, we find that ∼70% of the sample satisfies the criteria with respect to the median values of R_{Fe II} and $\dot{{ \mathcal M }}$. Meanwhile, for the $\dot{{ \mathcal M }}$–F_var correlation, ∼72% of our sources are within the median values considered. This is graphically shown in Figure 8, where we depict the correlation $\dot{{ \mathcal M }}$–R_{Fe II} (left panel) and the anticorrelation $\dot{{ \mathcal M }}$–F_var (right panel) for the studied sample of Mg ii sources. Therefore, we can conclude that the division into low and high accretors based on $\dot{{ \mathcal M }}$ is analogous to a division based on the independent parameters R_{Fe II} and F_var. To illustrate this, in Figure 4 we use R_{Fe II} for color-coding individual sources in the RL relation for the whole sample (left panel) and the high-accretor subsample (right panel).

Zoom In Zoom Out Reset image size

The fractional variability of the Mg ii line of HE 0435-4312 is ${F}_{\mathrm{var}}^{\mathrm{line}}\simeq 5.4 \%$ (excluding observation 19). The continuum fractional variability is larger, ${F}_{\mathrm{var}}^{\mathrm{cont}}\simeq 8.8 \%$ over 14.74 yr. However, when the continuum light curve is separated into the first 7 yr (CATALINA), ${F}_{\mathrm{var}}^{\mathrm{cont}1}\simeq 4.8 \%$, and the last 6 yr (SALT monitoring), ${F}_{\mathrm{var}}^{\mathrm{cont}2}\simeq 4.8 \%$, then both values are comparable to the fractional variability of Mg ii emission. This implies that Mg ii-emitting clouds reprocess the continuum emission very well for our source. In other words, the triggering continuum and a line echo have similar amplitudes, suggesting that sharp echoes are present and that the Mg ii-emitting region lies on (or close to) an isodelay surface, which is an important geometrical condition. In addition, the variability timescale of the continuum must be long enough, i.e., longer than the light-travel time through the locally optimally emitting cloud (LOC) model of the BLR. It seems that both conditions are fulfilled for our source.

This is in contrast with the study of Guo et al. (2020), who used the CLOUDY code and the LOC model to study the response of Mg ii emission to the variable continuum. They found that at the Eddington ratio of ∼0.4, the Mg ii emission saturates and does not further increase with the rise in the continuum luminosity. Observationally, Yang et al. (2020) found for extreme-variability quasars that Mg ii emission responds to the continuum but with a smaller amplitude, ${\rm{\Delta }}\mathrm{log}L(\mathrm{Mg}\,{\rm\small{II}})=(0.39\pm 0.07){\rm{\Delta }}\mathrm{log}L(3000\,\mathring{\rm A} )$. Although our source has a high Eddington ratio comparable to that studied in Guo et al. (2020), η ∼ 0.3–0.6, its Mg ii emission responds very efficiently to the variable continuum. For the previous highly accreting luminous quasar HE 0413-4031 that we studied (Zajaček et al. 2020), we also showed that the Mg ii line can respond strongly to the increase in the continuum, ${\rm{\Delta }}\mathrm{log}L(\mathrm{Mg}\,{\rm\small{II}})=(0.82\pm 0.26){\rm{\Delta }}\mathrm{log}L(3000\,\mathring{\rm A} )$, as shown by the intrinsic Baldwin effect.

A possible interpretation of the discrepancy between the saturation of Mg ii emission at larger Eddington ratios modeled nominally by Guo et al. (2020) and our two highly accreting sources, HE 0435-4312 and HE 0413-4031, which exhibit a strong response of the Mg ii emission to the continuum, is an order-of-magnitude difference in the black hole mass as well as in the studied luminosity at 3000 Å. In Guo et al. (2020), they used M_• = 10⁸ M_⊙ and L₃₀₀₀ = 10⁴⁴–10⁴⁵ erg s⁻¹ as well as (R_out, Γ) = (10^17.5, −2) for the outer radius and the slope of the LOC model, respectively. For our source, all of the characteristic parameters are increased: M_• ≃ 2 × 10⁹ M_⊙, L₃₀₀₀ ≃ 10^46.36 erg s⁻¹, and R_{Mg II} = 10^17.9 cm. Mainly the larger outer radius implies that not all of the Mg ii-emitting gas is fully ionized at these scales and can exhibit "partial breathing," as is also shown by Guo et al. (2020) for the case with R_out = 10¹⁸ cm, when the Mg ii line luminosity continues to rise with the increase in the continuum luminosity.

Since the RL relation for a highly accreting Mg ii subsample showed a relatively small scatter of ∼0.2 dex, we attempt to use it for cosmological purposes. We adopt the general approach outlined in Martínez-Aldama et al. (2019), that is we assume a perfect relation between the absolute luminosity and the measured time delay, and having absolute luminosity and the observed flux we can determine the luminosity distance for each source. We do not use here the best-fit relation given in Figure 4, right panel, since this relation was calibrated for a specific assumed cosmology. Therefore we assume a relation in a general form,

$$\begin{eqnarray}&&\mathrm{log}{L}_{3000}=\alpha \mathrm{log}\tau -\beta +44,\end{eqnarray} \tag{ 5 }$$

and we treat α and β as free parameters. We consider first the case of a flat cosmology, and we assume the value of the Hubble constant H₀ = 67.36 km s⁻¹ Mpc⁻¹ adopted from Planck Collaboration et al. (2020), and we minimize the fits to the predicted luminosity distance by varying Ω_m, α, and β. The preliminary fit was highly unsatisfactory, and we used a sigma-clipping method to remove the outliers (eight sources). The final sample was well fitted with the standard ΛCDM model for ${{\rm{\Omega }}}_{{\rm{m}}}={0.252}_{-0.044}^{+0.051}$ (1σ error). This result is fully consistent with the data from Planck Collaboration et al. (2020) with ${{\rm{\Omega }}}_{{\rm{m}}}^{\mathrm{Planck}}=0.3153\pm 0.0073$. The values obtained for the remaining two parameters were α = 1.8 and β = 2.1, which would correspond to a slightly different RL relation, $\mathrm{log}\tau =0.56\mathrm{log}{L}_{3000}+1.18$. The results are illustrated in Figure 9, left panel. We also attempted to constrain both Ω_Λ and Ω_m, waiving the assumption of a flat space. The parameters α and β in Equation (5) were kept fixed to the values inferred from the flat-cosmology fit, i.e., α = 1.8 and β = 2.1. The result is given in the right panel of Figure 9. The best fit in this case implies a slightly lower Ω_m, but the contour error is large and covers the best fit from Planck Collaboration et al. (2020) within the 1σ confidence level. The data do not yet tightly constrain the cosmological parameters, but the Mg ii highly accreting subsample does not imply any departures from the standard model, and the data quality is much higher, and constrains the parameters better than the larger sample based on the Hβ line and discussed in Martínez-Aldama et al. (2019), and somewhat better than the mixed sample discussed in Czerny et al. (2020).

Zoom In Zoom Out Reset image size

We first analyzed the continuum and Mg ii line-emission light curves using the interpolated cross-correlation function, which is a standard and well-tested method for assessing the time delay between the continuum and the line emission flux density (Peterson et al. 1998). Light curves are typically unevenly sampled, while the ICCF works with the continuum–line emission pairs with a certain time step of Δt = t_i+1 − t_i . The cross-correlation function for two light curves, x_i and y_i , with the same time step of Δt achieved by interpolation, evaluated for a time shift of τ_k = kΔt (k = 1, ..., N − 1), is defined as

$$\begin{eqnarray}&&\begin{array}{l}\mathrm{CCF}({\tau }_{k})\\ \,=\displaystyle \frac{(1/N){\displaystyle \sum }_{i=1}^{N-k}({x}_{i}-\overline{x})({y}_{i+k}-\overline{y})}{{\left[(1/N){\displaystyle \sum }_{i=1}^{N}{\left({x}_{i}-\overline{x}\right)}^{2}\right]}^{1/2}{\left[(1/N){\displaystyle \sum }_{i=1}^{N}{\left({y}_{i}-\overline{y}\right)}^{2}\right]}^{1/2}}.\end{array}\end{eqnarray} \tag{ C1 }$$

The same time step Δt can be achieved by interpolating the continuum light curve with respect to the line-emission light curve and vice versa (asymmetric ICCF). Typically, both interpolations are averaged to obtain the symmetric ICCF.

For the time-delay analysis, we use the Python implementation of ICCF in code PYCCF by Sun et al. (2018), which is based on the algorithm by Peterson et al. (1998). The code allows one to perform both asymmetric and symmetric interpolation. Based on the Monte Carlo techniques of random subset selection and flux randomization, one can obtain the centroid and the peak distributions and their corresponding uncertainties.

We studied the time delay between the continuum light curve consisting of 81 measurements, with eight Catalina-survey averaged detections, 16 SALTICAM measurement, 27 BMT data, and 30 OGLE data. Two SALTICAM measurements were excluded based on their poor quality. The emission-line light curve consists of 25 SALT measurements, where the 19th measurement has a poor quality due to weather conditions and is excluded from the further analysis. Hence, we have 79 continuum points with a mean cadence of 69.0 days and 24 Mg ii flux-density measurements with a mean cadence of 121.6 days. We set the interpolation interval to one day. For the redshift of z = 1.2231 and the d11_mod template, we display the ICCF as a function of time delay in the observer's frame in Figure 12 (left panel) for both asymmetric and symmetric interpolation. In the middle and right panels, we show the centroid and the peak distributions for the symmetric interpolation based on 3000 Monte Carlo realizations of random subset selection and flux randomization. The centroid and the peak are generally not well defined. The peak value of the correlation function for the symmetric interpolation is 0.32 for a time delay of 635 days, which is less than for our previous quasars CTS C30.10 (peak CCF of ∼0.65, Czerny et al. 2019; Zajaček et al. 2019) and HE 0413-4031 (peak CCF of 0.8, Zajaček et al. 2020). In the next step, we focus on the surroundings of this peak and we analyze the CCF centroid and peak distributions in the interval between 500 and 1000 days. The results for all interpolations are summarized in Table 5. For the symmetric interpolation, we obtain a centroid time delay of ${\tau }_{\mathrm{cent}}={663}_{-40}^{+66}$ days and a peak time delay of ${\tau }_{\mathrm{peak}}={672}_{-37}^{+49}$ days.

Zoom In Zoom Out Reset image size

Table 5. Summary of the Time-delay Determination for the Quasar HE 0435-4312 Using the ICCF

Interpolation Method	Time Delay (z = 1.2231)
Interpolated continuum—centroid [days]	${635}_{-34}^{+31}$
Interpolated continuum—peak [days]	${630}_{-46}^{+42}$
Interpolated line—centroid [days]	${673}_{-50}^{+60}$
Interpolated line—peak [days]	${683}_{-48}^{+38}$
Symmetric—centroid [days]	${663}_{-40}^{+66}$
Symmetric—peak [days]	${672}_{-37}^{+49}$

Note. We list the centroids as well as the peaks for the interpolated continuum light curve, the interpolated line emission light curve, and the symmetric interpolation. Time delays are expressed in days in the observer's frame of reference.

Download table as:  ASCIITypeset image

We also performed the time-delay analysis using the ICCF for the Mg ii light curve derived using the KDP15 template (case f in Table 4). For the symmetric ICCF, the CCF peak value is 692 days with CCF = 0.36, which is comparable to the analysis performed for the d11_mod template. The CCF peak and centroid distributions also contain the same subpeaks, with the peak close to 700 days being the most prominent; see Figure 13.

Zoom In Zoom Out Reset image size

For unevenly sampled and sparse light curves, the ICCF can distort the best time delay by adding additional data points to one or both light curves. In that case, the discrete correlation function (DCF) introduced by Edelson & Krolik (1988) is better suited to determine the best time lag.

In general, the DCF value is determined for any pair of light curves (x_i , y_j ) whose time difference Δt_ij falls into the time-delay bin of size δ τ, τ − δ τ/2 < Δt_ij < τ + δ τ/2, where τ is a given time delay. First, one calculates the unbinned DCF,

$$\begin{eqnarray}&&{\mathrm{UDCF}}_{{ij}}=\displaystyle \frac{({x}_{i}-\overline{x})({y}_{j}-\overline{y})}{\sqrt{{s}_{x}-{\sigma }_{x}^{2}}\sqrt{{s}_{y}-{\sigma }_{y}^{2}}},\end{eqnarray} \tag{ C2 }$$

where $\overline{x}$ and $\overline{y}$ are light-curve means in a given time-delay bin, s_x and s_y are corresponding light curve variances, and σ_x and σ_y are the mean measurement errors of the light-curve points in the given time-delay interval.

The DCF is then calculated as a mean over M light-curve points that are located in a given bin,

$$\begin{eqnarray}&&\mathrm{DCF}(\tau )=\displaystyle \frac{1}{M}\displaystyle \sum _{{ij}}{\mathrm{UDCF}}_{{ij}}.\end{eqnarray} \tag{ C3 }$$

The uncertainty can be estimated using the relation

$$\begin{eqnarray}&&{\sigma }_{\mathrm{DCF}}=\displaystyle \frac{1}{M-1}\sqrt{\sum {[{\mathrm{UDCF}}_{{ij}}-\mathrm{DCF}(\tau )]}^{2}}.\end{eqnarray} \tag{ C4 }$$

We made use of the Python script pyDCF (Robertson et al. 2015), where the general procedure described using Equations (C2), (C3), and (C4) is implemented. We searched for the DCF peak in the time interval from 100 to 1000 days, with a time step of 50 days. We obtained a better defined DCF peak using the Gaussian weighting scheme than using the slot weighting. In Figure 14 (left panel), we show the DCF versus time delay in the observer's frame. The time delay with the highest DCF value of 0.41 ± 0.16 is at 675 days. To determine the uncertainty of the peak, we ran 1000 bootstrap simulations, where at each step we constructed a pair of new light curves by randomly selecting a light-curve subset. The histogram of peak DCF time delays constructed from 1000 bootstrap realizations is shown in Figure 14 (right panel). The peak of the distribution is at ${\tau }_{\mathrm{DCF}}={656}_{-73}^{+18}$ days in the observer's frame. The left and right uncertainties represent standard deviations, where we included the area within 30% of the main peak.

Zoom In Zoom Out Reset image size

In the following, we applied the z-transformed discrete correlation function (Alexander 1997), which improves the classical DCF by replacing equal time binning by equal population binning. This is achieved by applying Fisher's z-transform. With this property, the zDCF processes satisfactorily especially undersampled, sparse, and heterogeneous light curves, which to some extent is the case here, since the continuum light curve originates from three different instruments and the Mg ii light curve is relatively sparse with respect to the continuum (24 versus 79 points). In our analysis, we applied the zDCF method several times, changing the minimum number of light-curve pairs per population bin. Finally, we set the minimum number of light-curve pairs to eight and used 5000 Monte Carlo-generated pairs of light curves to determine the errors in each bin. In Figure 15, we show the zDCF values as a function of time delay in the observer's frame. There are some peaks with a large DCF value within the first 200 days in the rest frame, e.g., 190.3 days (DCF = 0.69). These are, however, too short to correspond to the realistic rest-frame time delay for our highly luminous quasar. The most prominent peak at larger time delays is at 646 days (DCF = 0.46).

Zoom In Zoom Out Reset image size

Next we calculate the maximum likelihood (ML) using the zDCF values and we focus on the prominent peak in the interval between 500 and 800 days. From the ML analysis, we obtain the ML peak of ${\tau }_{\mathrm{zDCF}}={646}_{-57}^{+63}$ days in the observer's frame. This peak and the corresponding uncertainties are also highlighted in Figure 15 by vertical lines.

Another technique to evaluate the time-delay distribution is to model the continuum variability of AGNs as a stochastic process using the damped random walk process (Kelly et al. 2009; Kozłowski et al. 2010; MacLeod et al. 2010; Kozłowski 2016). Subsequently, the line emission is assumed to be time-delayed, smoothed, and a scaled version of the continuum emission. This method is implemented in the JAVELIN package (Just Another Vehicle for Estimating Lags in Nuclei; Zu et al. 2011, 2013, 2016) ¹⁵ . The package uses Markov Chain Monte Carlo (MCMC) methods to first determine the posterior probabilities for the continuum variability timescale and the amplitude. Based on that, the variability of the line emission is modeled to obtain the posterior probabilities for the time lag, smoothing width of the top-hat function, and the scaling ratio (the ratio between the line and the continuum amplitudes, A_line/A_cont).

First, we searched for the time delay in the longer interval between 0 and 2000 days in the observer's frame. There are four distinct peaks (see Figure 16, left panel) at ∼100, ∼650, ∼1300, and ∼2000 days. The first peak is too short, while the last two appear too long for our data set. In the zDCF analysis in Section C.3, the time-delay peaks at 1000 days and more also have large uncertainties. Therefore, in the next search, we perform a time-delay search in the narrower interval between 0 and 1000 days to focus on the intermediate peak at ∼650 days. We obtain a prominent peak at ∼652 days; see Figure 16 (right panel).

Zoom In Zoom Out Reset image size

To determine the precise position of the time-delay peak, we performed 200 bootstrap simulations—i.e., generating 200 subsets of the original continuum and Mg ii light curves. Then we applied the JAVELIN to determine the peak time delay for each individual pair of light curves. From 200 best time delays, we first construct a density plot in the plane of the time delay and the scaling factor; see Figure 17 (left panel). We see that there is a peak in the distribution close to 600 days in the observer's frame. In Figure 17 (right panel), we show the histogram of the best time delays with the peak value of ${\tau }_{\mathrm{peak}}={645}_{-41}^{+55}$ days, where the uncertainties were determined within 30% of the peak value.

Zoom In Zoom Out Reset image size

For the JAVELIN time-delay determination, we typically notice several narrow comparable peaks in the histogram of time delays for a single MCMC run as can be seen in Figure 16. These secondary peaks typically arise due to a limited duration of observational runs and a sparse and/or nonuniform sampling of continuum and line-emission light curves. However, they might also reflect the complex BLR geometry seen in advanced data analysis or models (e.g., Grier et al. 2017; Hu et al. 2020; Horne et al. 2021; Naddaf et al. 2021). As mentioned in Grier et al. (2017), these alias peaks can be of a comparable significance to the true peak in the time-delay distribution, which is also the case here; see Figure 18 (right panel, black line), where time delays at ∼105, ∼1300, and ∼1980 days appear to have a larger significance than the fiducial peak of ∼650 days in the observer's frame.

Zoom In Zoom Out Reset image size

Following Grier et al. (2017), we weight the posterior probability distribution of time delays by the number of overlapping light-curve points, where the prior probability distribution is given by P(τ) = [N(τ)/N(0)]², where N(τ) denotes the number of overlapping light-curve points for the time delay τ, and N(0) is the number of data points in the overlap region for zero time lag. Hence, the weight is P(τ) = 1 by definition for zero time lag and it will decrease for larger positive and negative time lags, being zero when there is no overlap. For our light curves, we construct the prior probability distribution P(τ), which is shown in Figure 18 (left panel).

Using P(τ), we weight the posterior time-delay distribution based on the number of overlapping light-curve points. This is demonstrated in Figure 18 (right panel, red line). The alias mitigation based on P(τ) helps to effectively suppress longer time-delays, which drop in significance below the fiducial time lag at ∼655 days. However, this technique is not efficient to suppress aliases at shorter time delays of ∼100 days, where the decrease is negligible due to the definition of P(τ). This peak at smaller time delays is effectively mitigated by the bootstrap or random subset selection technique. As we showed in Figure 17, this peak is not further represented in the posterior time-delay distribution, which leaves us with the fiducial peak close to ∼650 days. In addition, in Appendix D we show that the artifact peaks at ≲200 days arise in the probability density distributions (see Figures 22 and 23) constructed from mock light curves that are sampled with the same cadence as the original continuum and line-emission light curves.

We also analyzed the time- delay between the continuum and Mg ii light curves using a novel technique of the measures of the data regularity or randomness (Chelouche et al. 2017), which were previously applied extensively in cryptography and data compression. The advantage of these regularity measures is that they do not require the interpolation as for the ICCF or χ² technique and neither do they require binning in the correlation space like the DCF and zDCF methods. In addition, no AGN variability model is needed for the continuum light curve as for the JAVELIN.

One of the suitable estimators to analyze the time delay between two light curves is the optimized von Neumann scheme, which works with the combined light curve $F(t,\tau )={({t}_{i},{f}_{i})}_{i\,=\,1}^{N}={F}_{1}\cup {F}_{2}^{\tau }$, where F₁ is the continuum light curve and ${F}_{2}^{\tau }$ is the time-shifted line-emission light curve. The optimized von Neumann estimator for a given time delay τ is defined as the mean of the successive differences of F(t, τ),

$$\begin{eqnarray}&&E(\tau )=\displaystyle \frac{1}{N-1}\displaystyle \sum _{i=1}^{N-1}{\left[F({t}_{i})-F({t}_{i+1})\right]}^{2}.\end{eqnarray} \tag{ C5 }$$

The aim of the von Neumann scheme is to find the minimum $E(\tau ^{\prime} )$, where $\tau ^{\prime} \sim {\tau }_{0}$ is supposed to correspond to the actual time delay.

In Figure 19 (left panel), we calculate E(τ) for HE 0435-4312 using the Python script based on Chelouche et al. (2017). ¹⁶ The minimum is reached for a time delay between 300 and 400 days in the observer's frame (333–334 days and 391–392 days), and then for 690 days, which is consistent with the results of previous methods. To determine the peak value and its uncertainty close to this minimum, we performed ∼10,000 bootstrap realizations, from which we constructed the histogram of the best von Neumann time delays; see Figure 19 (right panel). In this histogram, we detect two prominent peaks at ∼27 days and ∼1200 days and a smaller peak at ∼635 days. After focusing on the interval between 150 and 1100 days, the best peak is at ${635}_{-66}^{+32}$ days in the observer's frame, where the uncertainties correspond to standard deviations within 30% of the best peak.

Zoom In Zoom Out Reset image size

For comparison, we also applied the Bartels estimator to both light curves; it uses the ranked version of the combined light curve F_R(t, τ). In Figure 20 (left panel), we show the Bartels estimator as a function of the time delay in the observer's frame. The minimum is clearer than for the von Neumann estimator and is located at 690 days. After running 10,000 bootstrap realizations, we detect three peaks as for the von Neumann estimator. Between 100 and 1100 days, there is a peak of ${644}_{-45}^{+27}$ days (see Figure 20, right panel), where the uncertainties were calculated within 30% of this time-delay peak.

Zoom In Zoom Out Reset image size

Inspired by the time-delay analysis in quasar lensing studies, we applied the χ² method to our set of light curves. The χ² method performed well and consistently in comparison with other time-delay analysis techniques for the previous two SALT quasars—CTS C30.10 (Czerny et al. 2019) and HE 0413-4031 (Zajaček et al. 2020). It also performs better than the classical ICCF for the case when the AGN variability can be interpreted as a red-noise process (Czerny et al. 2013). We subtracted the mean from both the continuum and the line-emission light curves, and subsequently they were normalized by their corresponding variances. The similarity between the continuum and the time-shifted line-emission light curves is evaluated using χ². We make use of symmetric interpolation.

In Figure 21 (left panel), we show χ² as a function of the time delay in the observer's frame. There is a global minimum at 696 days. The exact position and the uncertainty of this peak are inferred from the histogram of the best time delays constructed from 10,000 bootstrap realizations. We obtain ${706}_{-61}^{+70}$ days, where the uncertainties were determined as the left and the right standard deviations from the time delays that are positioned within 30% of the main peak; see Figure 21 (right panel).

Zoom In Zoom Out Reset image size