H0LiCOW VII. Cosmic evolution of the correlation between black hole mass and host galaxy luminosity

TL;DR

This study uses gravitational lensing of active galactic nuclei to investigate how the relationship between supermassive black hole mass and host galaxy luminosity evolves up to redshift 4.5, revealing black hole growth precedes galaxy growth.

Contribution

First to analyze high-redshift lensed AGN to measure black hole-host galaxy luminosity evolution, incorporating lens modeling and new calibrations for a comprehensive redshift range.

Findings

Black holes in distant galaxies are less luminous than local counterparts when accounting for passive evolution.

The black hole to host luminosity ratio increases with redshift, following a power-law with gamma around 0.6-0.8.

Results support black hole growth occurring before significant galaxy assembly.

Abstract

Strongly lensed active galactic nuclei (AGN) provide a unique opportunity to make progress in the study of the evolution of the correlation between the mass of supermassive black holes ($\mathcal M_{BH}$) and their host galaxy luminosity ($L_{host}$). We demonstrate the power of lensing by analyzing two systems for which state-of-the-art lens modelling techniques have been applied to Hubble Space Telescope imaging data. We use i) the reconstructed images to infer the total and bulge luminosity of the host and ii) published broad-line spectroscopy to estimate $\mathcal M_{BH}$ using the so-called virial method. We then enlarge our sample with new calibration of previously published measurements to study the evolution of the correlation out to z~4.5. Consistent with previous work, we find that without taking into account passive luminosity evolution, the data points lie on the local relation. Once passive luminosity evolution is taken into account, we find that BHs in the more distant Universe reside in less luminous galaxies than today. Fitting this offset as $\mathcal M_{BH}$/$L_{host}$ $\propto$ (1+z)$^{\gamma}$, and taking into account selection effects, we obtain $\gamma$ = 0.6 $\pm$ 0.1 and 0.8$\pm$ 0.1 for the case of $\mathcal M_{BH}$-$L_{bulge}$ and $\mathcal M_{BH}$-$L_{total}$, respectively. To test for systematic uncertainties and selection effects we also consider a reduced sample that is homogeneous in data quality. We find consistent results but with considerably larger uncertainty due to the more limited sample size and redshift coverage ($\gamma$ = 0.7 $\pm$ 0.4 and 0.2$\pm$ 0.5 for $\mathcal M_{BH}$-$L_{bulge}$ and $\mathcal M_{BH}$-$L_{total}$, respectively), highlighting the need to gather more high-quality data for high-redshift lensed quasar hosts. Our result is consistent with a scenario where the growth of the black hole predates that of the host galaxy.

Tables

Table 1. Table 3: The inference of HE0435 and RXJ1131.

Object	magnitude	$R_{\mathrm{eff}}$	Sérsic index ($n$)
		(arcsec)
HE0435	$21.75\pm 0.13$	$0.82\pm 0.14$	$3.94\pm 0.14$
RXJ1131disk	$20.07\pm 0.06$	0.84$\pm 0.09$	fixed 1
RXJ1131bulge	$21.81\pm 0.28$	0.20$\pm 0.08$	fixed 4

Table 2. Table 4: The summary for the different inference of γ 𝛾 \gamma .

Sample	Account for	${ℳ}_{\mathrm{BH}}$-$L_{\mathrm{bulge}}$	${ℳ}_{\mathrm{BH}}$-$L_{\mathrm{total}}$
	selection effects?
Entire	No	$0.75\pm 0.11$	$0.95\pm 0.11$
Entire	Yesa	$0.6\pm 0.1$	0.8$\pm 0.1$
Exclude P06, SS13	No	$0.88\pm 0.28$	0.51$\pm 0.28$
Exclude P06, SS13	Yes	$0.7\pm 0.4$	0.2$\pm 0.5$

Table 3. Table 1: Properties of AGNs in the distant sample.

Object	Line(s) Used	redshift	log ${ℳ}_{\mathrm{BH}}$	log $L_{\mathrm{host}}$,R
			$(M_{\odot }$)	$(L_{\odot }$)
			$\pm 0.4$ dex	$\pm 0.2$ dex
HE 0435	MgII	1.693	8.61	10.96
RXJ 1131bulge	MgII/H$\beta$	0.654	8.26/8.41	10.58
RXJ 1131disk	MgII/H$\beta$	0.654	8.26/8.41	11.12
Lensed AGNs from P06
RXJ 1131	H$\beta$	0.66	7.90	11.02
SDSS 1226-0006	MgII	1.12	8.41	10.74
FBQ 0951+2635	MgII	1.24	8.57	10.25
CTQ 414	MgII/CIV	1.29	7.78/8.19	10.87
B 0712+472	MgII	1.34	7.44	10.90
SBS 0909+532	MgII	1.38	9.13	10.54
Q 0957+561	MgII/CIV	1.41	9.06/8.97	11.79
FBQ 1633+3134	MgII/CIV	1.52	8.84/8.91	11.08
SDSS 0924+0219	MgII	1.54	7.61	11.09
B 1030+071	MgII	1.54	8.13	11.06
SDSS 1335+0118	MgII	1.57	8.77	10.97
B 1600+434	MgII	1.59	7.56	10.98
HE 0047-1756	CIV	1.66	8.83	11.24
HE 0435	CIV	1.69	8.36	11.12
PG 1115+080	MgII/CIV	1.72	8.67/8.63	11.08
SBS 1520+530	CIV	1.86	8.60	10.82
HE 2149-2745	CIV	2.03	9.48	11.47
HE 1104-1805	CIV	2.32	9.03	11.48
Q 1017-207	CIV	2.55	8.88	11.69
H 1413+117	CIV	2.55	8.08	11.48
MG 0414+0534	H$\beta$	2.64	9.07	11.41
J 1004+1229	CIV	2.65	8.97	11.70
Q 0142-110	CIV	2.72	9.01	11.37
LBQS 1009-0252	MgII/CIV	2.74	8.51/8.70	11.48
RXJ 0911+0551	CIV	2.80	8.57	10.93
PMNJ 1632-0033	CIV	3.42	8.25	11.35
B 1422+231	MgII/CIV	3.62	8.93/9.34	11.60
BRI 0952-0115	CIV	4.5	8.80	11.95
Non-lensed AGNs from P06
PKS 0440-00	MgII	0.844	8.09	11.29
MGC 2214+3550A	MgII	0.88	8.76	11.13
MGC 2214+3550B	MgII	0.88	8.26	10.59
3C 422	MgII	0.942	9.04	11.64
PKS 0938+18	MgII	0.943	8.53	11.16
SGP5:46	MgII	0.955	8.03	10.92
LBQS 1009-0252c	CIV	1.63	8.88	11.34
RXJ 0921+4528A	MgII/CIV	1.65	9.09/8.71	11.12
RXJ 0921+4528B	MgII/CIV	1.65	8.82/8.4	11.38
SGP4:39	CIV	1.716	8.07	10.49
MZZ 11408	CIV	1.735	8.00	10.82
SGP2:36	CIV	1.756	8.87	11.26
MZZ 1558	CIV	1.829	8.80	10.94
SGP2:25	CIV	1.868	8.45	11.27
MZZ 4935	CIV	1.876	8.02	10.43
SGP3:39	CIV	1.964	8.69	11.38
SGP2:11	CIV	1.976	8.69	11.04
4C 45.51	MgII/CIV	1.992	8.60/8.48	12.18
MZZ 9592	CIV	2.71	8.47	11.49
MZZ 9744	CIV	2.735	8.52	11.09
Non-lensed AGNs from P15				($L_{\mathrm{bulge}}$/ $L_{\mathrm{total}}$)
S09a	H$\beta$	0.3545	7.99	11.03/11.05
S10	H$\beta$	0.3513	8.44	10.60/11.03
S12	H$\beta$	0.3583	8.78	10.11/11.02
S21	H$\beta$	0.3546	8.93	10.38/11.30
S16	H$\beta$	0.3702	8.02	9.89/10.53
S23	H$\beta$	0.3513	8.82	10.27/11.04
S24	H$\beta$	0.3619	8.20	10.95/11.14
S26	H$\beta$	0.3692	8.00	10.47/10.76
S27	H$\beta$	0.3669	7.85	11.00/11.07
S01	H$\beta$	0.3594	8.15	10.41/10.87
S02	H$\beta$	0.3544	8.02	10.43/10.68
S03	H$\beta$	0.3584	8.17	9.64/11.00
S04	H$\beta$	0.3579	8.11	10.51/10.98
S05	H$\beta$	0.3535	8.73	9.75/10.94
S06	H$\beta$	0.3688	7.72	9.53/10.95
S07	H$\beta$	0.352	8.53	10.28/11.02
S08	H$\beta$	0.3586	7.74	10.02/10.89
S11	H$\beta$	0.3559	7.76	10.50/10.90
SS1	H$\beta$	0.3566	7.73	10.13/10.94
SS2	H$\beta$	0.3671	7.56	10.76/10.76
SS5	H$\beta$	0.3735	7.98	10.02/10.66

Table 4. Table 1: —continued.

Non-lensed AGNs from P15				log $L_{\mathrm{bulge}}$/log $L_{\mathrm{total}}$
S31	H$\beta$	0.3506	8.19	10.76/11.00
SS6	H$\beta$	0.3588	7.39	9.67/10.47
SS7	H$\beta$	0.3613	7.67	10.08/10.74
SS8	H$\beta$	0.3655	7.85	9.96/10.92
SS9	H$\beta$	0.3702	7.95	10.82/10.82
SS10	H$\beta$	0.3658	8.10	10.61/10.81
SS11	H$\beta$	0.3731	7.83	10.04/10.83
SS12	H$\beta$	0.3629	8.15	10.54/10.66
SS13	H$\beta$	0.3743	7.69	10.60/10.60
S28	H$\beta$	0.3678	8.12	10.70/10.91
SS14	H$\beta$	0.3706	7.45	10.48/10.48
S29	H$\beta$	0.3574	7.95	10.03/10.69
SS18	H$\beta$	0.3585	7.51	9.81/10.62
W11	H$\beta$	0.565	7.95	10.73/10.73
W22	H$\beta$	0.5652	8.68	11.27/11.27
W12	H$\beta$	0.5623	8.94	10.49/10.98
W20	H$\beta$	0.5761	8.60	11.03/11.03
W16	H$\beta$	0.578	7.86	10.81/10.81
W8	H$\beta$	0.5712	8.74	11.01/11.01
W3	H$\beta$	0.576	8.76	10.28/10.88
SS15	H$\beta$	0.3593	7.44	10.36/10.36
W1	H$\beta$	0.5736	8.84	10.71/10.95
W4	H$\beta$	0.5766	8.28	11.01/11.01
W5	H$\beta$	0.5767	8.29	10.96/10.96
SS3	H$\beta$	0.3566	7.51	10.10/10.78
SS4	H$\beta$	0.3629	7.85	10.81/10.81
W17	H$\beta$	0.5617	8.31	10.12/10.80
W2	H$\beta$	0.572	9.07	10.91/10.91
W10	H$\beta$	0.5711	7.94	10.51/10.77
W14	H$\beta$	0.5617	8.68	10.97/10.97
W9	H$\beta$	0.5654	8.70	10.89/10.96
J033252-275119	MgII	1.227	8.87	10.53/11.28
J033243-274914	MgII	1.900	9.17	11.73/11.73
J033239-274601	MgII	1.220	8.24	11.22/11.22
J033226-274035	MgII	1.031	7.85	10.37/11.46
J033225-274218	MgII	1.617	8.08	11.48/11.48
J033210-274414	MgII	1.615	8.30	11.52/11.52
J033200-274319	MgII	1.037	7.75	10.47/10.47
J033229-274529	MgII	1.218	8.37	11.39/11.39
J123553+621037	MgII	1.371	8.27	10.73/11.66
J123618+621115	MgII	1.021	8.35	10.11/11.37
J123618+621115	MgII	1.450	8.77	11.45/11.45
158b	MgII	0.717	7.28	10.71/10.89
170	MgII	1.065	7.07	10.16/10.46
271	MgII	0.960	7.43	10.46/11.36
273	MgII	0.970	8.23	10.19/10.45
305	MgII	0.544	8.61	11.02/11.14
333	MgII	1.044	7.90	10.45/10.91
339	MgII	0.675	7.95	10.83/11.00
348	MgII	0.569	8.11	10.64/11.17
379	MgII	0.737	9.14	11.03/11.70
413	MgII	0.664	7.05	10.19/10.74
417	MgII	0.837	8.37	10.18/10.87
465	MgII	0.740	8.02	10.67/11.47
516	MgII	0.733	7.93	11.02/11.38
540	MgII	0.622	7.61	10.84/11.30
597	MgII	1.034	8.12	10.80/10.91
712	MgII	0.841	8.65	11.22/11.59

Table 5. Table 2: Properties of local AGNs.

Object	redshift	log ${ℳ}_{\mathrm{BH}}$	log $L_{\mathrm{bulge}}$,R/log $L_{\mathrm{total}}$,R
		$(M_{\odot }$)	$(L_{\odot }$)
			$\pm 0.2$ dex
3C120	0.03301	7.71$\pm$0.21	10.51/10.51
3C390.3	0.0561	8.43$\pm$0.10	10.55/10.55
Ark120	0.03271	8.15$\pm$0.06	10.23/10.69
Mrk79	0.02219	7.69$\pm$0.12	9.62/10.01
Mrk110	0.03529	7.37$\pm$0.11	9.42/9.93
Mrk279	0.03045	7.51$\pm$0.11	9.87/10.32
Mrk335	0.02579	7.12$\pm$0.11	9.61/10.00
Mrk590	0.02639	7.65$\pm$0.07	9.98/10.45
Mrk817	0.03146	7.66$\pm$0.07	9.34/10.32
PG0052+251	0.155	8.54$\pm$0.09	11.34/11.34
PG0804+761	0.1	8.81$\pm$0.05	10.81/10.81
PG0844+349	0.064	7.94$\pm$0.18	10.54/10.54
PG1211+143	0.0809	8.13$\pm$0.13	10.43/10.43
PG1226+023	0.15834	8.92$\pm$0.09	11.64/11.64
PG1229+204	0.06301	7.83$\pm$0.21	10.30/10.70
PG1411+442	0.0896	8.62$\pm$0.14	10.58/10.58
PG1613+658	0.129	8.42$\pm$0.20	11.47/11.47
PG1700+518	0.292	8.86$\pm$0.10	11.51/11.51
PG2130+099	0.06298	7.55$\pm$0.17	9.82/10.42

Equations

$$\log \left(\frac{{\mathcal{M}}_{\mathrm{BH}}}{M_{\odot }}\right)$$

$$\displaystyle\log\big{(}\frac{\mathcal{M}_{\rm BH}}{10^{7}M_{\odot}})=\alpha+\beta\log(\frac{L_{R}}{10^{10}L_{\odot}}).$$

$$d\log L_{\mathrm{R}}/d\log (1+z)=\delta ,$$

$$\log (L_{R,0})=\log (L_R)-\delta \log (1+z).$$

$$\displaystyle\log\big{(}\frac{\mathcal{M}_{\rm BH}}{10^{7}M_{\odot}})$$

$$\Delta \log {\mathcal{M}}_{\mathrm{BH}}=\gamma \log (1+z),$$

Full text

††pubyear: 2017††pagerange: H0LiCOW VII. Cosmic evolution of the correlation between black hole mass and host galaxy luminosity–A

1 Introduction

It is commonly accepted that almost all the galaxies have a supermassive black hole (BH) in their center, whose mass ($\mathcal{M}_{\rm BH}$) is known to be correlated with the host properties. The tight correlations are usually, but not uniquely, explained as the results of their co-evolution (e.g., Magorrian et al., 1998; Ferrarese & Merritt, 2000; Gebhardt et al., 2001; Marconi & Hunt, 2003; Häring & Rix, 2004; Gültekin et al., 2009; Graham et al., 2011; Beifiori et al., 2012; Park et al., 2015; Kormendy & Ho, 2013) (see, however, Peng, 2007; Jahnke & Macciò, 2011, for a different view). A powerful way to explore the origin of this physical coupling and understand the role of active galactic nuclei (AGN) feedback in galaxy formation is to measure the correlations directly at high redshift and determine how and when they emerged and evolved over cosmic time (e.g., Treu et al., 2004; Salviander et al., 2006; Woo et al., 2006; Jahnke et al., 2009; Schramm & Silverman, 2013; DeGraf et al., 2015).

The most common technique used to estimate $\mathcal{M}_{\rm BH}$ beyond the local Universe ($z>0.1$) is the so-called virial method, based on the properties of broad emission lines in type 1 AGN (Shen, 2013; Peterson, 2014). However, the bright source associated with the AGN makes the study of its host galaxy very difficult. Strong gravitational lensing (see, e.g., Courbin et al., 2002; Schneider et al., 2006; Treu, 2010; Treu & Ellis, 2015, for reviews) stretches the host galaxy out from the wings of the bright point source as point spread function (PSF), providing a unique opportunity to infer its magnitude robustly (Peng et al., 2006). However, in order to measure host luminosity ($L_{\rm host}$) and construct the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ correlation from strongly lensed AGN, it is necessary to ensure that any systematic uncertainties associated with the gravitational lens model can be controlled to the desired level of accuracy.

Recently, Ding et al. (2017) studied the fidelity of the measurement of lensed AGN host brightness through a set of extensive and realistic simulations of Hubble Space Telescope observation and lens modeling. First, the mock images of the lensed AGNs in our sample (see Ding et al. (2017), Table 1) were generated as realistically as possible. Second, the simulated AGN host galaxy images were reconstructed with the state-of-the-art lens modelling tool (glee111Developed by Suyu & Halkola (2010) based on Suyu et al. (2006) and Halkola et al. (2008).). Third, by fitting the host magnitude with the software galfit (Peng et al., 2002) and comparing the inference to the input value, Ding et al. (2017) found that the $L_{\rm host}$ can be recovered with better accuracy and precision than the uncertainty on single epoch $\mathcal{M}_{\rm BH}$ estimates ($\sim 0.5$ dex) for hosts as faint as 2$-$4 magnitudes dimmer than the AGN itself.

In this paper, we apply our advanced techniques to two strongly lensed systems (i.e. HE0435$-$1223 and RXJ1131$-$1231), with excellent imaging data. The host galaxy luminosity is inferred from the lens detailed model developed as part of the H0LiCOW collaboration222$H_{0}$ Lenses in COSMOGRAIL’s Wellspring, http://www.h0licow.org/. with the goal of measuring cosmological parameters from gravitational time delays (Suyu et al., 2017; Bonvin et al., 2017). $\mathcal{M}_{\rm BH}$ is inferred by applying a set of self-consistent calibrations of the virial method to the broad emission line properties measured by Sluse et al. (2012). In addition, we combine our new measurements with a large sample of AGNs taken from the literature and consistently recalibrated, and study the evolution of the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation for $146$ objects in the redshift range $0<z<4.5$. It is still unclear whether the bulge or the total luminosity provides the tightest correlation with $\mathcal{M}_{\rm BH}$ (Jahnke et al., 2009; Bennert et al., 2011b; Park et al., 2015). Thus, we consider both of them in this study.

This paper is organized as follows. We briefly describe the sample selection in Section 2. The host galaxy surface photometry and the $\mathcal{M}_{\rm BH}$ are inferred in Section 3 and 4, respectively. In Section 5, we present our main result. Discussion and conclusion are presented in Section 6 and 7. Throughout this paper, we adopt a standard concordance cosmology $H_{0}=70$ km s*-1* Mpc*-1*, $\Omega{{}_{m}}=0.30$, and $\Omega{{}_{\Lambda}}=0.70$. Magnitudes are given in the AB system.

2 Sample selection

First, we analyze the two quadruply-imaged AGN HE0435$-$1223 and RXJ1131$-$1231 (hereafter HE0435 and RXJ1131) with source redshifts at $1.693$ and $0.654$, respectively. Detailed information for these two systems is given by Suyu et al. (2017). Accurate lens models have been derived in an effort to measure cosmological parameters from gravitational time delays as described by Wong et al. (2017) and Suyu et al. (2013). These models provide the reconstructions of AGN hosts, from which in turn we estimate $L_{\rm host}$.

Second, we combine and compare our new measurements with those by Peng et al. (2006) (hereafter, P06). P06 explored the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ based on 20 non-lensed AGNs and 31 gravitationally lensed AGNs (including HE0435 and RXJ1131). P06 is so far the only paper in which the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation has been comprehensively investigated using lensed AGNs observed with HST. We note that for the two systems in common, the HST images used in our work are much deeper than those used by P06, and the lens models are much more detailed. Also, P06 was based on NIC2 images, as opposed to the much more powerful more modern cameras used in our work. Therefore, our measurements supercede those by P06 for these two systems. Furthermore, we exclude MG 2016+112 because it is a type II AGN (Koopmans et al., 2002) and the black hole mass using the virial method cannot be considered reliable. We also exclude the lens system B2045+265 used by P06 because of the incorrect redshift identification of the AGN spectrum by Fassnacht et al. (1999) (Nierenberg et al. 2017, in preparation).

Third, we combine our new measurements with samples of non-lensed AGN that have been measured by members of our team using the same techniques as those applied here. The samples consist of 52 intermediate redshift AGNs ($0.36<z<0.57$) summarized by Park et al. (2015) (hereafter P15), 27 distant AGNs ($0.5<z<1.9$) measured by Bennert et al. (2011b) and Schramm & Silverman (2013) (hereafter, B11 and SS13), and 19 local AGNs measurements (Bennert et al., 2010; Peterson et al., 2004). It is worth noticing that they are so far the largest HST imaging samples which are carefully selected as moderate-luminosity AGN, for which the contrast between nucleus and host galaxies is much more favorable for the inference of $L_{\rm host}$ than for high luminosity lensed quasars. Thus, their host luminosities are measured with high accuracy even without lensing.

Overall, our sample consists of two new lensed systems and active galaxies from the literature, including elliptical and spiral hosts with redshift up to 4.5. This total sample of 146 objects is the largest compilation of AGNs from HST which are cross-calibrated to study the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation. The objects and their basic properties are listed in Tab. 1 and 2.

3 Surface photometry

In this section, we describe the measurement of host luminosity. For HE0435 and RXJ1131, we first derived their host magnitude from the reconstructed surface brightness maps in the source plane. Then, we inferred the rest-frame R-band luminosities based on their spectral energy distribution (SED). For the other AGNs, we collected and homogenized their luminosities from the literature.

3.1 Surface photometry of HE0435 and RXJ1131

We used the software galfit to model the reconstructions from Wong et al. (2017) and Suyu et al. (2013). The reconstruction of HE0435 was fitted as the Sérsic profile with $n$ limited between $1-4$. It has been tested that this prior on $n$ does not bias the inference of magnitude (Ding et al., 2017). In the case of RXJ1131 a clearly visible residual image was present and the resulting parameters were physically acceptable when fitted with an additional profile, we concluded that the host galaxy is composed of a disk and a bulge. In this case, we fixed the reconstruction as two-component Sérsic profiles with $n$ equals to 1 and 4, corresponding to exponential disk profile and de Vaucouleurs (1948) profile, respectively. Although the luminosities of lens systems are corrected from lensing magnification using a lens model, small differences exist between models of different groups. The derived magnification rarely differs by more than $20\%$. According to detailed simulations presented by Ding et al. (2017), the inferred values of $L_{\rm host}$ can be recovered with sufficient accuracy and precision to study the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation using our approach. Finally, we derived the rest-frame R-band luminosity using a standard K-correction. These steps are described below in more detail for each system.

3.1.1 HE0435

HE0435 was imaged with HST/WFC3-IR through filter F160W from program HST-GO-12889 (PI: S. H. Suyu). Wong et al. (2017) produced a set of twelve reconstructions for this system, based on different assumptions, in order to estimate the amplitude of systematic errors associated with these choices. In 9/12 cases the source plane resolution was set to 40$\times$40 pixels. For the other three cases a higher resolution of 50$\times$50 pixels was adopted. The reconstructions were based on an image plane size of $\sim 1\aas@@fstack{\prime\prime}9$ square.

By fitting each of the twelve reconstructions with single Sérsic profile, we summarized the inference and found the mean value and the scatter of the host magnitudes are $m_{\rm host}=21.75\pm 0.13$; the inferred effective radius and Sérsic index are R_{\mathrm{eff}}$$=0.82\pm 0.14 arcsecond; $n=3.94\pm 0.14$, as shown in Tab. 3. Furthermore, to test the type of the host galaxy, we fitted the reconstructions as two-component Sérsic profile. However, we obtained unphysical results and no improvements in the fit indicating that the host galaxy of HE0435 is consistent with being a pure elliptical. One example of the reconstruction and its corresponding galfit best-fit are shown in Fig. 1, panel (a). We also note that there is a small structure at the lower left of the host. However, its brightness is negligible compared to the host which do not affect the inference of the $L_{\rm host}$. Interestingly, this could correspond tidal features in the host galaxies. If true, the mergers could be related to triggered AGN activity. It is beyond the scope of this work to pursue this further, but it would be intriguing to simulate the hosts with merger signature and to see if they can be recovered in the source reconstruction.

We can verify the accuracy of our result by carrying out simulations as described in our previous paper (Ding et al., 2017), using our inferred parameters as input. The observed and simulated HE0435 images are shown in Fig. 1, panel (b). By repeating the analysis on the simulated image, we recover the input value (input: $m_{\rm host}=21.75$ mag; output 21.88 mag) showing an accuracy much better than our target 1.25 mag (0.5 dex). We note that while in the simulations the PSF is assumed to be perfectly known, for the real data the PSF is inferred from the data using an iterative correction procedure (see Chen et al., 2016; Wong et al., 2017, Suyu et al. in preparation).

Following P06, we made no corrections for dust extinction of the host galaxy because they are likely to be small for a pure elliptical. The observed magnitudes were then transformed to rest-frame R-band by applying K-correction with Sbc template spectrum, using Coleman et al. (1980) templates. We used the Sbc template because the stellar populations cannot be older than a few Gyrs at this redshift and the local elliptical template would be too red. Nevertheless, since the HE0435 is observed through the F160W filter, which roughly corresponds to the rest-frame R-band at $z\sim 1.5$, the K-correction are only weakly dependent on the assumed SED (see Fig. 7 in P06), and do not contribute significantly to the error budget. Finally, the best inferred value of $L_{\rm host}$ of HE0435 in rest-frame R-band is \log$$L_{\rm host}$$=10.96 which is very close to the one inferred in P06 (i.e. \log$$L_{\rm host}$$=11.12).

Lens models based on archival HST Advanced Camera for Surveys (ACS) in the filter F555W and F814W are also available from Wong et al. (2017). Unfortunately, due to the short exposure time, the signal to noise ratio of the reconstructed host images in these bands is insufficient to infer the luminosity robustly in these bands and study the colors of the host. Thus they are not considered in this study.

3.1.2 RXJ1131

RXJ1131 is imaged with HST/ACS through filter F814W. A set of seven source resolutions including 50$\times$50, 52$\times$52, 54$\times$54, 56$\times$56, 58$\times$58, 60$\times$60, and 64$\times$64 pixels were selected when modelling the host image into source plane (Suyu et al., 2013), with a frame size of $\sim 2\aas@@fstack{\prime\prime}9$ square333Suyu et al. (2014) updated the model of RXJ1131. Given the similarity in the composite and power-law model by Suyu et al. (2014), a similar time delay distance is obtained (within $\sim 2\%$, and hence spatial scaling of the source due to mass-sheet degeneracy). This means the inference of total flux of the host should be unchanged to within $\sim 4\%$..

As noted by Suyu et al. (2013), all the reconstruction of the host show a compact peak near the center (see Fig. 2, panel (a), left panel), exhibiting the boundary line between the dominated area of bulge and disk which indicates the host galaxy is a spiral galaxy. Similarly, Claeskens et al. (2006) reconstructed the host of RXJ1131 and found it to be a spiral, disk dominated galaxy. Thus, we fitted the reconstructions as two-component Sérsic profiles, and the inferred properties of the disk are $m_{\rm disk}=20.07\pm 0.06$ mag; R_{\mathrm{eff}}$$\__{\rm disk}=0.84\pm 0.09 arcsecond and the properties of the bulge are $m_{\rm bulge}=21.81\pm 0.28$ mag; R_{\mathrm{eff}}$$\__{\rm bulge}=0.20\pm 0.08 arcsecond, as summarized in Tab. 3. An example of the reconstruction and the best-fit image are shown in Fig. 2, panel (a).

In the simulations of Ding et al. (2017), the host of RXJ1131 was assumed to be a single Sérsic profile with the magnitude between 19.0 and 20.5. In this work, we simulate a more realistic two-component profile, with key parameters (i.e. $m_{\rm host}$ and $R_{\mathrm{eff}}$) equal to the inferred values. The real and mock RXJ1131 image are shown in Fig. 2, panel (b). We first use a single Sérsic profile to fit the reconstruction, but applying this model is a poor representation with an obvious residual in the central image (i.e. Fig. 2, panel (c), left). This result suggests the lens model of RXJ1131 reconstructs the host with sufficiently high resolution to distinguish a bulge+disk model from a single component. Fitting with two-component Sérsic profile, we find that the residual map is much improved and both components can be reconstructed accurately with our data and analysis techniques: input $m_{\rm disk}=20.07$ mag and $m_{\rm bulge}=21.80$ mag; inferred values are $m_{\rm disk}=20.37$ mag and $m_{\rm bulge}=22.07$ mag.

As for HE0435, we derived the rest-frame R-band magnitude using a standard K-correction. At the redshift of RXJ1131, the conversion to R-band magnitude depends significantly on the adopted SED. Therefore, we determined the K-correction directly from the color of lensed host arc, based on the multi-band SED fitting available in the archive (GO-9744; PI: C. S. Kochanek). The final estimations are $\Delta{\rm mag}_{\rm disk}$(R$-$F814W)$\approx-0.3$ and $\Delta{\rm mag}_{\rm bulge}$(R$-$F814W)$\approx-0.7$. For detail, see Appendix A.

3.2 Surface photometry for the literature samples

In this section we describe our inference of the rest-frame R-band luminosity for the P06 and P15 samples.

P06 used the galfit (for non-lensed source) and lensfit444lensfit is a version of galfit that has been extended to fit lensed host galaxies while optimizing the mass model for the lens galaxy. For details, see P06. (for lensed source) softwares to infer the brightness of the AGN hosts, describing the host galaxy as single Sérsic profile. In P06, they reported a single value of luminosity for each object, suggesting that the host galaxies are ellipticals. However, in our analysis, we find the RXJ1131 is a spiral galaxy which suggests the approach in P06 may not be accurate for all the host galaxies. We will return to this issue in Section 6. Their measurements of absolute magnitude are presented by P06 (Tab 3 and 4 therein) in rest-frame R-band, Vega system. Thus, we transfer to AB system using $m_{\rm AB,R}-m_{\rm Vega,R}=0.21$ (Blanton & Roweis, 2007).

Similarly, for the P15 sample, which includes the samples from B11 and SS13, the host galaxy was fitted as a $n=4$ profile to model the bulge component; an exponential disk profile was added if deemed necessary. The rest-frame V-band luminosity is derived (see P15, Tab. 4, Column 3) by applying the K-correction with an early-type galaxy template spectrum. The same template is taken, we converted their results to rest-frame R-band. As the scatter in V-R colors is small, the associated uncertainty is estimated to be 0.16 mag (i.e. 0.06 dex in luminosity). Likewise, the luminosities for 19 local active galaxies are converted to rest-frame R-band.

Having obtained the R-band mag, the luminosity is derived by $\log L_{R}/L_{R,\odot}=0.4(M_{R,\odot}-M_{R})$, where $M_{R,\odot}=4.61$ (Blanton & Roweis, 2007). We summarized the homogenized R-band luminosities in Tab. 1 and 2.

4 Black hole mass

Assuming that the dynamics of the broad-line region (BLR) is dominated by the gravity of the central supermassive black hole, $\mathcal{M}_{\rm BH}$ can be derived by applying the so-called virial method, based on the size of the BLR ($R_{BLR}$) and the line-of-sight velocity width ($\Delta V$) which can be inferred in turn from continuum luminosity and emission line width, respectively. Usually, the CIV($\lambda 1549$), MgII($\lambda 2798$) and H$\beta$($\lambda 4861$) emission lines width and their local continuum luminosities \lambda\rm L_{\lambda}$$(1300\text{\AA}), \lambda\rm L_{\lambda}$$(3000\text{\AA}) and \lambda\rm L_{\lambda}$$(5100\text{\AA}) are used, respectively.

Sluse et al. (2012), P06, and P15 used different lines and different calibrations of the virial method. Thus we need to cross-calibrate them in order to avoid any systematic bias between the samples.

We first choose the recipe of P15 as the baseline:

[TABLE]

Then we align the self-consistent recipes (including emission lines using H$\beta$ and MgII) from McGill et al. (2008) with this baseline, by adding a small constant to the intercept (i.e. $-0.144$). In order to cross-calibrate the CIV-based estimator, we exploit the nine AGNs in our sample for which both MgII and CIV are available. We take the CIV recipe from P06 and add a small constant intercept (i.e. $-0.331$) to match on average the value inferred from MgII. Overall, we adopt the following virial formalism:

[TABLE]

with a{CIV, MgII, H$\beta$}={6.322, 6.623, 6.882}, b{CIV, MgII, H$\beta$}={0.53, 0.47, 0.518}, $\lambda_{line}${CIV, MgII, H$\beta$}={1350, 3000, 5100}. Having achieved a consistent cross-calibration, the $\mathcal{M}_{\rm BH}$ is estimated by adopting the emission line properties measured by Sluse et al. (2012), P06 and P15.

For the 19 local AGNs, rather than using continuum luminosity, $R_{BLR}$ was derived from time lags between continuum and emission-line variations (Peterson et al., 2004). Thus, same as P15, we adopt the reverberation-mapping $\mathcal{M}_{\rm BH}$ measurements with virial factor ($\log f=0.71$, Park et al., 2012; Bentz et al., 2009), noting that they are the anchor for the virial method and thus are inherently self-consistent.

$\mathcal{M}_{\rm BH}$ estimates are listed in Tab. 1 and 2, together with details on the emission line used. For RXJ1131, since the estimated $\mathcal{M}_{\rm BH}$ using MgII and H$\beta$ are very similar, we adopt their average. Moreover, we note that the values of $\mathcal{M}_{\rm BH}$ for HE0435 and RXJ1131 inferred in this paper are larger than the estimates by P06 (\Delta\log$$\mathcal{M}_{\rm BH} $=0.25$ and $0.44$ for HE0435 and RXJ1131, respectively), due to the fact that the properties of the emission lines of these two systems have been revised upwards by Sluse et al. (2012) based on data of superior quality.

5 Results

Following P15 and Ding et al. (2017), for the distant objects, we adopt total uncertainty for $L_{\rm host}$ and $\mathcal{M}_{\rm BH}$ of 0.2 dex ($\sim$0.5 mag) and 0.4 dex, respectively.

5.1 The observed $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation

The $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\mathcal{M}_{\rm BH}$-$L_{\rm total}$ relation defined by our samples are shown in Fig. 3, panel (a) and (b). There is a clearly positive correlation between $\mathcal{M}_{\rm BH}$ and $L_{\rm host}$ as in local samples. For a direct comparison to local samples, we fit the local $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ as:

[TABLE]

Using a Markov Chain Monte Carlo (MCMC) process we derive $\alpha=0.68\pm 0.18$; $\beta=0.74\pm 0.09$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\alpha=0.33\pm 0.22$; $\beta=0.95\pm 0.15$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$, with intrinsic scatter $\sigma_{\rm int}\sim 0.25$ for both of them. Consistent with previous work (e.g. P06, P15), the observed correlation at high redshift is nearly identical to the local. This is perhaps surprising, considering that both the black hole mass and host galaxy luminosity are expected to evolve over cosmic time. For example, in a minimal evolution toy model, the elliptical galaxies and their black hole are formed at high redshift and evolve passively thereafter. Thus, we expect $L_{\rm host}$ to fade over time, owing to aging stellar populations. To allow a direct comparison to the local samples, we considered this scenario in the next section.

5.2 The passive evolution-corrected $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation

In order to test the passive evolution scenario, we correct the observed $L_{\rm host}$ to its expected value at $z=0$ by accounting for the aging of the stellar populations. It has been shown that the evolution of the mass-to-light ratio of early-type galaxies can be effectively described as that of a single burst stellar population formed at appropriate redshifts (e.g., Treu et al., 2005). In order to represent the uncertainty in the star formation history we consider a range of single burst models formed at $z_{f}$ equals to 2, 3 and 5555 Stellar evolution is calculated with Galaxev (Bruzual & Charlot, 2003), based on Padova-1994 stellar evolutionary tracks, assuming Salpeter IMF, Solar metallicity, and no dust attenuation.. We choose to parametrize the evolution with the functional form $d{\rm mag}_{\rm R}=\delta_{m}d\log(1+z)$, i.e.

[TABLE]

with $\delta=-\delta_{m}/2.5$, so that

[TABLE]

For this parametrization, we derive that $\delta_{m}\simeq-3.7\pm 0.2$ (i.e. $\delta=1.48\pm 0.08$) provides a good representation of typical star formation histories.

This formalism is more accurate when considering a broad range in redshift with respect to adopting a single slope as a function of $d{\rm mag}/dz$ as done by P06 and P15. For a direct comparison, we also plot the passive evolution correction as a function of redshift in Fig. 4. Note that our chosen functional form describes well the P15 form at $z<1$ and the P06 form at $z\sim 3$ redshift, while avoiding the extreme corrections at very high-$z$ implied by previous parametrizations. Furthermore, our chosen functional form facilitate the analysis of the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ evolution in the following way. Combining Eq. 3 with the passive evolving correction, i.e. Eq. 5, and adding $\gamma^{\prime}$ term which describes the evolution of the correlation between $\mathcal{M}_{\rm BH}$ and observed $L_{\rm host}$, leads to the following formalism:

[TABLE]

In this equation, $\beta\delta$ represents the effects of passive evolution. The evolution at fixed present-day luminosity is given by $\gamma=\gamma^{\prime}+\beta\delta$. In this way the effects of the passive evolution correction can be easily separated and a different passive evolution model can be applied to the data, if desired. In our specific case, since we derived $\beta=0.74\pm 0.09$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ relation, the passive evolution term corresponds to approximately $\beta\delta=0.74\times 1.48\sim 1.0$, neglecting the effects of scatter and errors. Likewise, the passive evolution term is $\beta\delta=0.95\times 1.48\sim 1.4$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$ relation.

The resulting $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation after applying the passive evolution correction is shown in Fig. 3, panels (c) and (d). Clearly, after the correction, the high redshift samples are offset with respect to the local samples, indicating a tendency of BH in the more distant Universe to reside in less luminous hosts at fixed $\mathcal{M}_{\rm BH}$. This tendency is consistent with previous work, and also consistent with the studies of the $\mathcal{M}_{\rm BH}$-$\sigma_{*}$ (stellar velocity dispersion) and $\mathcal{M}_{\rm BH}$-$\mathcal{M}_{*}$ (stellar mass) correlations, which do not require correction for passive evolution (Treu et al., 2004; Woo et al., 2006, 2008; Bennert et al., 2011a).

We fit the offset in black hole mass at fixed passively evolved luminosity as a function of redshift in the form:

[TABLE]

where \Delta\log\mathcal{M}_{\rm BH}=\log\big{(}\frac{\mathcal{M}_{\rm BH}}{10^{7}M_{\odot}})-\alpha-\beta\log(\frac{L_{R,0}}{10^{10}L_{\odot}}), and obtain $\gamma=0.75\pm 0.11$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\gamma=0.95\pm 0.11$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$, as shown in Fig. 5, panels (a), (b). We also obtain $\gamma^{\prime}=-0.14\pm 0.11$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\gamma^{\prime}=-0.26\pm 0.12$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$, when not taking into account the passive evolving correction. As expected, the difference $\gamma-\gamma^{\prime}$ is consistent with the effects of the passive evolving correction, i.e. $\beta\delta\sim 1.0$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\beta\delta\sim 1.4$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$.

We conclude by noting that this fit does not take into account selection effects, which are discussed in the next section.

5.3 Selection effects

From Fig. 3, we can see that at high redshift we preferentially study systems with the larger $\mathcal{M}_{\rm BH}$ and $L_{\rm total}$. This is expected as observational samples tend to be flux limited and thus favor the high luminosity tail (and hence typically high $\mathcal{M}_{\rm BH}$) of the distribution. Like many other instances in astronomy, it is essential to take into account the selection function when estimating the evolution of the black hole mass host galaxy correlations (Treu et al., 2007; Lauer et al., 2007; Bennert et al., 2011a; Schulze & Wisotzki, 2014; Park et al., 2015).

Following P15, we take selection effects into account by using a Monte Carlo simulation method based upon the methodology introduced by Treu et al. (2007) and Bennert et al. (2010). The simulated samples are generated from a combination of the local active BH mass function from Schulze & Wisotzki (2010) and the local $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ relation from Bennert et al. (2010) with Gaussian random noise added as a function of the two free parameters $\gamma$ and intrinsic scatter of the correlation $\sigma_{\rm int}$. Note that the scatter is assumed to be independent of redshift in our description. For each object, the likelihood of the observed $\mathcal{M}_{\rm BH}$ with a given $L_{\rm host}$ is calculated from the simulated sample at the given $\gamma$ and $\sigma_{\rm int}$, and taking into account whether the object would be selected or not based on our sensitivity. Finally, by adopting uninformative uniform (flat) prior or lognormal prior from Bennert et al. (2010, $\sigma_{\rm int}=0.21\pm 0.08$), the posterior distribution function of $\gamma$ and $\sigma_{\rm int}$ is evaluated. Selection effects are modelled in the same way for the lensed-quasar sample, neglecting any second order effects related to lensing magnification. We note however, that these effects are small (Collett & Cunnington, 2016) and magification-related biases should affect the quasar and host galaxy in a similar manner, thus moving objects mostly along the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ correlation and not away from it.

Taking into account selection effects, the results of the inference are shown in Fig. 6. The fitted values of $\gamma$ are $0.6\pm 0.1$ ($\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$) and $0.8\pm 0.2$ ($\mathcal{M}_{\rm BH}$-$L_{\rm total}$), almost independent of the choice of prior. These values are consistent with the previous inference in Section 5.1 and 5.2.

Interestingly, the intrinsic scatter of the correlations is found to be consistent with typical values inferred for local samples ($0.3-0.4$ dex). This result is consistent with the hypothesis that well defined correlations exist at the redshifts probed by our sample, and indicates that we have not significantly underestimated our errors at high-$z$. It would be beneficial to study how the selection bias changes as a function of some key factors such as the values of $L_{\rm host}$ and $\mathcal{M}_{\rm BH}$, the level of the uncertainties and the redshift distribution of the samples. However, this topic is trivial in this study as we obtained consistent inference by either or not talking selection effects into account. Moreover, to study this relation quantitively requires considerate tests and simulations. Thus, we leave it in the future study.

6 Discussion

In this section, we first estimate the importance of potential systematic errors in § 6.1. Then, we carry out a detailed comparison with previous observational work in § 6.2. Finally we discuss how our measurements fit into our understanding of galaxies and BHs co-evolution in § 6.3.

6.1 Systematic errors

We have combined our new measurements with ones taken from the literature in order to increase the sample size and reduce statistical uncertainties. Even though we have restricted our analysis to the samples that have been analyzed in the most similar manner to our new data and we have cross-calibrated the black hole mass estimators, there are still some residual differences.

First, P06 obtained the luminosity of one galaxy by combined the fluxes together, even though some of them may include a disk component (e.g. RXJ1131). According to morphological studies of AGN host galaxies, the fraction of spiral/elliptical hosts of AGN is approximately one third (Kocevski et al., 2012), with the exact value depending on $\mathcal{M}_{\rm BH}$ and luminosity. Thus it is possible that P06 overestimates the bulge component of some of the host galaxies. The total luminosity should be less affected by this bias, even though not completely immune.

Furthermore, the subsample by SS13 included in the compilation by P15, was X-ray selected as opposed to optically selected like the rest of the non-lens sample (some of the lenses are radio-selected). This difference in selection could potentially lead to a systematic difference between the two samples.

In order to estimate these systematic uncertainties, we repeat the analysis by excluding the P06 and SS13 samples. This reduced sample will have significantly less statistical power, owing to the reduced size and redshift coverage, but should be more robust with respect to the systematic uncertainties discussed above. Given this reduced sample, we obtain $\gamma=0.88\pm 0.28$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\gamma=0.51\pm 0.28$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$, as shown in Fig. 5, panel (c), (d). Moreover, we use the same approach to study the selection effects and obtain the consistent inference, as illustrated in Fig. 7. Even though as expected the uncertainties are larger than for the full sample, the results are statistically mutually consistent at 1-$\sigma$ level. To facilitate the comparison between different $\gamma$, we summarize our inference in Tab. 4. We conclude that our inferred trends are not dominated by systematic differences between the samples, and systematic uncertainties of this kind are smaller than the random ones.

In this work, $\mathcal{M}_{\rm BH}$ estimates are derived using the CIV, MgII and H$\beta$ emission lines. However, the CIV and MgII lines are usually in outflow (Baskin & Laor, 2005; Richards et al., 2011; Denney, 2012) and therefore may not be dominated by the gravity of the central $\mathcal{M}_{\rm BH}$ and result in biased $\mathcal{M}_{\rm BH}$ estimates, especially for the CIV line (Trakhtenbrot & Netzer, 2012). Following McGill et al. (2008) the potential bias has been mitigated by cross-calibrating the $\mathcal{M}_{\rm BH}$ estimates based on the different lines. As a further sanity check, we fitted the $\gamma$ using only H$\beta$-based samples. We note that this H$\beta$ sample is very smilar to the subsample excluding P06 and SS13, and in fact the results are similar ($\gamma=1.10\pm 0.36$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\gamma=0.7\pm 0.37$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$). We conclude that any potential residual bias related to the use of lines other than H$\beta$ is smaller than statistical uncertainties or biases related to sample selection.

6.2 Comparison with previous work

P15, using a sample of 79 active galaxies, inferred the following evolutionary trends: $\gamma=0.9\pm 0.7$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\gamma=0.4\pm 0.5$ for the $\mathcal{M}_{\rm BH}$-$L_{\rm total}$. These are consistent with our inference, although their uncertainties are much larger, owing to the smaller sample size and reduced high redshift coverage. A similar result was obtained by P06, where they found that the ratio between $\mathcal{M}_{\rm BH}$ and $\mathcal{M}_{*}$ was $\sim 4$ times larger at $z\sim 2-4$ than today (i.e. $\gamma\sim 0.8-1.2$). The consistency between their measurements and ours are expected since the overall samples in this work are mostly composed of the samples by P15 and P06, even though there are some differences in the rest-frame bands chosen for photometry (we and P06 adopt rest-frame R, while P15 adopts rest-frame V), in the passive evolution correction, and in the black hole mass calibration.

The cosmic evolution of the $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ relation is a topic of intense debate in the literature. Many works have reported an evolutionary signal based on different relations including the $\mathcal{M}_{\rm BH}$-$L_{\rm host}$ (e.g. Treu et al., 2007; Bennert et al., 2010), the $\mathcal{M}_{\rm BH}$-$\mathcal{M}_{*}$ (e.g. McLure et al., 2006; Jahnke et al., 2009; Decarli et al., 2010; Cisternas et al., 2011; Bennert et al., 2011b; Trakhtenbrot et al., 2015) and the $\mathcal{M}_{\rm BH}$-$\sigma_{*}$ (e.g. Woo et al., 2006, 2008) correlations. Nevertheless, other observational studies (e.g. Shields et al., 2003; Greene & Ho, 2005; Komossa & Xu, 2007; Shen et al., 2008) found no evidence for evolution. In Shankar et al. (2016), they find serious biases in the $\mathcal{M}_{\rm BH}$-$\mathcal{M}_{*}$ relation and prove that $\sigma_{*}$ is more fundamental than any other variable. However, in Shankar et al. (2009), they show that there is no evolution in the $\mathcal{M}_{\rm BH}$-$\sigma_{*}$ relation once one accounts for the ages of local galaxies and the Sołtan argument. Moreover, Schulze & Wisotzki (2011, 2014) concluded that there is no statistically significant evidence for evolution once these selection effects are taken into account and corrected. Taking a different approach, DeGraf et al. (2015) used the results of the high-resolution numerical simulation MassiveBlackII to compare the observed and intrinsic evolution of the black hole mass host galaxy correlations and reproduced the evolutionary trend of the relation. Consistent with other considerations, they also found that the observed samples display steeper slopes than random ones, suggesting the selecting effects can exhibit faster evolution than a random sample. Similarly, by generating Monte Carlo realizations of the $\mathcal{M}_{\rm BH}$-$\sigma_{*}$ relation at $z=6$, Volonteri & Stark (2011) also found that due to selection bias the ’observable’ subsample would suggest an average positive evolution even when the intrinsic correlations is characterized by no or negative evolution at high redshift. These studies highlight once again the importance of taking selection effects into account.

Clearly, absence of evidence does not imply evidence of absence, and one way to make progress is to improve the precision and accuracy of the measurement. In our work, we attain much higher precision than previous work owing to the enlarged sample, including lensed quasars. Thanks to the large sample size, even when selection effects are taken into account, the evolutionary trend is detected at high significance ($\gamma\neq 0$ at more than 5-$\sigma$). However, using a reduced sample by excluding the subsamples from P06 and SS13, we obtain a smaller evolutionary trend, with larger uncertainties. These results are consistent at 1-$\sigma$ level (see Tab. 4), and highlight the importance of studying larger sample of high redshift lenses with the state-of-the-art data.

6.3 Implication for the co-evolution of black holes and their host galaxies

Our results are consistent with a scenario in which BHs in the distant Universe typically reside in lower stellar mass galaxies than today, assuming that the passively evolved luminosity tracks approximately stellar mass (see Bennert et al., 2011a; Schramm & Silverman, 2013, for a consistent direct measurement based on stellar mass determination). In order to end up on the local final relation, the stellar mass of the host galaxy would have to grow faster than $\mathcal{M}_{\rm BH}$.

An interesting clue to the physical mechanism driving the evolution could perhaps be found by comparing the inferred evolution for the correlation between $\mathcal{M}_{\rm BH}$ and the total host galaxy luminosity, and that with the bulge luminosity. We found those two to be comparable within the uncertainties. Previous work found the $\mathcal{M}_{\rm BH}$-$\mathcal{M}_{*}$bulge correlation to evolve somewhat faster than $\mathcal{M}_{\rm BH}$-$\mathcal{M}_{*}$total, albeit at low statistical significance, suggesting that one of the mechanisms at work is the build up of the bulge component from stars in the disk (Croton, 2006). It is difficult to perform a direct comparison, because of the fact that the P06 sample did not attempt bulge-disk decomposition, and we have also assumed a single passive evolution trend for the entire galaxy. Both effects could potentially suppress the differences between the evolution of the bulge and total luminosity with respect to the black hole mass.

Also, our sample extends to much larger redshift than that of B11. One possible explanation of this possible tension is that the dominant evolutionary mechanisms changes with redshift. At low redshift ($z\lesssim 1$), the growth of the bulge is dominated by the secular evolution with the redistribution of disk stars while at high redshift ($z\gtrsim 1$), the growth is dominated by major mergers (see Bennert et al., 2010, for a similar conjecture).

To settle this issue, it is crucial to obtain high-quality data and model large samples of lens systems, so that robust bulge to total luminosity decompositions can be carried out. It would also be beneficial to obtain multi-color data to estimate directly stellar mass, and ideally stellar kinematic information to distinguish pressure supported systems from rotationally supported ones.

Since the $L_{\rm host}$ of ellipticals would not change when considering the bulge and the total, we examine the offset using the sample limited to spiral galaxies. In our sample, there are 9 local spirals and 41 distant spirals, excluding SS13. Fitting the offset with this subsample, we obtain $\gamma=2.15\pm 0.41$ and $1.18\pm 0.41$ for $\mathcal{M}_{\rm BH}$-$L_{\rm bulge}$ and $\mathcal{M}_{\rm BH}$-$L_{\rm total}$, respectively, which are larger than the previous inference listed in Tab. 4 for the entire sample. This difference could suggest the spiral galaxies are undergoing a more rapid evolution than the ellipticals in order to end up on the local final relation. However, we caution that this result should be taken with a grain salt, given the small sample size of the local disk comparison sample.

7 summary

We presented a new measurement of the co-evolution of supermassive black holes and their host galaxies. First, we carried out a new analysis of two strongly lensed quasars, HE0435$-$1223 and RXJ1131$-$1231. By using the state-of-the-art lens models by Wong et al. (2017) and Suyu et al. (2013), we found that the host galaxies of HE0435 and RXJ1131 are well described by an elliptical and spiral surface brightness density profile, respectively. Then, we measured the host galaxy magnitude and tested for potential biases by carrying out realistic simulations following the procedure outlined by Ding et al. (2017). We found that the bias of our inference of $L_{\rm host}$ is small (0.1$-$0.2 mags) and that we can recover the host image precisely even if the host has multiple components (see Fig. 2, panel (c)). We estimated $\mathcal{M}_{\rm BH}$ by using a set of self-consistent single epoch estimators based on the quasar emission line properties as measured by Sluse et al. (2012).

Second, we combined our measurements with the published ones from the literature (Peng et al., 2006; Park et al., 2015), thus expanding our sample to $146$ active galaxies up to $z=4.5$. We have taken care of using self-consistent recipes to re-derive the black hole mass estimates and convert all the luminosities self-consistently to the rest-frame R-band.

Our main findings can be summarized as follows:

The observed correlations - without correction for evolution - are consistent with those observed in the local Universe. 2. 2.

The data are inconsistent with a passive evolution scenario. By correcting the host galaxy rest-frame luminosity to $z=0$, we find that galaxies are underluminous for a given $\mathcal{M}_{\rm BH}$, even neglecting growth by accretion. 3. 3.

The passively evolved correlations are well described by a relationship of the form \Delta\log$$\mathcal{M}_{\rm BH}$$=\gamma\log(1+z) with $\gamma=0.6\pm 0.1$ and $\gamma=0.8\pm 0.2$, respectively at fixed bulge and total host luminosity, taking into account selection effects.

Considering that stellar populations must fade as they get older, and considering that similar results have been found when studying the correlations between $\mathcal{M}_{\rm BH}$ and host galaxy velocity dispersion (Treu et al., 2004; Woo et al., 2006, 2008) and stellar mass (Jahnke et al., 2009; Bennert et al., 2011b; Schramm & Silverman, 2013), we are forced to conclude that the co-evolution of galaxies and black holes is non-trivial, in the sense that systems do not stay on the correlation as they evolve. At least for active galaxies in the range of black hole and stellar masses that can be analyzed with current technology, it appears that the growth of the black hole predates that of the bulge (Croton, 2006). However, given the complexity and variety of processes involved, direct comparisons with detailed numerical simulations are needed to further our understanding of the co-evolution of black holes and their hosts. Recent cosmological simulations including some prescriptions for black hole growth and feedback have been shown to reproduce the observations at least at $z<1$ (DeGraf et al., 2015). It will be interesting to carry out similar detailed comparisons, taking into account errors and observational selection functions, for a variety of models (e.g. Sijacki et al., 2015; Taylor & Kobayashi, 2016; Volonteri et al., 2016) and extending to higher redshifts. These comparisons will provide a powerful test of the various recipes that have been adopted to describe accretion and star formation physics at sub-grid level in numerical simulations.

Looking at the future, the sample of lensed quasars that can be analyzed with high fidelity is going to grow. Currently, ultra deep HST imaging data have been obtained for six additional strongly lensed systems666WFI2033$-$4723 and HE1104$-$1805 from H0LiCOW program HST-GO-12889 (PI: Suyu); SDSS1206$+$4332, HE0047$-$1756, SDSS0246$-$0825 and HS2209$+$1914 as part of Program HST-GO-14254 (PI: Treu). and their analysis will be described in a forthcoming paper. The sample of lensed quasars and their hosts that can be studied at high fidelity is likely to continue to grow as more such systems are discovered in wide field imaging and spectroscopic surveys (e.g. Agnello et al., 2015; More et al., 2016; Schechter et al., 2017; Ostrovski et al., 2017).

Acknowledgements

Based in part on observations made with the NASA/ESA Hubble Space Telescope, obtained at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. These observations are associated with programs # 9744, 12889, 14254. Financial support was provided by NASA through grants from the Space Telescope Science Institute.

We are grateful to Vivien Bonvin, Geoff C.-F. Chen, Frederic Courbin, Matthew A. Malkan, Cristian E. Rusu, Jong-Hak Woo, and Andreas Schulze for useful comments and suggestions that improved this manuscript. We thank Chien Peng for his help with the estimates of black hole mass. X.D. is supported by the China Scholarship Council. T.T. acknowledges support by the Packard Foundations through a Packard Research Fellowship and by the NSF through grants AST-1450141 and AST-1412315. S.H.S. gratefully acknowledges support from the Max Planck Society through the Max Planck Research Group. C.E.R. acknowledges support from the NSF grant AST-1312329. D.S. acknowledges funding support from a Back to Belgium grant from the Belgian Federal Science Policy (BELSPO). K.C.W. and D.P. is supported by an EACOA Fellowship awarded by the East Asia Core Observatories Association, which consists of the Academia Sinica Institute of Astronomy and Astrophysics, the National Astronomical Observatory of Japan, the National Astronomical Observatories of the Chinese Academy of Sciences, and the Korea Astronomy and Space Science Institute. V.N.B. gratefully acknowledges assistance from a National Science Foundation (NSF) Research at Undergraduate Institutions (RUI) grant AST-1312296. Note that findings and conclusions do not necessarily represent views of the NSF. V.B. acknowledge the support of the Swiss National Science Foundation (SNSF).

Appendix A The K-correction for the RXJ1131 host

We apply the K-correction to the observed magnitudes to obtain the rest-frame R-band magnitude. At the redshift of RXJ1131 ($z_{s}=0.654$), the conversion from F814W to R-band magnitude depends on the adopted SED templates, as shown in Fig. 8, panel(a). Therefore, we determine directly the K-correction of the disk component through SED fitting with the multi-band images (F555W, F814W, F160W), available in the archive (GO-9744; PI: C. S. Kochanek).

Briefly, we fit the SED for each pixel by using Fast (Kriek et al., 2009) based on Galaxev stellar evolution track (Bruzual & Charlot, 2003), assuming the solar metallicity, exponentially declining star formation history and Calzetti et al. (2000) dust extinction law, while the redshift is fixed to the spectroscopic one. The error for each pixel is calculated based on empty regions of images. We then derive the rest-frame R-band magnitude by using the best-fit template for each pixel, and see the offset from the observed F814W magnitude. As shown in Fig. 8, panel (b), $\Delta{\rm mag}_{\rm disk}$(R$-$F814W)$\approx-0.3$ is an appropriate estimate for the disk region. For the bulge, a direct measurement of SED is affected by the residual AGN contamination, and hence the arc corresponding to the bulge is half blue and half red. Thus, we assume that the age of stellar populations are $>3$ Gyr, where the correction hardly change (see Fig. 8, panel (a)), and adopt $\Delta{\rm mag}_{\rm bulge}$(R$-$F814W)$\approx-0.7$.