Subaru High-z Exploration of Low-Luminosity Quasars (SHELLQs) VIII. A less biased view of the early co-evolution of black holes and host galaxies

TL;DR

This study uses ALMA observations of low-luminosity quasars at z > 6 to explore early black hole and galaxy co-evolution, revealing they often align with local mass relations, unlike more luminous counterparts.

Contribution

It provides new insights into the co-evolution of black holes and host galaxies at high redshift by focusing on less biased, low-luminosity quasars.

Findings

Low-luminosity quasars have smaller [CII] luminosities than luminous quasars.

Most quasars show [CII]/FIR ratios similar to local star-forming galaxies.

Many low-luminosity quasars are on or below the local M_BH - M_bulge relation.

Abstract

We present ALMA [CII] line and far-infrared (FIR) continuum observations of three $z > 6$ low-luminosity quasars ($M_{\rm 1450} > -25$) discovered by our Subaru Hyper Suprime-Cam (HSC) survey. The [CII] line was detected in all three targets with luminosities of $(2.4 - 9.5) \times 10^8~L_\odot$, about one order of magnitude smaller than optically luminous ($M_{\rm 1450} \lesssim -25$) quasars. The FIR continuum luminosities range from $< 9 \times 10^{10}~L_\odot$ (3$\sigma$ limit) to $\sim 2 \times 10^{12}~L_\odot$, indicating a wide range in star formation rates in these galaxies. Most of the HSC quasars studied thus far show [CII]/FIR luminosity ratios similar to local star-forming galaxies. Using the [CII]-based dynamical mass ($M_{\rm dyn}$) as a surrogate for bulge stellar mass ($M_{\rm bulge}$), we find that a significant fraction of low-luminosity quasars are located on or even below the local $M_{\rm BH} - M_{\rm bulge}$ relation, particularly at the massive end of the galaxy mass distribution. In contrast, previous studies of optically luminous quasars have found that black holes are overmassive relative to the local relation. Given the low luminosities of our targets, we are exploring the nature of the early co-evolution of supermassive black holes and their hosts in a less biased way. Almost all of the quasars presented in this work are growing their black hole mass at much higher pace at $z \sim 6$ than the parallel growth model, in which supermassive black holes and their hosts grow simultaneously to match the local $M_{\rm BH} - M_{\rm bulge}$ relation at all redshifts. As the low-luminosity quasars appear to realize the local co-evolutionary relation even at $z \sim 6$, they should have experienced vigorous starbursts prior to the currently observed quasar phase to catch up with the relation.

Equations

$$\log \left(\frac{L_{[\mathrm{C}\text{emissiontype}\mathrm{II}]}}{L_{\odot }}\right)=(2.21\pm 1.13)+(0.55\pm 0.09)\times \log \left(\frac{\lambda L_{1450}}{L_{\odot }}\right).$$

$$\log \left(\frac{L_{[\mathrm{C}\text{emissiontype}\mathrm{II}]}}{L_{\mathrm{FIR}}}\right)=(3.05\pm 0.67)+(-0.50\pm 0.05)\log \left(\frac{L_{\mathrm{FIR}}}{L_{\odot }}\right),$$

$$\frac{{\mathrm{EW}}_{[\mathrm{C}\text{emissiontype}\mathrm{II}]}}{\mu \mathrm{m}}=0.527\times \frac{S_{[\mathrm{C}\text{emissiontype}\mathrm{II}]}[\mathrm{Jy}\mathrm{km}{\mathrm{s}}^{-1}]}{f_{\mathrm{cont}}[\mathrm{mJy}]},$$

$$\log \left(\frac{\mathrm{E}{\mathrm{W}}_{[\mathrm{C}\text{emissiontype}\mathrm{II}]}}{\mu \mathrm{m}}\right)=(1.76\pm 0.64)+(0.07\pm 0.02)\left(\frac{M_{1450}}{\mathrm{mag}}\right).$$

$$\log \left(\frac{L_{[\mathrm{C}\text{emissiontype}\mathrm{II}]}}{L_{\mathrm{FIR}}}\right)=(2.11\pm 1.00)+(-0.48\pm 0.09)\log \left(\frac{{\Sigma }_{\mathrm{FIR}}}{L_{\odot }\mathrm{kp}{\mathrm{c}}^2}\right),$$

$$M_{\mathrm{dyn}}/M_{\odot }=1.16\times 10^5{\left(\frac{v_{\mathrm{circ}}}{\mathrm{km}{\mathrm{s}}^{-1}}\right)}^2\left(\frac{D}{\mathrm{kpc}}\right)$$

$$\log \left(\frac{L_{\mathrm{FIR}}}{L_{\odot }}\right)=(2.94\pm 1.47)+(0.69\pm 0.11)\log \left(\frac{L_{\mathrm{Bol}}}{L_{\odot }}\right),$$

Full text

1 Introduction

The discovery of a tight correlation between the masses of central supermassive black holes ($M_{\rm BH}$) and those of galactic bulges ($M_{\rm bulge}$) or the stellar velocity dispersion in the local universe (e.g., Magorrian et al., 1998; Ferrarese & Merritt, 2000; Marconi & Hunt, 2003; Kormendy & Ho, 2013) strongly suggests that the formation and growth of the supermassive black holes (SMBHs) and their host galaxies are intimately linked, and the two undergo a co-evolution. Although the detailed mechanism by which the correlation arises unclear, theoretical models suggest that radiative and kinetic feedback of active galactic nuclei (AGNs) connected to the merger histories of galaxies play a pivotal role (e.g., Granato et al., 2004; Di Matteo et al., 2005; Hopkins et al., 2006; Li et al., 2007). A recent high resolution simulation work based on this scheme suggests that even a quasar at $z=7$ would follow the local co-evolution relation once we properly assess the mass of the host galaxy (Lupi et al., 2019). Detections of galaxy-scale massive AGN-driven outflows (e.g., Nesvadba et al., 2008; Maiolino et al., 2012; Cicone et al., 2014; Carniani et al., 2016; Toba et al., 2017), as well as the remarkable similarity of global star formation and SMBH accretion histories (Madau & Dickinson, 2014, for a review) would support this evolutionary scheme. As theoretical models usually make specific predictions for the time evolution of the systems, observations of the physical properties of both SMBHs and their host galaxies over cosmic time are essential to test and/or refine our current understanding of their build-up (Gallerani et al., 2017; Valiante et al., 2017).

From this perspective, high redshift quasars are a unique beacon of the early formation of SMBHs and their host galaxies, even in the first billion years of the universe (e.g., Mortlock et al., 2011; Bañados et al., 2018; Matsuoka et al., 2019). The last two decades have witnessed the discovery of $>200$ quasars at $z>5.7$ owing to wide-field optical and near-infrared (NIR) surveys, including the Sloan Digital Sky Survey (SDSS, e.g., Fan et al., 2003, 2006; Jiang et al., 2016), the Canada-France High-z Quasar Survey (CFHQS, e.g., Willott et al., 2007, 2010a), the Visible and Infrared Survey Telescope for Astronomy (VISTA) Kilo-degree Infrared Galaxy (VIKING, Venemans et al., 2013, 2015b), the United Kingdom Infrared Telescope (UKIRT) Infrared Deep Sky Survey (UKIDSS, Mortlock et al., 2009, 2011; Bañados et al., 2018), the Panoramic Survey Telescope & Rapid Response System (Pan-STARRS1, Bañados et al., 2014, 2016; Mazzucchelli et al., 2017), and several other projects (e.g., Kashikawa et al., 2015; Carnall et al., 2015; Kim et al., 2015; Tang et al., 2017; Reed et al., 2017). These surveys have found that luminous (absolute UV magnitude $M_{\rm 1450}\lesssim-25$ mag) quasars at $z\gtrsim 6$ are typically powered by SMBHs heavier than one billion solar masses and appear metal-enriched (e.g., Mortlock et al., 2011; De Rosa et al., 2014; Wu et al., 2015; Shen et al., 2019).

As host galaxies of $z>4$ quasars are hard to detect at rest-frame ultraviolet (UV) to optical wavelengths (Mechtley et al., 2012), longer wavelengths (i.e., far-infrared = FIR and sub/millimeter = sub/mm) cold gas and dust emission from star-forming regions have been used to probe such galaxies instead. These host galaxies possess copious amount of dust ($\sim 10^{8}~{}M_{\odot}$) and gas ($\sim 10^{10}~{}M_{\odot}$) with high FIR luminosities likely due to intense starburst (star formation rate $\gtrsim 100-1000$ $M_{\odot}$ yr*-1*) at $z>6$ (e.g., Bertoldi et al., 2003a, b; Petric et al., 2003; Priddey et al., 2003, 2008; Walter et al., 2004; Wang et al., 2007, 2008, 2011a, 2011b). These gaseous and dusty starburst regions appear to be spatially compact, with sizes of a few kpc or less (e.g., Wang et al., 2013; Venemans et al., 2017a), corresponding to the typical size of nearby bulges.

Among various emission lines of atoms and molecules in the cold interstellar medium (ISM), the fine structure line of singly ionized carbon, the 157.74 $\mu$m [C \emissiontypeII] ${}^{2}P_{3/2}\rightarrow{}^{2}P_{1/2}$ emission line (rest frequency $\nu_{\rm rest}=1900.5369$ GHz), of $z\gtrsim 6$ objects can be conveniently observed with ground-based sub/mm telescopes owing to an atmospheric window of $\sim 250$ GHz. The [C \emissiontypeII] line is the main coolant of the cold ISM, particularly of photodissociation regions (Hollenbach & Tielens, 1999), which makes this line an important tracer of star-forming activity. Sub/mm interferometers such as the Atacama Large Millimeter/submillimeter Array (ALMA) and the IRAM Plateau de Bure Interferometer (PdBI, now NOEMA) have sufficient resolution and sensitivity to resolve the gas dynamics of galaxies hosting not only optically luminous quasars (e.g., Wang et al., 2013, 2016; Bañados et al., 2015; Cicone et al., 2015; Venemans et al., 2016, 2017a, 2017b; Decarli et al., 2017, 2018; Mazzucchelli et al., 2017; Shao et al., 2017) but also low-luminosity ($M_{\rm 1450}\gtrsim-25$ mag) quasars (Willott et al., 2013, 2015, 2017; Izumi et al., 2018a).

Those high resolution studies have provided dynamical masses of the host galaxies. They have found that $z\gtrsim 6$ luminous quasars have ratios of SMBH mass to host galaxy mass $\sim 10$ times larger than the $z\sim 0$ relation, implying that these SMBHs formed significantly earlier than their hosts. However, there may be an observational bias, whereby more luminous quasars are powered by more massive SMBHs at high redshifts. This affects how closely these observations trace the underlying SMBH mass function at $z\gtrsim 6$, if there is a large scatter in $M_{\rm BH}$ for a given galaxy mass (e.g., Lauer et al., 2007; Schulze & Wisotzki, 2014). Indeed, early observations of low-luminosity CFHQS quasars showed that they are powered by less massive SMBHs ($\sim 10^{8}~{}M_{\odot}$), and show SMBH mass to bulge mass ratios roughly consistent with local galaxies (Willott et al., 2013, 2015, 2017). Thus, to achieve a comprehensive view on the early co-evolution of SMBHs and galaxies, observations of less-luminous (or smaller $M_{\rm BH}$) quasars are needed.

With this in mind, we have conducted ALMA pilot observations of several optically low-luminosity quasars discovered in an on-going deep multi-band ($g,r,i,z,y$), wide area imaging survey (Aihara et al., 2018a) with the Hyper Suprime-Cam (HSC, Miyazaki et al., 2018; Komiyama et al., 2018; Kawanomoto et al., 2018; Furusawa et al., 2018) mounted on the 8.2 m Subaru telescope. We have discovered more than 80 low-luminosity quasars at $z\gtrsim 6$ down to $M_{\rm 1450}\sim-22$ mag in this survey thus far (Matsuoka et al., 2016, 2018a, 2018b), including one $z>7$ object (Matsuoka et al., 2019). Most of these quasars constitute the break of the $z\sim 6$ quasar luminosity function (Matsuoka et al., 2018c), indicating that they represent the bulk of the quasar population at that high redshift. We then organized an intensive multi-wavelength follow-up consortium: Subaru High-z Exploration of Low-Luminosity Quasars (SHELLQs). In Izumi et al. (2018a), we presented Cycle 4 ALMA observations toward four $z\gtrsim 6$ HSC quasars in the [C \emissiontypeII] line and 1.2 mm continuum. The [C \emissiontypeII] and continuum luminosities of those HSC quasars are both comparable to local luminous infrared galaxy (LIRG)-class objects, suggesting that most of the quasar-host galaxies are less extreme starburst objects than has been found for luminous quasars. Like the CFHQS quasars, our HSC quasars tend to show SMBH mass to galaxy mass ratios similar to, or even lower than, the local co-evolution relation, which is also a clear contrast to their luminous counterparts.

In this paper, we report our ALMA Cycle 5 observations of [C \emissiontypeII] and underlying FIR continuum emission towards another three low-luminosity HSC quasars. We describe the observations in § 2. The basic observed properties of both the [C \emissiontypeII] line and the underlying FIR continuum emission are given in § 3. We then discuss [C \emissiontypeII]/FIR luminosity ratio, as an ISM diagnostic, and the less-biased early co-evolution of SMBHs and galaxies in § 4. Our findings are basically consistent with our previous Cycle 4 work (Izumi et al., 2018a): we summarize them in § 5. Throughout the paper, we adopt the cosmological parameters $H_{0}=70$ km s*-1* Mpc*-1*, $\Omega_{\rm M}=0.3$, and $\Omega_{\rm\Lambda}=0.7$.

2 Observations and data reduction

Three $z>6$ HSC quasars were observed during ALMA Cycle 5 (ID = 2017.1.00541.S, PI: T. Izumi) at band 6 between 2018 March 20 and 26. Our observations are summarized in Table 2, along with the basic target information. These observations were each conducted in a single pointing with $\sim 25\arcsec$ diameter field of view (FoV), which corresponds to $\sim 140$ kpc at the source redshifts (1$\arcsec$ $\sim$ 5.6 kpc). The phase tracking centers were set to the optical quasar locations (Matsuoka et al., 2018a). The absolute positional uncertainty is $\sim 0^{\prime\prime}.1$ according to the ALMA Knowledgebase111https://help.almascience.org/index.php?/Knowledgebase/List. With the minimum baseline length (15.1 m), the maximum recoverable scales of our observations are $\sim 9.5\arcsec$.

The receivers were tuned to cover the redshifted [C \emissiontypeII] line emission, whose frequencies were estimated from the measured redshifts of Ly$\alpha$. For the J2228$+$0152 observations, the total bandwidth was $\sim 7.5$ GHz, divided into four spectral windows of 1.875 GHz width. For the J1208$-$0200 and J2239$+$0207 observations, we set three spectral windows (i.e., 1.875 $\times$ 3 $\sim$ 5.6 GHz width in total) in one sideband, given the large uncertainties of their Ly$\alpha$ redshifts. The native spectral resolution was 7.813 MHz (8.7–8.9 km s*-1*), but 11–12 channels were binned to improve the signal-to-noise ratio ($S/N$), resulting in a final common velocity resolution of $\simeq 100$ km s*-1*.

Reduction and calibration of the data were performed with the Common Astronomy Software Applications package (CASA, McMullin et al., 2007) version 5.1.1, in the standard manner. All images were reconstructed with the CASA task clean with the Briggs weighting (gain = 0.1, robust = 0.5). The achieved synthesized beams and rms sensitivities at a velocity resolution of $\sim$100 km s*-1* are summarized in Table 2. All channels free of line emission were averaged to generate a continuum map for each source (see also Table 2), which was subtracted in the $uv$ plane before making the line cube. Throughout the paper, only statistical errors are displayed unless otherwise mentioned. The systematic uncertainty of the absolute flux calibration at ALMA band 6 is 10%, according to the ALMA Cycle 5 Proposer’s Guide.

{longtable}

*4c Description of our sample and the ALMA observations

J1208$-$0200 J2228$+$0152 J2239$+$0207

\endhead\endfootNote. Rest-frame UV properties are quoted from Matsuoka et al. (2018a) and Onoue et al. (2019).

The coordinates are updated after tying astrometric calibrations to the Gaia database.

*†*The Mg \emissiontypeII redshifts are 6.148 (J1208$-$0200) and 6.246 (J2239$+$0207), respectively (Onoue et al., 2019).

We did not obtain Mg \emissiontypeII these measurements in hand at the time of our ALMA observations.

\endlastfootRA (J2000.0) \timeform12h08m59s.22 \timeform22h28m47s.71 \timeform22h39m47s.48

Dec (J2000.0) $-$\timeform02D00’34”.9 $+$\timeform01D52’40”.4 $+$\timeform02D07’47”.4

$z_{\rm Ly\alpha}$ 6.2*†* 6.08 6.26*†*

$M_{\rm 1450}$ (mag) $-$24.3 $-$24.0 $-$24.6

Number of antennas 44–46 44–47 46–47

Baseline (m) 15.1–783.5 15.1–783.5 15.1–783.5

On-source time (minutes) 87 87 65

Bandpass calibrator J1229$+$0203 J2148$+$0657 J2148$+$0657

Complex gain calibrator J1218$-$0119 J2226+0052 J2226+0052

Flux calibrator J1229+0203 J2148+0657 J2148+0657

[C \emissiontypeII] cube

Beam size 0$\arcsec$.48 $\times$ 0$\arcsec$.39 0$\arcsec$.44 $\times$ 0$\arcsec$.40 0$\arcsec$.45 $\times$ 0$\arcsec$.38

Position Angle (East of North) $-$59$\arcdeg$.6 $-$74$\arcdeg$.2 $-$83$\arcdeg$.2

rms noise per 100 km s*-1* 0.10 0.10 0.11

(mJy beam*-1*)

rms noise per 100 km s*-1* 0.15 0.15 0.17

(mJy; 1$\arcsec$.0 aperture)

Continuum map

Observed continuum frequency (GHz) 266.0 260.1 261.4

Beam size 0$\arcsec$.48 $\times$ 0$\arcsec$.38 0$\arcsec$.45 $\times$ 0$\arcsec$.41 0$\arcsec$.45 $\times$ 0$\arcsec$.38

Position Angle (East of North) $-$60$\arcdeg$.7 $-$81$\arcdeg$.1 $-$82$\arcdeg$.4

rms noise 16.3 11.2 19.1

($\mu$Jy beam*-1*)

rms noise 19.5 15.7 26.0

($\mu$Jy; 1$\arcsec$.0 aperture)

3 Results

Figure 1 displays the spatial distribution of the velocity-integrated [C \emissiontypeII] line emission (0th moment) as well as the underlying rest-frame FIR continuum emission ($\lambda_{\rm rest}$ $\simeq$ 158 $\mu$m) of the three objects observed in Cycle 5. Those moment 0 maps were made with the CASA task immoments, integrating over the full velocity range containing the line emission 222The [C \emissiontypeII] emission line of J2239+0207 spans over two spectral windows (see Figure 2). We thus integrated the emission between 261.50 – 262.25 GHz in one window and 262.25 – 262.75 GHz in the other, and combined them to generate the moment 0 map.. The [C \emissiontypeII] emission was detected in all sources, with no apparent spatial offset of their peak locations from the optical centroids given the astrometric uncertainties. However, FIR continuum emission was only significantly detected in J2239+0207 at the original angular resolution. Given that the [C \emissiontypeII] emission seems to be somewhat extended relative to the synthesized beams, we measured FIR properties of these quasar host galaxies within a common 1\arcsec.0 aperture. With this treatment, FIR continuum emission was marginally detected in J1208$-$0200 as well. The rms sensitivities for these 1\arcsec.0 aperture measurements are listed in Table 2. The relevant FIR properties of the targets are shown in Table 3.1, which will be explained in detail in the following.

3.1 [C \emissiontypeII] line properties

Figure 2 displays the [C \emissiontypeII] line spectra measured with the 1\arcsec.0 aperture. Given the modest signal-to-noise ratio (S/N) achieved, we simply fitted each spectrum with a single Gaussian profile, which delivered the velocity centroid (or redshift = $z_{\rm[C\,\emissiontype{II}]}$), line width (full width at half maximum = FWHM${}_{\rm[C\,\emissiontype{II}]}$), and the velocity-integrated line flux ($S_{\rm[C\,\emissiontype{II}]}$) of the quasar host galaxy (Table 3.1).

The FWHM${}_{\rm[C\,\emissiontype{II}]}$ of these three HSC quasar host galaxies presented here, as well as those measured in our Cycle 4 work (Izumi et al., 2018a) 333The mean and standard deviation of the distribution is 327 $\pm$ 135 km s*-1* for the full sample of seven objects (Cycle 4 + 5)., are consistent with the distribution constructed from a large sample of $z\gtrsim 6$ quasar host galaxies (Decarli et al., 2018), but J1208$-$0200 and J2239+0207 lie at the lower and higher extremes of the distribution, respectively (Figure 3). Thus, there seems to be no clear correlation of FWHM${}_{\rm[C\,\emissiontype{II}]}$ and quasar luminosity: indeed, the Spearman rank correlation coefficient for the relation in Figure 3 is only $-$0.21 with a null hypothesis probability of 0.16. The line profile of J2239+0207 is clearly flat at the peak, in the velocity range of $-$300 to +300 km s*-1*, although we fit the profile with a single Gaussian for simplicity. A flat line profile was also found in a $z=4.6$ quasar (Kimball et al., 2015) and a $z=6.13$ quasar (Shao et al., 2017). Such a profile suggests that the [C \emissiontypeII] line emission originates from a rotating disk. On the other hand, we cannot discuss the detailed dynamical nature (e.g., rotation-dominant or dispersion-dominant) of J1208$-$0200 and J2228$+$0152 as these are barely resolved.

The velocity centroid ($z_{\rm[C\,\emissiontype{II}]}$) of J2228+0152 agrees well with $z_{\rm Ly\alpha}$ (Table 2), but $z_{\rm[C\,\emissiontype{II}]}$ is offset significantly blueward from $z_{\rm Ly\alpha}$ for J1208$-$0200 and J2239+0207. These offsets could simply be the consequence of the considerable uncertainties in $z_{\rm Ly\alpha}$ due to severe intergalactic absorption ($\Delta z\sim 1000$ km s*-1*, e.g., Eilers et al., 2017). Among the ionized lines predominantly emerging from the broad line region of quasars, C \emissiontypeIV $\lambda$1549 is usually found to be blueshifted with respect to the host galaxies (Shen et al., 2016; Mazzucchelli et al., 2017). In the case of our 7 HSC quasars with [C \emissiontypeII] measurements, the measured blueshifts of C \emissiontypeIV are $\sim 400-600$ km s*-1* with respect to $z_{\rm[C\,\emissiontype{II}]}$ (Onoue et al., 2019). The blueshifted nature of C \emissiontypeIV indicates the existence of outflowing gas close to the central SMBHs.

Regarding Mg \emissiontypeII $\lambda$2798, however, we do not see significant blueshifts ($\gtrsim 500$ km s*-1*) in these HSC quasars (Onoue et al., 2019) 444J1208$-$0200 even shows a large redshift in Mg \emissiontypeII with respect to [C \emissiontypeII] by $\sim 1260$ km s*-1* (out of the displayed range in Figure 2), although its Mg \emissiontypeII profile would be affected by OH sky emission (Onoue et al., 2019)., while some other quasars show such blueshifts (e.g., Venemans et al., 2016; Wang et al., 2016; Decarli et al., 2018), as shown in the compilation of the velocity shift ($\Delta v=v_{\rm Mg\,\emissiontype{II}}-v_{\rm[C\,\emissiontype{II}]}$, Figure 4). The mean and standard deviation of the distribution is $-$284 $\pm$ 607 km s*-1*. Although nuclear outflows can be a possible origin of such shifts, only a marginal and no clear correlations are found for $\Delta v$, between quasar nuclear luminosity or Eddington ratio 555We used the bolometric correction factor of 4.4 from 1450 Å luminosity (Richards et al., 2006) and Mg \emissiontypeII-based $M_{\rm BH}$ data to compute the Eddington ratios., with the corresponding Spearman rank coefficient of 0.35 ($p$-value = 0.05) and 0.04 ($p$-value = 0.82), respectively, for the samples shown in Figure 4. Hence the cause of the Mg \emissiontypeII blueshifts remains unclear. It is nevertheless advisable to use the Mg \emissiontypeII line with caution as a good indicator of the systemic redshift, given the wide distribution of $\Delta v$. This is particularly an issue for ALMA observations, as the width of its single baseband is only $\sim 2250$ km s*-1* at 250 GHz.

The [C \emissiontypeII] line luminosities of our sources were calculated with the standard equation, $L_{\rm[C\,\emissiontype{II}]}=1.04\times 10^{-3}~{}S_{\rm[C\,\emissiontype{II}]}~{}\nu_{\rm rest}~{}(1+z_{\rm[C\,\emissiontype{II}]})^{-1}~{}D^{2}_{L}$ (Solomon & Vanden Bout, 2005), where $L_{\rm[C\,\emissiontype{II}]}$ is the [C \emissiontypeII] line luminosity in units of $L_{\odot}$ and $D_{L}$ is the luminosity distance in units of Mpc. The resultant luminosities lie in the range $L_{\rm[C\,\emissiontype{II}]}=(2.4-9.5)\times 10^{8}~{}L_{\odot}$, consistent with our Cycle 4 measurements for the other four HSC quasars ($L_{\rm[C\,\emissiontype{II}]}\simeq(4-10)\times 10^{8}~{}L_{\odot}$, Izumi et al., 2018a). These values are also comparable to other low-luminosity quasars at $z\gtrsim 6$ (Willott et al., 2013, 2015, 2017). In contrast, optically-luminous quasars at $z\gtrsim 6$ show higher $L_{\rm[C\,\emissiontype{II}]}$, typically in the range $\simeq(1-10)\times 10^{9}~{}L_{\odot}$ (e.g., Wang et al., 2013, 2016; Venemans et al., 2016; Decarli et al., 2018). Indeed there is a positive correlation between the quasar nuclear luminosity (represented by the 1450 Å monochromatic luminosity) and $L_{\rm[C\,\emissiontype{II}]}$ as shown in Figure 5. The derived Spearman rank coefficient is 0.59, with the associated $p$-value of 2.0 $\times$ 10*-5*. Our linear regression analysis returns the best fit line for the correlation as

[TABLE]

As the FIR emission is likely to be dominated by cold dust emission powered by star formation (e.g., Leipski et al., 2014), one plausible and simple explanation for the above correlation is that both the black hole accretion and the star formation is driven by a common reservoir of gas.

The SFR of these quasars can then be estimated by attributing the [C \emissiontypeII] emission heating solely to young stars: ${\rm SFR}_{\rm[C\,\emissiontype{II}]}/M_{\odot}~{}{\rm yr^{-1}}=1.0\times 10^{-7}(L_{\rm[C\,\emissiontype{II}]}/L_{\odot})^{0.98}$ (De Looze et al., 2011). This relation has an intrinsic scatter of 0.3 dex and is based on the Kroupa initial mass function (IMF, Kroupa, 2001). The derived values range from 16 to 63 $M_{\odot}$ yr*-1* (Table 3.1), well within the SFR-range of local LIRG-class systems (e.g., Díaz-Santos et al., 2013; Sargsyan et al., 2014). Note that the De Looze et al. (2011) relation was derived for objects with $L_{\rm FIR}\lesssim 10^{12}~{}L_{\odot}$, thus it may not be appropriate to apply this for J2239+0207 ($L_{\rm FIR}\simeq 2\times 10^{12}~{}L_{\odot}$).

{longtable}

*4c Rest-frame FIR properties of the HSC quasars J1208$-$0200 J2228+0152 J2239+0207

\endhead\endfootNote. These were measured with a common 1$\arcsec$.0 aperture.

The FIR luminosities were estimated with a graybody spectrum model.

The upper limits are the 3$\sigma$ values.

${\rm SFR}_{\rm[C\,\emissiontype{II}]}$/$M_{\odot}$ yr*-1* = 1.0 $\times$ 10*-7* ($L_{\rm[C\,\emissiontype{II}]}$/$L_{\odot}$)0.98 (De Looze et al., 2011).

${\rm SFR}_{\rm TIR}/M_{\odot}~{}{\rm yr^{-1}}=1.49\times 10^{-10}L_{\rm TIR}/L_{\odot}$ (Murphy et al., 2011).

\endlastfoot$z_{\rm[C\,\emissiontype{II}]}$ 6.1165 $\pm$ 0.0002 6.0805 $\pm$ 0.0004 6.2497 $\pm$ 0.0004

FWHM${}_{\rm[C\,\emissiontype{II}]}$ (km s*-1*) 184 $\pm$ 28 266 $\pm$ 48 607 $\pm$ 44

$S_{\rm[C\,\emissiontype{II}]}$ (Jy km s*-1*) 0.280 $\pm$ 0.056 0.253 $\pm$ 0.059 0.955 $\pm$ 0.085

$L_{\rm{[C\,\emissiontype{II}]}}$ (108 $L_{\odot}$) 2.71 $\pm$ 0.54 2.43 $\pm$ 0.57 9.53 $\pm$ 0.85

$f_{\rm 1.2mm}$ ($\mu$Jy) 85 $\pm$ 20 $<$47 1110 $\pm$ 26

EW${}_{\rm{[C\,\emissiontype{II}]}}$ ($\mu$m) 1.73 $\pm$ 0.52 $>$2.82 0.45 $\pm$ 0.04

SFR${}_{\rm[C\,\emissiontype{II}]}$ ($M_{\odot}$ yr*-1*) 18 $\pm$ 4 16 $\pm$ 5 63 $\pm$ 8

Assumption: $T_{\rm d}$ = 47 K, $\beta$ = 1.6, $\kappa_{\rm 250}$ = 0.4 cm2 g*-1*

$L_{\rm FIR}$ (1011 $L_{\odot}$) 1.62 $\pm$ 0.37 $<$0.94 21.74 $\pm$ 0.51

$L_{\rm TIR}$ (1011 $L_{\odot}$) 2.29 $\pm$ 0.52 $<$1.34 30.65 $\pm$ 0.72

SFRTIR ($M_{\odot}$ yr*-1*) 34 $\pm$ 8 $<$20 453 $\pm$ 10

$M_{\rm dust}$ ($10^{7}~{}M_{\odot}$) 1.2 $\pm$ 0.3 $<$0.7 15 $\pm$ 1

$L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ($10^{-3}$) 1.67 $\pm$ 0.51 $>$2.58 0.44 $\pm$ 0.04

We used the CASA task imfit to fit a two-dimensional Gaussian profile to the [C \emissiontypeII] integrated intensity (0th moment) maps, and estimated their beam-deconvolved spatial extents. The maps made with the original resolutions (Figure 1) were used for this purpose. This image-plane fitting method has been widely used in previous submm studies of $z\gtrsim 6$ quasar host galaxies (e.g., Willott et al., 2015; Venemans et al., 2016), which enables a direct comparison with these earlier studies. The estimated values are listed in Table 3.1: their FWHM sizes are $\sim 2.1-4.0$ kpc (major axis). Although the associated uncertainties are admittedly large, these sizes are comparable to those found in our previous work on the other four HSC quasars, as well as to many optically luminous quasars (Izumi et al., 2018a).

{longtable}

*3c Spatial extent of the star-forming region Name Size ([C \emissiontypeII] FWHM) Size (FIR continuum FWHM)

\endhead\endfootNote. The original resolution data was used for the measurements.

\endlastfootJ1208$-$0200 (0\arcsec.63 $\pm$ 0\arcsec.06) $\times$ (0\arcsec.35 $\pm$ 0\arcsec.07) -

(3.6 $\pm$ 0.3) kpc $\times$ (1.9 $\pm$ 0.4) kpc

J2228+0152 (0\arcsec.38 $\pm$ 0\arcsec.05) $\times$ (0\arcsec.18 $\pm$ 0\arcsec.07) -

(2.1 $\pm$ 0.3) kpc $\times$ (1.0 $\pm$ 0.4) kpc

J2239+0207 (0\arcsec.72 $\pm$ 0\arcsec.02) $\times$ (0\arcsec.46 $\pm$ 0\arcsec.02) (0\arcsec.22 $\pm$ 0\arcsec.03) $\times$ (0\arcsec.11 $\pm$ 0\arcsec.05)

(4.0 $\pm$ 0.1) kpc $\times$ (2.6 $\pm$ 0.1) kpc (1.2 $\pm$ 0.2) kpc $\times$ (0.7 $\pm$ 0.5) kpc

3.2 FIR continuum properties

The observed 1.2 mm continuum flux densities ($f_{\rm 1.2mm}$) are used to determine their FIR luminosities ($L_{\rm FIR}$) integrated over the rest-frame wavelengths of $\lambda_{\rm rest}=42.5-122.5$ $\mu$m (Helou et al., 1988). Here we assume a graybody spectrum with dust temperature of $T_{\rm d}=47$ K and emissivity index666Emissivity $\propto\nu^{\beta}$ of $\beta=1.6$ based on the mean spectral energy distribution of high redshift optically and FIR luminous quasars (Beelen et al., 2006; Leipski et al., 2014), to be consistent with previous work (e.g., Wang et al., 2013; Willott et al., 2013; Venemans et al., 2016). However, these values are likely to vary significantly from source to source (Venemans et al., 2018; Liang et al., 2019). If our HSC quasars instead have $T_{\rm d}$ close to the value found for nearby LIRG-class systems ($\sim 35$ K, U et al., 2012), the resultant inferred $L_{\rm FIR}$ values would be $\sim 3\times$ lower. We hereafter only consider the uncertainties of flux measurements, not that of the $T_{\rm d}$, which should be constrained further with future multi-wavelength observations. Note that the influence of the cosmic microwave background (CMB) radiation on the submm observations at high redshifts (da Cunha et al., 2013) is not considered, as that effect is negligible as long as we adopt $T_{\rm d}=47$ K.

The $f_{\rm 1.2mm}$ and $L_{\rm FIR}$ measured with the common 1\arcsec.0 aperture are listed in Table 3.1. J1208$-$0200 was marginally detected ($\sim 4\sigma$) with $L_{\rm FIR}=(1.6\pm 0.4)\times 10^{11}~{}L_{\odot}$, which is slightly smaller than those of the four Cycle 4 HSC quasars ($L_{\rm FIR}\sim(3-5)\times 10^{11}~{}L_{\odot}$, Izumi et al., 2018a). J2228+0152 is undetected, with a $3\sigma$ upper limit of $f_{\rm 1.2mm}<47$ $\mu$Jy and $L_{\rm FIR}<9\times 10^{10}~{}L_{\odot}$ (i.e., below the luminosity range of a LIRG), making it one of the lowest $L_{\rm FIR}$ quasar host galaxies ever studied at $z\gtrsim 6$. The $L_{\rm FIR}$ of these HSC quasars are thus smaller by factors of $\sim 10-100$ than most of the $z\gtrsim 6$ optically-luminous quasars studied thus far (e.g., Wang et al., 2007, 2008). On the other hand, for J2239+0207 we found $L_{\rm FIR}\simeq 2\times 10^{12}~{}L_{\odot}$, showing that there is a broad distribution in $L_{\rm FIR}$ even among HSC quasars of comparable UV/optical luminosities. The relationship between quasar luminosity and $L_{\rm FIR}$ is further discussed in § 4.2.2.

We measured the size of the FIR continuum-emitting region of J2239+0207 with the imfit task, finding a significantly smaller size than that of the [C \emissiontypeII] emitting region (Table 3.1; see also Figure 1). The [C \emissiontypeII] emitting region is often larger than the continuum-emitting region in high-$z$ quasars (e.g., Wang et al., 2013; Venemans et al., 2016), although the cause for this remains unclear.

The total infrared luminosity ($L_{\rm TIR}$) integrated over the 8–1000 $\mu$m range is supposed to be powered by star formation, and thus gives an independent estimate of SFR. We use the conversion, ${\rm SFR_{\rm TIR}}/M_{\odot}~{}{\rm yr}^{-1}=1.49\times 10^{-10}~{}L_{\rm TIR}/L_{\odot}$ (Murphy et al., 2011), which is based on the Kroupa IMF (Kroupa, 2001). The SFRs based on this method are also listed in Table 3.1. We also derived dust mass $M_{\rm dust}$ from $L_{\rm FIR}$ with $M_{\rm dust}=L_{\rm FIR}/(4\pi\int\kappa_{\nu}B_{\nu}d\nu)$, where $\kappa_{\nu}$ is the mass absorption coefficient, taken to be $\kappa_{\nu}=\kappa_{0}(\nu/250~{}{\rm GHz})^{\beta}$ with $\kappa_{0}=0.4$ cm2 g*-1* (Alton et al., 2004), and $B_{\nu}$ is the Planck function. The derived values span a wide range from $<7\times 10^{6}~{}M_{\odot}$ to $1.5\times 10^{8}~{}M_{\odot}$ (Table 3.1).

3.3 Further details of J2228+0152 and J2239+0207: interactions and companions?

The [C \emissiontypeII] emission of J2228+0152 seems to be extended to the east, although the statistical significance of the extended component is only modest ($\sim 3-3.5\sigma$). However, there is also weak ($\sim 2-3\sigma$) FIR continuum emission around the eastern extension, which motivates us to further investigate its structure. To this end, we constructed [C \emissiontypeII] velocity channel maps of J2228+0152 using the MIRIAD software (Sault et al., 1995). We found that there is one [C \emissiontypeII] emitter candidate at the location of the eastern extension (Figure 6), which is detected at $5.5\sigma$. If this is a [C \emissiontypeII] emitter, it is located at $z_{\rm[C\,\emissiontype{II}]}=6.068$ with $S_{\rm[C\,\emissiontype{II}]}=0.055$ Jy km s*-1* or $L_{\rm[C\,\emissiontype{II}]}=5.3\times 10^{7}~{}L_{\odot}$. The velocity offset and the projected separation on the sky of this emitter candidate, measured from the central quasar, are $\sim-540$ km s*-1* and $\sim 7$ kpc, respectively. Thus, the eastern extension may be related to the interaction of J2228+0152 and this companion emitter. Although we need higher sensitivity observations to further study the nature of this emitter as it was only detected in one channel, companion galaxies and extended (or interacting) morphologies have been identified in [C \emissiontypeII] emission around some quasars at $z\gtrsim 6$ (Decarli et al., 2017).

Another interesting object is J2239+0207, as it shows the highest $L_{\rm FIR}$ of the seven HSC quasars thus far studied with ALMA. One possible origin of the high FIR luminosity is a merger of two or more galaxies, as is often observed in nearby ULIRGs (e.g., Sanders & Mirabel, 1996). Thus, we searched for features indicative of galaxy mergers or interactions in [C \emissiontypeII] velocity channel maps of J2239+0207 (Figure 7). As we remarked in § 2, the [C \emissiontypeII] line spans two spectral windows, allowing us to test the robustness of the features we see. Note, however, that the velocity spacings of these two data sets are not perfectly matched.

The [C \emissiontypeII] emission of J2239+0207 is slightly extended relative to the synthesized beam and shows filamentary structures in some channels. Although these components are of low statistical significance ($\lesssim 3-5\sigma$) and do not seem to be well matched in the two windows, we found a [C \emissiontypeII] emission candidate at exactly the same location in the two datasets, at similar velocities (+87 km s*-1* in one spectral window and +66 km s*-1* in the other). If this is indeed a [C \emissiontypeII] emitter, an interaction with J2239+0207 (projected separation $\sim 1\arcsec$ $\sim$ 5.6 kpc) may have triggered the observed starburst activity. Meanwhile, if it is real, this [C \emissiontypeII] emitter ($z_{\rm[C\,\emissiontype{II}]}=6.248$) has $S_{\rm[C\,\emissiontype{II}]}\sim 0.063$ Jy km s*-1* or $L_{\rm[C\,\emissiontype{II}]}\sim 6.2\times 10^{7}~{}L_{\odot}$ (corresponding to ${\rm SFR}_{\rm[C\,\emissiontype{II}]}\sim 4~{}M_{\odot}$ yr*-1*), which are rather modest values. Given the compact size and narrow velocity width, a dynamical mass of this emitter should be small. Thus a possible future merger of this emitter and J2239+0207 would be rather minor. J2239+0207 itself must thus have copious amount of cold material to support the observe starburst activity. However, our limited sensitivity (the significance of the emission is $\sim 5\sigma$ and $\sim 3\sigma$ in each window, respectively; FIR continuum emission is below 3$\sigma$ at this position) means that we cannot call this a robust detection. Again, we need much deeper observations to reveal its nature.

3.4 Continuum emitters

We also searched for companion continuum emitters in these three HSC quasar fields. Here we conservatively define an emitter as one that shows $\geq 5\sigma$ significance. In the FoV of J1208$-$0200 ($5\sigma=82$ $\mu$Jy beam*-1*) and J2239+0207 ($5\sigma=96$ $\mu$Jy beam*-1*), no such emitter was found. In contrast, one emitter was found slightly outside the nominal FoV of J2228+0152 ($5\sigma=56$ $\mu$Jy beam*-1*), located 13\arcsec.7 away from the quasar position (Figure 8). This emitter is bright with a peak flux density of 0.25 mJy beam*-1*, and is clearly more extended than the synthesized beam. The deconvolved size by simply fitting a two-dimensional Gaussian with the imfit task is (0\arcsec.417 $\pm$ 0\arcsec.078) $\times$ (0\arcsec.306 $\pm$ 0\arcsec.080). We found that this object is also detected in the HSC optical bands as $g=25.12\pm 0.21$ mag, $r=24.57\pm 0.18$ mag, $i=25.04\pm 0.37$ mag, $z=23.91\pm 0.25$ mag, and $y=24.32\pm 0.82$, respectively. The HSC photometric redshift catalog (Tanaka et al., 2018) from the first data release (Aihara et al., 2018b), along with the Bayesian-based Mizuki code (Tanaka, 2015), suggests that the redshift of this source is $z_{\rm photo}=2.26\pm 0.56$. If we rely on this $z_{\rm photo}$, the above source size is equivalent to $\sim 3.4$ kpc $\times$ 2.5 kpc. Such faint (e.g., a few 100 $\mu$Jy at $\sim 1$ mm) continuum sources have been uncovered by recent unbiased surveys of, e.g., HUDF/GOODS-S using ALMA (e.g., Aravena et al., 2016; Dunlop et al., 2017; Hatsukade et al., 2018; Franco et al., 2018).

Combining with Cycle 4 data (Izumi et al., 2018a), we have observed seven HSC quasars in total with ALMA Band 6, and this object is the only one identified as a continuum emitter in the field, except for the target quasars. This detection rate (one continuum emitter in $7\times 0.135$ arcmin2 fields) seems to be smaller than recent measurements of the 1.2 mm number counts in the field (e.g., Aravena et al., 2016; Fujimoto et al., 2016). Indeed, the best-fit cumulative number count in Fujimoto et al. (2016) 777We adopted $S_{*}=2.35$ mJy, $\phi_{*}=1.54\times 10^{3}$ deg*-2*, $\alpha=-2.12$, and effective area = 80% of each field (see Table 5 of Fujimoto et al., 2016). predicts that we should have detected $\sim 11-16$ continuum sources (when corrected for errors) over our seven fields, or $\sim 1-3$ in each 0.135 arcmin2 field, given the depth of each one. This may be an overestimate for our fields because our high resolution observations are less sensitive to extended emission (see discussion in Fujimoto et al., 2017). We may also be subject to strong cosmic variance. Even so, our results will not support that these HSC quasars reside in overdense regions of submm sources. In this context, Champagne et al. (2018) reported no submm overdensity around a sample of 35 $z>6$ optically-luminous quasars (total effective area = 4.3 arcmin2), although their observations are considerably shallower than ours. It is noteworthy, however, that some studies at $z\lesssim 5$ suggest that luminous quasars tend to reside in overdense regions of emitters (Silva et al., 2015; Trakhtenbrot et al., 2017). Further observations are required to reconcile this discrepancy.

4 Discussion

4.1 The [C \emissiontypeII]-FIR luminosity relation

We discuss the [C \emissiontypeII] to FIR luminosity ratio of our quasars and various comparison samples here. The ratio quantifies the contribution of [C \emissiontypeII] line emission to the cooling of the cold ISM (the Milky Way value of $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ is $3\times 10^{-3}$, Carilli & Walter, 2013), but it has long been known that the $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ratio is an order of magnitude smaller in ULIRG-like FIR-luminous systems (e.g., Malhotra et al., 1997; Brauher et al., 2008; Stacey et al., 2010; Graciá-Carpio et al., 2011; Farrah et al., 2013; Díaz-Santos et al., 2013, 2017). This [C \emissiontypeII]-deficit trend has also been found in high-redshift quasars (e.g., Wang et al., 2013; Venemans et al., 2016; Decarli et al., 2018). Several processes may contribute to the deficit in quasars, including an AGN contribution to $L_{\rm FIR}$ (Sargsyan et al., 2014), reduction of C+ abundance due to AGN irradiation (Langer & Pineda, 2015), high gas surface densities of individual clouds (giving a high molecular-to-atomic gas fraction, Narayanan & Krumholz, 2017), etc., but an overwhelming factor appears to be the presence of a high FIR luminosity density region and/or high-temperature dust-emitting region, as has been invoked for local ULIRGs (e.g., Malhotra et al., 1997; Díaz-Santos et al., 2013, 2017, see also the discussion in the last part of this subsection).

Figure 9 shows the [C \emissiontypeII]-deficit trend with a compilation of galaxies at various redshifts, including $z>5.7$ quasars with available data (the 7 HSC quasars + 43 previously-studied objects). By adding optically-selected lower-luminosity ($L_{\rm FIR}\sim 10^{11}~{}L_{\odot}$) objects like our HSC quasars, we have increased the dynamic range of the plot, making it easier to see any correlation between $L_{\rm FIR}$ and $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$. Here the $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ratios of the HSC quasars are not drastically different from those of low redshift galaxies with comparable $L_{\rm FIR}$, except for the ULIRG-class object J2239+0207. The modest ratios found for the HSC quasars were also observed in comparably optically- and FIR-faint CFHQS quasars (Willott et al., 2013, 2015, 2017). Note that some of the previously-studied quasars may be biased toward low $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ values as they were originally selected based on their high $L_{\rm FIR}$. However, as recent ALMA follow-up studies have been performed for a number of optically-selected quasars (e.g., Decarli et al., 2018), the trend in Figure 9 becomes less biased. The Spearman rank correlation coefficient for the whole quasar sample in Figure 9, excluding objects with upper or lower limits, is fairly high ($\rho=-0.78$) with a null-hypothesis probability of $p=1.29\times 10^{-10}$. The relationship can be expressed as

[TABLE]

which is consistent with previous analysis (Willott et al., 2017; Izumi et al., 2018a).

In order to investigate the physical origin of the [C \emissiontypeII]-deficit in high-redshift quasars, we measured their rest frame [C \emissiontypeII] equivalent width (EW${}_{\rm[C\,\emissiontype{II}]}$) to explore the quasar contribution to the FIR light, defined as

[TABLE]

where the continuum emission is measured at frequencies close to that of the line. The use of EW${}_{\rm[C\,\emissiontype{II}]}$ has the advantage that it does not require any assumption about the shape of the IR SED. With quasar samples in Figure 9, we found a correlation between EW${}_{\rm[C\,\emissiontype{II}]}$ and the quasar absolute UV magnitude (Figure 10) as

[TABLE]

However, this correlation is marginal with $\rho=0.32$ and $p=0.03$; i.e., its significance is much weaker than the $L_{\rm FIR}$–$L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ correlation. Indeed, there are a number of optically-luminous quasars that show comparable EW${}_{\rm[C\,\emissiontype{II}]}$ to the much fainter HSC quasars (Figure 10). Thus, while the marginal correlation coefficient (0.32) implies a certain level of quasar contribution to the 1.2 mm continuum flux density (hence $L_{\rm FIR}$), it is not a prime driver of the [C \emissiontypeII]-deficit. Furthermore, given the positive correlation observed between quasar nuclear luminosity and $L_{\rm[C\,\emissiontype{II}]}$ (Figure 5), AGN irradiation does not play a primary role in the [C \emissiontypeII]-deficit.

Another plausible factor is the high FIR luminosity density and/or the existence of a high temperature dust-emitting region (e.g., Díaz-Santos et al., 2013, 2017; Smith et al., 2017; Herrera-Camus et al., 2018). In this context, Figure 11 investigates the dependence of the $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ratio on the FIR luminosity density ($\Sigma_{\rm FIR}=L_{\rm FIR}/2\pi R^{2}_{\rm FIR}$), computed using the two-dimensional Gaussian fit size in the FIR continuum map. In addition to the HSC quasars, we plot $z\gtrsim 6$ quasars compiled from Decarli et al. (2018) and Venemans et al. (2018). These quasars were selected from optical surveys, and thus are free from a $L_{\rm FIR}$-based selection bias. We then found a strong anti-correlation between $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ and $\Sigma_{\rm FIR}$:

[TABLE]

with $\rho=-0.78$ and $p=4.9\times 10^{-5}$, respectively. This strong trend is consistent with previous results for $z\gtrsim 6$ quasars (Decarli et al., 2018). It is noteworthy that the local ULIRGs also fit to this scenario well, as they typically have even smaller FIR emitting regions ($\lesssim 500$ pc, e.g., Soifer et al., 2000; Sakamoto et al., 2008; Imanishi et al., 2011) than high-$z$ submillimeter galaxies and quasar hosts having comparable $L_{\rm FIR}$ (e.g., Hodge et al., 2016; Fujimoto et al., 2017), which will lead to the their small $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ratios seen in Figure 9. Note that recent very high resolution ALMA observations enabled spatially-resolved measurements of $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$, which revealed similarly small ratios in the central ($r\lesssim 1$ kpc) high $\Sigma_{\rm FIR}$ regions of some submillimeter galaxies (SMGs, e.g., Gullberg et al., 2018; Rybak et al., 2019).

There are three intimately linked physical processes which may drive the [C \emissiontypeII]-deficit in high $\Sigma_{\rm FIR}$ regions. First process is driven by the increased radiation field, under which dust particles have more positive charge. This results in a reduction of the number of free electrons released from the dust particles, which contribute to [C \emissiontypeII] excitation (e.g., Malhotra et al., 1997; Negishi et al., 2001). Second, increased ionized-to-atomic hydrogen ratio will reduce the fraction of UV photons absorbed by gas, which then leads to the reduction of the $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ratio (Abel et al., 2009). Third, the temperature of the dust ($T_{\rm d}$) itself matters. Dust grains are heated to higher $T_{\rm d}$ in higher $\Sigma_{\rm FIR}$ regions, in which a larger number of ionizing photons is available, as seen in local ULIRGs. This greatly enhances $L_{\rm FIR}$ and thus reduces $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ (e.g., Díaz-Santos et al., 2013, 2017). To further test this possibility, however, shorter wavelength continuum observations are necessary to constrain $T_{\rm d}$ directly as we now assumed $T_{\rm d}=47$ K (Beelen et al., 2006) for all quasars in Figure 11 for consistency.

4.2 Early SMBH–host galaxy co-evolution

We now investigate the early co-evolution of SMBHs and their host galaxies, in both integrated and differential forms. The relevant properties of the seven HSC quasars with ALMA data are compared with those of previously observed $z\gtrsim 6$ quasars, to give a less-biased view on early mass assembly. We recall that the HSC quasars constitute the break ($M^{\star}_{\rm 1450}=-24.9$ mag) or further lower luminosity regime of the quasar luminosity function at $z\sim 6$ (Matsuoka et al., 2018c), which then represent the bulk of the quasar population at that era. Subsequent NIR follow-up observations started to reveal that they possess a wide range of $M_{\rm BH}$ ($\sim 10^{7.5}-10^{9}~{}M_{\odot}$, Onoue et al., 2019), which are therefore characterized by a wide range of Eddington ratio ($\sim 0.1-1$). Among the six HSC quasars reported in Onoue et al. (2019), four objects indeed have Eddington ratios of $0.16-0.24$. Note that dust obscuration does not play a major role in shaping the low-luminosity nature of most of the HSC quasars, as judged from the rest-UV SED modeling (Onoue et al., 2019).

4.2.1 Integrated form: $M_{\rm BH}-M_{\rm dyn}$

We first computed $M_{\rm BH}$ of J1208$-$0200 and J2239+0207 (Table 4.2.1) using the broad Mg \emissiontypeII-based virial mass calibration (the so-called single epoch method, Vestergaard & Osmer, 2009). Details of the procedure are described in Onoue et al. (2019). J2228+0152 has no Mg \emissiontypeII spectroscopy, so we assumed Eddington-limited accretion to derive its $M_{\rm BH}$. Note that this assumption is often made in $z\gtrsim 6$ quasar studies (e.g., Wang et al., 2013; Decarli et al., 2018) to compute the lower limit of $M_{\rm BH}$. The bolometric luminosity was calculated from the 1450 Å monochromatic luminosity with a correction factor of 4.4 (Richards et al., 2006). The estimated $M_{\rm BH}$ of the three HSC quasars presented here fall in the range $(1.1-11)\times 10^{8}~{}M_{\odot}$. Along with our Cycle 4 measurements, these HSC quasars populate the middle to lower regime of the $z\gtrsim 6$ quasar mass distribution observed thus far (Onoue et al., 2019). Hence we are indeed probing a quasar population less biased in terms of $M_{\rm BH}$ by observing these low-luminosity quasars.

Next, we estimate their host galaxy dynamical masses ($M_{\rm dyn}$) from the [C \emissiontypeII] spatial extents and line widths, by following the standard procedure used in $z\gtrsim 6$ quasar studies (Wang et al., 2013; Willott et al., 2015, 2017; Venemans et al., 2016; Izumi et al., 2018a). Here the [C \emissiontypeII] emission is assumed to originate from a thin rotating circular disk. The inclination angle ($i$) of the disk is determined from the axis ratio of the deconvolved Gaussian fit to the [C \emissiontypeII] emitting region. The circular velocity is expressed as $v_{\rm circ}=0.75{\rm FWHM}_{\rm[C\,\emissiontype{II}]}/\sin i$. The disk size is given by $D=1.5\times a_{\rm maj}$, where $a_{\rm maj}$ is the deconvolved size of the Gaussian major axis, and the factor 1.5 is used to account for spatially extended low level emission (e.g., Wang et al., 2010). The $M_{\rm dyn}$ within $D$ is then,

[TABLE]

The resultant values for both $M_{\rm dyn}\sin^{2}i$ and $M_{\rm dyn}$ are listed in Table 4.2.1. Note that formal errors on $M_{\rm dyn}$ are not given due to multiple unconstrained uncertainties including the inclination angles and the true geometry of the [C \emissiontypeII] emitting regions. Moreover, if the host galaxies have dispersion-dominant gas dynamics, their $M_{\rm dyn}$ will be significantly smaller than the ones derived with the rotating disk assumption (e.g., Venemans et al., 2017a), which will affect our discussion in the following. Unfortunately, it is difficult at this moment to investigate if the host galaxies discussed here are indeed rotation-dominant systems, as the spatial resolutions obtained thus far are only modest in most cases. Further higher resolution ALMA observations will reveal the dynamical nature of $z\gtrsim 6$ quasars (see recent examples in, Venemans et al., 2019; Wang et al., 2019), and give a better insight on the early co-evolution. Despite these large uncertainties, we hereafter use $M_{\rm dyn}$ as a surrogate for the stellar mass ($M_{\star}$) of the quasar hosts, as is often done in high redshift quasar studies (e.g., Wang et al., 2013; Willott et al., 2015; Venemans et al., 2016; Izumi et al., 2018a). The derived values of $M_{\rm dyn}$ exceed $10^{10}~{}M_{\odot}$, or even $10^{11}~{}M_{\odot}$ for the case of J2239+0207, which lie at the massive end of the $M_{\star}$-distribution for $z\sim 6$ galaxies in general (e.g., Grazian et al., 2015). Thus, the host galaxies of the HSC quasars are among the most evolved systems known at $z\sim 6$ in terms of their galaxy-masses.

We also compiled $M_{\rm BH}$ and $M_{\rm 1450}$ values of other $z\geq 5.7$ quasars from the literature (Willott et al., 2010b; De Rosa et al., 2014; Kashikawa et al., 2015; Venemans et al., 2015a; Bañados et al., 2016; Jiang et al., 2016; Shao et al., 2017; Mazzucchelli et al., 2017; Decarli et al., 2018). Most of the $M_{\rm BH}$ values from the literature were indeed measured with the same Mg \emissiontypeII-based single epoch method that we have used above (Vestergaard & Osmer, 2009). Exceptions are J1319+0950 (Shao et al., 2017), as well as J0842+1218, J1207+0630, and J2310+1855 (Shen et al., 2019), for which we recalculated $M_{\rm BH}$ with the Vestergaard & Osmer (2009) calibration. For those without available $M_{\rm BH}$ data, we assumed the Eddington-limited accretion. Although this assumption may not be valid for our HSC quasar J2228+0152 as it shows a rather quiescent nature (non-detection of FIR continuum emission), we nevertheless keep this assumption as we would like to take a statistical approach. The $M_{\rm dyn}$ of these quasars from the literature were also computed in the same way we adopted for our HSC quasars, i.e., the thin disk approximation with a 2D Gaussian decomposition (Wang et al., 2013; Willott et al., 2013, 2015, 2017; Venemans et al., 2016; Shao et al., 2017; Decarli et al., 2018). Quasars that do not have such a Gaussian decomposition were excluded from our sample. One exception is the second highest redshift quasar known, J1120+0641 ($z=7.08$, Mortlock et al., 2011). Previous high resolution [C \emissiontypeII] observations revealed that its host galaxy does not show ordered rotation, and an upper limit on $M_{\rm dyn}$ was provided by applying the virial theorem (Venemans et al., 2017a).

In Figure 12 we display the relation between $M_{\rm BH}$ and $M_{\rm dyn}$ for the 40 above-mentioned quasars at $z\gtrsim 6$, overlaid with the local $M_{\rm BH}-M_{\rm bulge}$ relation after equating $M_{\rm dyn}$ to $M_{\rm bulge}$ (Kormendy & Ho, 2013). These quasars are color-coded by their $M_{\rm 1450}$. Regarding the optically-luminous objects, this figure supports conclusions in previous works (e.g., Wang et al., 2013; Venemans et al., 2016). That is, the luminous quasars ($M_{\rm 1450}\lesssim-25$ mag) typically have overmassive SMBHs with respect to the local relation, while the discrepancy becomes less evident at $M_{\rm dyn}\gtrsim 10^{11}~{}M_{\odot}$: some luminous quasars even show comparable $M_{\rm BH}/M_{\rm dyn}$ ratios to the local relation. It is likely for such high-mass host galaxies that past multiple major and minor mergers already led to this convergence.

On the other hand, most of the low-luminosity quasars with $M_{\rm 1450}\gtrsim-25$ mag show comparable ratios to, or even lower ratios than, the local relation, particularly at a range of $M_{\rm dyn}\gtrsim 4\times 10^{10}~{}M_{\odot}$. The existence of the undermassive SMBHs even implies an evolutionary path, in which galaxies grow earlier than SMBHs, such as expected in a standard merger-induced evolution model (e.g., Hopkins et al., 2008). In this high $M_{\rm dyn}$ range, our result demonstrates that previous works on luminous quasars have been largely biased toward the most massive SMBHs (Lauer et al., 2007; Schulze & Wisotzki, 2014), easily resulting in objects lying above the local $M_{\rm BH}-M_{\rm dyn}$ relation. Therefore, our result highlights the power of the sensitive Subaru survey to probe the fainter part of the quasar luminosity- and mass-functions and to reveal the nature of early co-evolution of black holes and their host galaxies in a less biased way.

However, we also found that the low-luminosity quasars at $M_{\rm dyn}\lesssim 3\times 10^{10}~{}M_{\odot}$ start to show overmassive $M_{\rm BH}$ with respect to the local relation. Regarding the HSC quasars, we argue that this is still due to our selection bias: we basically selected objects with $M_{\rm 1450}<-24$ for the targets in our NIR follow-up observations (Onoue et al., 2019), which corresponds to $M_{\rm BH}\sim 10^{8}~{}M_{\odot}$ if it is accreting at the Eddington-limit. Meanwhile, the $M_{\rm BH}$ expected at $M_{\rm dyn}\sim(1-3)\times 10^{10}~{}M_{\odot}$ from the local relation is well smaller than $10^{8}~{}M_{\odot}$. Thus, given the large scatter of the Eddington ratio distribution at $z\sim 6$ ($\sim 0.01-1$, Onoue et al., 2019), we clearly need to observe much fainter objects to surely probe the $M_{\rm BH}<10^{8}~{}M_{\odot}$ region and to reveal the unbiased shape of the relation. Future sensitive observations with the James Webb Space Telescope (JWST) or ground-based extremely large telescopes will allow to probe down this very low mass region even at $z>6$.

Given their low nuclear luminosities (hence mass accretion rate onto SMBHs) and SFRs, the low-luminosity quasars studied here will not move drastically in the $M_{\rm BH}-M_{\rm dyn}$ plane even over the next $\sim 10$ Myr (typical life-time of high-redshift quasars expected with the transverse proximity effect, Borisova et al., 2016), unless there is a rich supply of gas which infalls to the galaxy. This is particularly true if the quasar duty-cycle is $\lesssim 10^{-2}$ as suggested for $z\sim 6$ objects recently (Chen & Gnedin, 2018). Note that Izumi et al. (2018a) argued that the four HSC quasar hosts studied in ALMA Cycle 4 are on or even below the star formation main sequence at $z\sim 6$ (e.g., Salmon et al., 2015). This trend also holds for the Cycle 5 samples. Therefore, these low-luminosity quasars may already be near the end of their SMBH growth and galaxy growth, and are transforming into a quiescent phase even at $z\gtrsim 6$. This is consistent with the conclusion of Onoue et al. (2019), who raised a similar argument based on the Eddington ratio measurements for HSC quasars. If this is true, there must be a quite rapid physical process to realize the $M_{\rm BH}-M_{\rm dyn}$ relation at that high redshift, such as expected in a merger-induced evolution (e.g., Di Matteo et al., 2005; Hopkins et al., 2006, 2008). Indeed, a recent very high resolution simulation based on this scheme suggests that even a quasar at $z=7$ follows the local co-evolution relation (Lupi et al., 2019).

4.2.2 Differential form: $L_{\rm Bol}-L_{\rm FIR}$

The $M_{\rm BH}-M_{\rm dyn}$ relationship studied in § 4.2.1 reflects the integrated history of past mass accumulation. We now investigate the on-going mass accumulation of $z\gtrsim 6$ quasars by comparing quasar bolometric luminosity ($L_{\rm Bol}$) and $L_{\rm FIR}$. As $L_{\rm Bol}$ and $L_{\rm FIR}$ can be converted to the growth rate of the central SMBH ($\equiv\dot{M}_{\rm BH}$) and the SFR of the host galaxy (see § 3.2), respectively, the $L_{\rm Bol}-L_{\rm FIR}$ relation indicates a differential form of the Maggorian relation (e.g., Willott et al., 2013; Drouart et al., 2014).

We again compiled literature data for $M_{\rm 1450}$ (Willott et al., 2003, 2010b; De Rosa et al., 2014; Kashikawa et al., 2015; Wu et al., 2015; Jiang et al., 2016; Mazzucchelli et al., 2017; Shao et al., 2017; Bañados et al., 2016, 2018; Decarli et al., 2018; Venemans et al., 2018) and $L_{\rm FIR}$ (Bertoldi et al., 2003a; Petric et al., 2003; Maiolino et al., 2005; Wang et al., 2007, 2008, 2011a, 2013, 2016; Omont et al., 2013; Willott et al., 2013, 2015, 2017; Bañados et al., 2015; Venemans et al., 2016; Shao et al., 2017; Venemans et al., 2017a, c, 2018; Mazzucchelli et al., 2017; Decarli et al., 2017, 2018) of $z\gtrsim 6$ quasars. In addition to the objects we compiled in § 4.2.1, we appended a number of quasars, mostly with single dish $L_{\rm FIR}$ measurements. Thus the total number of quasars on this analysis is increased to 97, including the seven HSC quasars. Their $L_{\rm Bol}$ and $L_{\rm FIR}$ values were computed in the same manner as described earlier.

Figure 13 shows the $L_{\rm Bol}-L_{\rm FIR}$ relationship for the $z\gtrsim 6$ quasars on a logarithmic scale. While a substantial fraction of the sample quasars only have upper limits on $L_{\rm FIR}$, there seems to be a marginal correlation between the two quantities. Fitting a linear regression, including the upper limits, we find

[TABLE]

using the IRAF STSDAS package888STSDAS is a product of the Space Telescope Science Institute, which is operated by AURA for NASA.. This relation is consistent with that found in Wang et al. (2011a) for stacked averages of quasars at $2<z<7$. Note that this correlation is strengthened by adding low-luminosity objects like the HSC quasars, as the Spearman rank coefficient of the $\log L_{\rm Bol}-\log L_{\rm FIR}$ relation decreases from 0.51 (full sample in Figure 13) to 0.34 when restricting the sample to optically luminous ($L_{\rm Bol}>10^{13}~{}L_{\odot}$) quasars only. This is qualitatively consistent with Venemans et al. (2018), who reported no significant correlation between these two quantities when focusing on optically luminous quasars ($L_{\rm Bol}\gtrsim 10^{13}~{}L_{\odot}$). The physical origin of this kind of correlation remains unclear, but one simple explanation is that both the quasar activity and star-forming activity are supported by a common source of gas, such as a $\sim 100$ pc scale circumnuclear gas disk seen in nearby Seyfert galaxies (e.g., Izumi et al., 2016), as also found in some simulations (e.g., Hopkins & Quataert, 2010). Note that, whatever the origin of the correlation is, time-variability of the quasar luminosity would weaken the correlation, given that star-formation occurs on a much longer time-scale than does nuclear mass accretion (Novak et al., 2011; Hickox et al., 2014).

In Figure 13, we can define a line of parallel growth (e.g., Willott et al., 2013; Drouart et al., 2014). This indicates an evolutionary path in which SMBHs grow in tandem with their host galaxies to give rise to the local $M_{\rm BH}-M_{\rm bulge}$ relation. Here we adopt the local relation of Kormendy & Ho (2013), i.e., $M_{\rm BH}\simeq 0.005\times M_{\rm bulge}$. The parallel growth line is thus equivalent to $\dot{M}_{\rm BH}=0.005\times SFR$ (Figure 13). The line moves upward in the figure if we adopt the McConnell & Ma (2013) relation instead ($M_{\rm BH}\simeq 0.003\times M_{\rm bulge}$). $L_{\rm Bol}$ and $\dot{M}_{\rm BH}$ are equated as $\dot{M}_{\rm BH}=((1-\eta)/\eta)\cdot(L_{\rm Bol}/c^{2})$, where $\eta$ is the radiative efficiency ($\sim 0.1$, Volonteri & Rees, 2005; Zhang & Lu, 2019) and $c$ is the speed of light. It is evident that most of the $z\gtrsim 6$ quasars discussed here, including the low-luminosity objects, are on or below the parallel growth line (quasar-dominant phase) in this differential Magorrian plot.

This result supports the evolutionary scenario that sequentially links vigorous starburst and SMBH growth (e.g., Sanders et al., 1988; Hopkins et al., 2006, 2008; Alexander et al., 2008; Alexander & Hickox, 2012; Simpson et al., 2012). In this scenario, the peak of the quasar luminosity comes after the peak of the dust-obscured starburst phase, in which AGN feedback would expel the surrounding dusty ISM and quenches the starburst. As we saw, the low-luminosity quasars, including the HSC quasars, particularly hosted in massive galaxies ($M_{\rm dyn}\gtrsim 4\times 10^{10}~{}M_{\odot}$) are located on or below the local Magorrian relation (Figure 12). Given the low $M_{\rm dust}$ of these low-luminosity quasars (mostly on the order of 107 $M_{\odot}$), an expected gas mass fraction in $M_{\rm dyn}$ after multiplying a gas-to-dust mass ratio (e.g., 100, Draine et al., 2007) is small, indicating that stellar masses dominate their $M_{\rm dyn}$. Meanwhile, the on-going star formation in their hosts are also weak; they are entering a quiescent phase (Izumi et al., 2018a). Note that this view is sensitive to the assumed $T_{\rm d}$ to compute $L_{\rm FIR}$ from the single ALMA photometry. We recall, however, that now we assume the same high $T_{\rm d}$ for the HSC quasars as the one canonically adopted for the optically and FIR much brighter quasars (= 47 K). It may be more plausible for the HSC quasars to have lower $T_{\rm d}$ such as 35 K (see also Section 3.2). If this is true, the HSC quasars will reside in a further quasar-dominant region on Figure 13. Future multi-band ALMA observations will enable to better constrain the $T_{\rm d}$ of these HSC quasars, as has been performed in some luminous quasars (e.g., Wang et al., 2019). In summary, we suggest that these HSC quasars would have experienced an earlier vigorous starburst phase to generate their host galaxy masses, prior to the currently observed quasar-dominant phase, in order for them to reach the local Magorrian relation.

Prime candidates for such starbursting predecessors are SMGs at even higher redshifts. So far, only three SMGs have been spectroscopically identified at $z>6$: SPT0311$-$58 at $z=6.90$ (Marrone et al., 2018), HFLS3 at $z=6.34$ (Riechers et al., 2013), and HATLAS G09 83808 at $z=6.03$ (Zavala et al., 2018). These are known to host extreme starbursts with SFR $\sim 3000$ $M_{\odot}$ yr*-1* (for the former two) and $\sim 380$ $M_{\odot}$ yr*-1* (for G09 83808, after correcting for gravitational magnification), which has allowed their host galaxies to grow as massive as $\sim 10^{11}~{}M_{\odot}$ within a fairly short time. This picture is consistent with recent simulation work of merger-induced galaxy evolution (Ginolfi et al., 2019), which suggested that $z\sim 7$ hyper-luminous infrared galaxies (HyLIRGs) are indeed the ancestors of $z\sim 6$ luminous quasars.

Whilst it is quite difficult to determine if there is an AGN in these gas-rich and dusty systems at $z>6$ using a current instrument, studies in the nearby universe have shown that IR-luminous or gas-rich systems tend to possess AGNs with high Eddington ratio (e.g., Hao et al., 2005; Xia et al., 2012; Izumi, 2018b). Such a high Eddington ratio (or even super-Eddington) phase may be essential for cosmic SMBH growth (e.g., Kawaguchi et al., 2004; Di Matteo et al., 2008). Furthermore, recent ALMA observations, supplemented by deep Chandra 7 Ms survey data, have revealed a high fraction of AGNs in moderate luminosity SMGs (i.e., 90${}^{+8}_{-19}$% of ULIRG-class objects and 57${}^{+23}_{-25}$% of LIRG-class objects are AGNs) at $z\sim 1-3$ (Ueda et al., 2018). About two-thirds of their sample (25 SMGs in total) are in the starburst-dominant phase in terms of the differential Magorrian relation (see also Wang et al., 2013), which is qualitatively consistent with some earlier works (e.g., Alexander et al., 2005, 2008). Our results, along with these previous findings, support the idea that SMBHs and their host galaxies do not actually co-evolve in a synchronized way. Rather, there seems to be an evolutionary path from the starburst-dominant phase to the quasar-dominant phase.

5 Summary

As a part of the SHELLQs project, a large optical survey of low-luminosity quasars ($M_{\rm 1450}\gtrsim-25$ mag) at $z\gtrsim 6$ with the Subaru Hyper Suprime-Cam (HSC), we performed Cycle 5 ALMA follow-up observations toward three HSC quasars, in addition to our previous Cycle 4 study of four quasars. We thus continue to grow the sample of hosts of low-luminosity, high redshift quasars, from the pioneering works of Willott et al. (2013, 2015, 2017). The main findings of this paper strengthen our previous arguments in Izumi et al. (2018a), and are summarized as follows.

We detected [C \emissiontypeII] line emission in all three target HSC quasars, with [C \emissiontypeII] luminosities between $(2.4-9.5)\times 10^{8}~{}L_{\odot}$. These are consistent with the Cycle 4 measurements for other HSC quasars, but are one order of magnitude smaller than those measured in optically luminous quasars. The [C \emissiontypeII] line width shows no clear dependence on the quasar luminosity.

• 2.

Within a common 1\arcsec.0 aperture, we detected underlying FIR continuum emission from two of the three target quasars. It is intriguing that one HSC quasar shows a ULIRG-like FIR luminosity ($L_{\rm FIR}\simeq 2\times 10^{12}~{}L_{\odot}$), while another object was not detected even with our sensitive ALMA observations (3$\sigma$ limit of $L_{\rm FIR}<9\times 10^{10}~{}L_{\odot}$). There is a wide spread in $L_{\rm FIR}$ among the HSC quasars, even though their quasar luminosities in rest-UV are comparable. These $L_{\rm FIR}$ are again an order of magnitude smaller than those typically found in optically luminous quasars.

• 3.

The spatial extents of the (barely resolved) [C \emissiontypeII]-emitting regions of these HSC quasars are $\sim 2-4$ kpc, roughly consistent with previous measurements of both optically luminous and low-luminosity quasars. Thus the difference in $L_{\rm FIR}$ between optically luminous quasars and our HSC quasars roughly translates to the difference in FIR luminosity density.

• 4.

The $L_{\rm[C\,\emissiontype{II}]}/L_{\rm FIR}$ ratios of the HSC quasars (except for the ULIRG-class object J2239+0207) are consistent with local star-forming galaxies. However, the ratio is an order of magnitude smaller in some optically and FIR-luminous quasars from the literature. We suggest that the high FIR luminosity density of the observed region and/or existence of a high temperature dust-emitting region may be the physical origin of this [C \emissiontypeII] deficit.

• 5.

From the dynamical measurements based on the thin rotating disk assumption, along with the $M_{\rm BH}$ data, we found that most of the HSC quasars and similarly low-luminosity $z\gtrsim 6$ CFHQS quasars ($M_{\rm 1450}\gtrsim-25$ mag) tend to lie on or even below the local $M_{\rm BH}-M_{\rm bulge}$ relation particularly at the high $M_{\rm dyn}$ range ($\gtrsim 4\times 10^{10}~{}M_{\odot}$). This would require a quite rapid mechanism of co-evolution, such as merger-induced galaxy evolution, to grow both the hosts and black holes, given the high redshifts ($z\gtrsim 6$) of our sample quasars. Optically luminous quasars ($M_{\rm 1450}\lesssim-25$ mag), on the other hand, host overmassive SMBHs with respect to the local relation, while the discrepancy becomes less evident at the massive-end of $M_{\rm dyn}\gtrsim 10^{11}~{}M_{\odot}$. However, we also implied our current limitation to fully probe down the low-$M_{\rm BH}$ range ($<10^{8}~{}M_{\odot}$), which will be the subject of future sensitive NIR observations.

• 6.

All $z\gtrsim 6$ quasars compiled in this work are located in the quasar-dominant region of the $L_{\rm Bol}-L_{\rm FIR}$ plane. As low-luminosity quasars seem to be consistent with the local $M_{\rm BH}-M_{\rm bulge}$ relation at least at the high galaxy-mass range, and are transforming into a quiescent population in terms of both SMBH growth and stellar mass accumulation, there should have been a starburst-dominant phase for them to gain their large host dynamical masses ($\sim 10^{10-11}~{}M_{\odot}$), prior to the currently observed quasar phase. Submillimeter galaxies (SMGs) at even higher redshifts ($z\gtrsim 7$) would be the prime candidate for such a starbursting ancestor.

In this work, we demonstrated the importance of investigating low-luminosity quasars to probe the nature of early co-evolution of black holes and host galaxies in an unbiased way. The trends of low-luminosity quasars shown above are indeed clearly different from those of optically-luminous quasars, although our conclusions are based on the small sample. As the number of low-luminosity HSC quasars at $z>6$ is growing dramatically, we can statistically confirm the trends we found thus far with sensitive ALMA surveys.

{ack}

We deeply appreciate the anonymous referee for his/her thorough reading and useful comments which greatly improved this paper. We also appreciate Toshihiro Kawaguchi for his intensive advice on this work. This paper makes use of the following ALMA data: ADS/JAO.ALMA 2017.1.00541.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/ NRAO and NAOJ.

The Hyper Suprime-Cam (HSC) collaboration includes the astronomical communities of Japan and Taiwan, and Princeton University. The HSC instrumentation and software were developed by the National Astronomical Observatory of Japan (NAOJ), the Kavli Institute for the Physics and Mathematics of the Universe (Kavli IPMU), the University of Tokyo, the High Energy Accelerator Research Organization (KEK), the Academia Sinica Institute for Astronomy and Astrophysics in Taiwan (ASIAA), and Princeton University. Funding was contributed by the FIRST program from Japanese Cabinet Office, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), the Japan Society for the Promotion of Science (JSPS), Japan Science and Technology Agency (JST), the Toray Science Foundation, NAOJ, Kavli IPMU, KEK, ASIAA, and Princeton University.

This paper makes use of software developed for the Large Synoptic Survey Telescope. We thank the LSST Project for making their code available as free software at http://dm.lsstcorp.org.

The Pan-STARRS1 Surveys (PS1) have been made possible through contributions of the Institute for Astronomy, the University of Hawaii, the Pan-STARRS Project Office, the Max-Planck Society and its participating institutes, the Max Planck Institute for Astronomy, Heidelberg and the Max Planck Institute for Extraterrestrial Physics, Garching, The Johns Hopkins University, Durham University, the University of Edinburgh, Queen’s University Belfast, the Harvard-Smithsonian Center for Astrophysics, the Las Cumbres Observatory Global Telescope Network Incorporated, the National Central University of Taiwan, the Space Telescope Science Institute, the National Aeronautics and Space Administration under Grant No. NNX08AR22G issued through the Planetary Science Division of the NASA Science Mission Directorate, the National Science Foundation under Grant No. AST-1238877, the University of Maryland, and Eotvos Lorand University (ELTE).

T.I., M.I., Y.M., N.K., T.N., K.K., and Y.T. are supported by Japan Society for the Promotion of Science (JSPS) KAKENHI grant numbers JP17K14247, JP15K05030, JP17H04830, JP15H03645, JP17H01114, and JP17H06130, respectively. Y.M. is also supported by the Mitsubishi Foundation grant No. 30140. T.G. acknowledges the support by the Ministry of Science and Technology of Taiwan through grant105-2112-M-007-003-MY3. T.I. is supported by the ALMA Japan Research Grant of NAOJ Chile Observatory, NAOJ-ALMA-217.