#  SILVERRUSH. XI. Intensity Mapping for Ly$\alpha$ Emission Extending over
100-1000 comoving kpc around $z\sim 2-7$ LAEs with Subaru HSC-SSP and CHORUS
Data

Shotaro Kikuchihara Institute for Cosmic Ray Research, The University of
Tokyo, 5-1-5 Kashiwanoha, Kashiwa, Chiba 277-8582, Japan  Department of
Astronomy, Graduate School of Science, The University of Tokyo, 7-3-1 Hongo,
Bunkyo-ku, Tokyo 113-0033, Japan<EMAIL_ADDRESS>Yuichi Harikane
Institute for Cosmic Ray Research, The University of Tokyo, 5-1-5 Kashiwanoha,
Kashiwa, Chiba 277-8582, Japan  Department of Physics and Astronomy,
University College London, Gower Street, London WC1E 6BT, UK  Masami Ouchi
Institute for Cosmic Ray Research, The University of Tokyo, 5-1-5 Kashiwanoha,
Kashiwa, Chiba 277-8582, Japan  National Astronomical Observatory of Japan,
2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan  Graduate University for Advanced
Studies (SOKENDAI), 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan  Kavli
Institute for the Physics and Mathematics of the Universe (Kavli IPMU, WPI),
The University of Tokyo, 5-1-5 Kashiwanoha, Kashiwa, Chiba 277-8583, Japan
Yoshiaki Ono Institute for Cosmic Ray Research, The University of Tokyo, 5-1-5
Kashiwanoha, Kashiwa, Chiba 277-8582, Japan  Takatoshi Shibuya  Kitami
Institute of Technology, 165 Koen-cho, Kitami, Hokkaido 090-8507, Japan
Ryohei Itoh  Institute for Cosmic Ray Research, The University of Tokyo, 5-1-5
Kashiwanoha, Kashiwa, Chiba 277-8582, Japan  Department of Physics, Graduate
School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo
113-0033, Japan  Ryota Kakuma  Institute for Cosmic Ray Research, The
University of Tokyo, 5-1-5 Kashiwanoha, Kashiwa, Chiba 277-8582, Japan
Department of Astronomy, Graduate School of Science, The University of Tokyo,
7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan  Akio K. Inoue Waseda Research
Institute for Science and Engineering, Faculty of Science and Engineering,
Waseda University, 3-4-1 Okubo, Shinjuku, Tokyo 169-8555, Japan  Department of
Physics, School of Advanced Science and Engineering, Faculty of Science and
Engineering, Waseda University, 3-4-1 Okubo, Shinjuku, Tokyo 169-8555, Japan
Haruka Kusakabe Observatoire de Genève, Université de Genève, 51 chemin de
Pégase, 1290 Versoix, Switzerland  Kazuhiro Shimasaku Department of Astronomy,
Graduate School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku,
Tokyo 113-0033, Japan  Research Center for the Early Universe, Graduate School
of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033,
Japan  Rieko Momose Department of Astronomy, Graduate School of Science, The
University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan  Yuma
Sugahara National Astronomical Observatory of Japan, 2-21-1 Osawa, Mitaka,
Tokyo 181-8588, Japan  Waseda Research Institute for Science and Engineering,
Faculty of Science and Engineering, Waseda University, 3-4-1 Okubo, Shinjuku,
Tokyo 169-8555, Japan  Satoshi Kikuta Center for Computational Sciences,
University of Tsukuba, Ten-nodai, 1-1-1 Tsukuba, Ibaraki 305-8577, Japan  Shun
Saito Institute for Multi-messenger Astrophysics and Cosmology, Department of
Physics, Missouri University of Science and Technology, 1315 N. Pine Street,
Rolla, MO 65409, USA  Nobunari Kashikawa Department of Astronomy, Graduate
School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo
113-0033, Japan  Research Center for the Early Universe, Graduate School of
Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033,
Japan  Haibin Zhang Department of Astronomy, Tsinghua University, No. 1
Qinghuayuan, Beijing 100084, China  Chien-Hsiu Lee NSF’s National Optical-
Infrared Astronomy Research Laboratory, Tucson, AZ 85719, USA

###### Abstract

We conduct intensity mapping to probe for extended diffuse Ly$\alpha$ emission
around Ly$\alpha$ emitters (LAEs) at $z\sim 2-7$, exploiting very deep ($\sim
26$ mag at $5\sigma$) and large-area ($\sim 4.5$ deg2) Subaru/Hyper Suprime-
Cam narrow-band (NB) images and large LAE catalogs consisting of a total of
1781 LAEs at $z=2.2$, $3.3$, $5.7$, and $6.6$ obtained by the HSC-SSP
SILVERRUSH and CHORUS projects. We calculate the spatial correlations of these
LAEs with $\sim 1-2$ billion pixel flux values of the NB images, deriving the
average Ly$\alpha$ surface brightness (${\rm SB_{Ly\alpha}}$) radial profiles
around the LAEs. By carefully estimating systematics such as fluctuations of
sky background and point spread functions, we detect diffuse Ly$\alpha$
emission ($\sim 10^{-20}-10^{-19}$ erg s-1 cm-2 arcsec-2) at $100-1000$
comoving kpc around $z=3.3$ LAEs at the $4.1\sigma$ level and tentatively
($\sim 2\sigma$) at the other redshifts, beyond the virial radius of a dark-
matter halo with a mass of $10^{11}\ M_{\odot}$. While the observed ${\rm
SB_{Ly\alpha}}$ profiles have similar amplitudes at $z=2.2-6.6$ within the
uncertainties, the intrinsic ${\rm SB_{Ly\alpha}}$ profiles (corrected for the
cosmological dimming effect) increase toward high redshifts. This trend may be
explained by increasing hydrogen gas density due to the evolution of the
cosmic volume. Comparisons with theoretical models suggest that extended
Ly$\alpha$ emission around a LAE is powered by resonantly scattered Ly$\alpha$
photons in the CGM and IGM that originates from the inner part of the LAE,
and/or neighboring galaxies around the LAE.

galaxies: formation — galaxies: evolution — galaxies: high-redshift —
galaxies: halos — intergalactic medium

††journal: ApJ

##  1 Introduction

The gas surrounding a galaxy, called the circumgalactic medium (CGM), falls
into the galaxy and triggers star formation activity, and is subsequently
ejected from the galaxy due to outflows (e.g., Tumlinson et al., 2017; Péroux
& Howk, 2020). The hydrogen gas inside the CGM can be traced by Lyman-alpha
(Ly$\alpha$) emission, which is observed as a Ly$\alpha$ halo (LAH).
Therefore, observing LAHs is key to understanding the properties and
kinematics of the CGM, and eventually provide information on galaxy formation
and evolution.

Many studies have detected LAHs around nearby galaxies (e.g., Östlin et al.,
2009; Hayes et al., 2013, 2014). At high-redshift ($z>2$), meanwhile, LAHs
have been identified mainly around massive galaxies, such as Lyman break
galaxies (e.g., Hayashino et al., 2004; Swinbank et al., 2007; Steidel et al.,
2011) and quasars (e.g., Goto et al., 2009; Cantalupo et al., 2014; Martin et
al., 2014; Borisova et al., 2016; Arrigoni Battaia et al., 2019; Kikuta et
al., 2019; Zhang et al., 2020). However, it remains difficult to detect
diffuse emission around less massive star-forming galaxies (SFGs), such as
Ly$\alpha$ emitters (LAEs), at high redshift, due to their faintness and
sensitivity limits.

To overcome this difficulty, Rauch et al. (2008), for example, performed a
very deep (92 hr) long-slit observation with the ESO Very Large Telescope
(VLT)/FOcal Reducer and low dispersion Spectrograph 2 (FORS2) that reached a
$1\sigma$ surface brightness (SB) detection limit of $8\times 10^{-20}$ erg
cm-2 s-1 arcsec-2. They investigated 27 LAEs at $z=2.67-3.75$, identifying
Ly$\alpha$ emission extending over 26 physical kpc (pkpc) around one of the
LAEs. Individual detections of many high redshift LAHs have been enabled by
the advent of the Multi-Unit Spectroscopic Explorer (MUSE) installed on the
VLT. Recently, Leclercq et al. (2017) identified individual LAHs around 145
LAEs at $z=3-6$ in the Hubble Ultra Deep Field (HUDF) with the VLT/MUSE (see
also Wisotzki et al., 2016). Their data reached a SB limit of $\lesssim
10^{-19}$ erg s-1 cm-2 arcsec-2 at radii of $>10$ pkpc. The magnification
effect caused by gravitational lensing was also utilized in conjunction with
the MUSE for studying individual LAHs (Patrício et al., 2016; Smit et al.,
2017; Claeyssens et al., 2019).

A stacking method has been widely used to obtain averaged LAH profiles with
high signal-to-noise (S/N) ratios (e.g., Matsuda et al., 2012; Momose et al.,
2014, 2016; Xue et al., 2017; Wisotzki et al., 2018; Wu et al., 2020). For
example, Momose et al. (2014) stacked $>100$ narrow-band (NB) images around
LAEs at $z=2.2-6.6$ using Subaru Telescope/Suprime-Cam (SC) data,
investigating LAH profiles up to $\sim 50$ pkpc radial scales. Matsuda et al.
(2012) and Momose et al. (2016) investigated the LAH size dependence on LAE
properties, such as Ly$\alpha$ luminosity, rest-frame ultra-violet (UV)
magnitude, and overdensity, at $z=3.1$ and 2.2, respectively. We note that
some studies (e.g. Bond et al., 2010; Feldmeier et al., 2013; Jiang et al.,
2013) reported no evidence of extended Ly$\alpha$ emission at $z>2$.

Another approach is the intensity mapping technique (Kovetz et al. 2017 for a
review; see also Carilli 2011; Gong et al. 2011; Silva et al. 2013; Pullen et
al. 2014; Comaschi & Ferrara 2016a, b; Li et al. 2016; Fonseca et al. 2017),
which utilizes cross-correlation functions between objects and their emission
or absorption spectra. This technique enables us to detect signals from
targeted galaxies with a high S/N ratio by efficiently estimating and removing
contaminating signals from foreground interlopers. Croft et al. (2016, 2018)
derived cross-correlation functions between the Ly$\alpha$ emission and quasar
positions at $z=2-3.5$ using data from the Sloan Digital Sky Survey (SDSS;
Eisenstein et al., 2011) Baryon Oscillation Spectroscopic Survey (BOSS; Dawson
et al., 2013). This work enabled the detection of positive signals up to a
$\sim 15\ h^{-1}$ comoving Mpc (cMpc) radial scale.

Recently, Kakuma et al. (2021) applied the intensity mapping technique to LAEs
at $z=5.7$ and 6.6 using Subaru/Hyper Suprime-Cam (HSC) data. They tentatively
identified very diffuse ($\sim 10^{-20}\ {\rm erg\ s^{-1}\ cm^{-2}\
arcsec^{-2}}$) Ly$\alpha$ emission extended around the LAEs over the radial
scale of the virial radius ($R_{\rm vir}$) of a dark-matter halo (DMH).
Although their finding has shed light on the potential existence of
Ly$\alpha$-emitting hydrogen gas beyond $R_{\rm vir}$ around LAEs in the
reionization epoch, it remains an open question whether LAEs host such
extended structures even at lower redshifts ($z<5$). In this paper, we exploit
NB images offered in two Subaru/HSC surveys, the Cosmic HydrOgen Reionization
Unveiled with Subaru (CHORUS; Inoue et al., 2020) and the Subaru Strategic
Program (HSC-SSP; Aihara et al., 2019), which enable us to trace Ly$\alpha$
emission from LAEs across $z=2.2-6.6$. Our goal is to systematically
investigate diffuse Ly$\alpha$ emission extended beyond $R_{\rm vir}$ around
$z=2.2-6.6$ LAEs, taking advantage of the intensity mapping technique and
ultra-deep images from the CHORUS and HSC-SSP projects.

Another open question about extended Ly$\alpha$ emission is its physical
origin (see review by Ouchi et al. 2020 and Figure 15 of Momose et al. 2016).
Theoretical studies have suggested several physical processes producing a LAH
around a galaxy, which can be attributed mainly to 1) resonant scattering and
2) in-situ production. 1) Resonant scattering: Ly$\alpha$ photons are produced
in the interstellar medium (ISM) of the galaxy, and then resonantly scattered
by neutral hydrogen gas while escaping the galaxy into the CGM and
intergalactic medium (IGM) (e.g., Laursen & Sommer-Larsen, 2007; Laursen et
al., 2011; Steidel et al., 2011; Zheng et al., 2011; Dijkstra & Kramer, 2012;
Jeeson-Daniel et al., 2012; Verhamme et al., 2012; Kakiichi & Dijkstra, 2018;
Smith et al., 2018, 2019; Garel et al., 2021). 2) In-situ production:
Ly$\alpha$ photons are produced not inside the galaxy but in the CGM. This can
be further classified into three processes: i) recombination, ii) collisional
excitation, and iii) satellite galaxies. i) Recombination: ionizing radiation
from the galaxy or extragalactic background (UVB) photoionizes the hydrogen
gas in the CGM, which in turn emit Ly$\alpha$ emission via recombination
(‘fluorescence’; e.g., Furlanetto et al., 2005; Cantalupo et al., 2005;
Kollmeier et al., 2010; Lake et al., 2015; Mas-Ribas & Dijkstra, 2016; Gallego
et al., 2018; Mas-Ribas et al., 2017b). ii) Collisional excitation: Hydrogen
gas in the CGM is compressively heated by shocks and then emit Ly$\alpha$
photons by converting its gravitational energy into Ly$\alpha$ emission when
they accrete onto the galaxy (‘gravitational cooling’ or ‘cold stream’; e.g.,
Haiman et al., 2000; Fardal et al., 2001; Goerdt et al., 2010; Faucher-Giguère
et al., 2010; Rosdahl & Blaizot, 2012; Lake et al., 2015). iii) Satellite
galaxies: Ly$\alpha$ emission is produced by star formation (SF) in unresolved
dwarfs surrounding the galaxy (‘satellite galaxies’; e.g., Mas-Ribas et al.,
2017a, b). To recognize which process plays a major role, Kakuma et al. (2021)
compared the observed Ly$\alpha$ SB profiles against the prediction by Zheng
et al. (2011) that considers resonant scattering, although large uncertainties
of the data prevented to draw a conclusion. Since their comparison was limited
to the case for resonant scattering at $z=5.7$, we need to systematically
investigate various origins by utilizing multiple models at multiple redshifts
to pin down the key origin(s), in the similar way as Byrohl et al. (2020) and
Mitchell et al. (2021).

This paper is organized as follows. Our data are described in Section 2. In
Section 3, we use the intensity mapping technique to derive the cross-
correlation SB of Ly$\alpha$ emission around the LAEs. We discuss the redshift
evolution and physical origins of extended Ly$\alpha$ emission in Section 4
and summarize our findings in Section 5.

Throughout this paper, magnitudes are given in the AB system (Oke & Gunn,
1983). We adopt the concordance cosmology with $\Omega_{\rm m,0}=0.7$,
$\Omega_{\Lambda,0}=0.3$, and $H_{0}=70\ {\rm km\ s^{-1}\ Mpc^{-1}}$, where
$1\arcsec$ corresponds to transverse sizes of (8.3, 7.5, 5.9, 5.4) pkpc and
(26, 32, 39, 41) comoving kpc (ckpc) at $z=(2.2,\ 3.3,\ 5.7,\ 6.6)$.

##  2 Data

In this Section, we describe the images and sample catalogs used for our
analyses. All the images and catalogs are based on Subaru/HSC data.

###  2.1 Images

We use NB and broad-band (BB) imaging data that were obtained in two
Subaru/HSC surveys, the HSC-SSP and the CHORUS. The HSC-SSP and CHORUS data
were obtained in March 2014-January 2018 and January 2017-December 2018,
respectively. We specifically use the internal data of the S18A release. The
HSC-SSP survey is a combination of three layers: Wide, Deep, and UltraDeep
(UD). We use the UD layer images in the fields of the Cosmological Evolution
Survey (UD-COSMOS; Scoville et al., 2007) and Subaru/XMM Deep Survey (UD-SXDS;
Sekiguchi et al., 2005), because the wide survey areas ($\sim 2$ deg2 for each
field) and deep imaging (the $5\sigma$ limiting magnitudes are $\sim 26$ mag
in a $2\arcsec$-diameter aperture) in these fields are advantageous for the
detection of very diffuse Ly$\alpha$ emission. The CHORUS data were obtained
over the UD-COSMOS field. The HSC-SSP and CHORUS data were reduced with the
HSC pipeline v6.7 (Bosch et al., 2018).

The HSC-SSP program of S18A provides the data of two NB (NB816 and NB921)
filters in the UD-COSMOS and UD-SXDS fields, while the CHORUS images are
offered in four NB (NB387, NB527, NB718, and NB973) filters in the UD-COSMOS
field. In this work, we present the results in the NB387, NB527, NB816, and
NB921 filters. The NB718 and NB973 filters are not used in the following
sections, because the number of LAEs and the image depths are not sufficient
to detect diffuse Ly$\alpha$ emission. The NB387, NB527, NB816, and NB921
filters are centered at 3863, 5260, 8177, and 9215 Å with the full widths at
half maximum (FWHMs) of 55, 79, 113, and 135 Å, respectively, which cover the
observed wavelengths of Ly$\alpha$ emission from $z=2.178\pm 0.023$, $3.327\pm
0.032$, $5.726\pm 0.046$, and $6.580\pm 0.056$, respectively. Five BB ($g$-,
$r2$-, $i2$-, $z$-, and $y$-band) filters are also available in both HSC-SSP
and CHORUS. Figure 1 shows the NB and BB filter throughputs, and Table 1
summarizes the images and filters.

Figure 1: Effective response curves of the HSC NB filters: NB387 (blue), NB527
(green), NB816 (orange), and NB921 (red), and BB filters: $g$, $r2$, $i2$,
$z$, and $y$ (black). These response curves include the quantum efficiency of
the HSC CCD, transmittance of the dewar window and of the Primary Focus Unit,
reflectivity of the Primary Mirror, and airmass at the Telescope site. Table
1: Summary of the NB filters and images.

|  |  |  | UD-COSMOS |  | UD-SXDS
---|---|---|---|---|---|---
NB |  $\lambda_{\rm c}$ (Å) | FWHM (Å) | $z_{\rm LAE}$ | Area (deg2) |  $m_{{\rm NB},5\sigma}$ (mag) |  | Area (deg2) |  $m_{{\rm NB},5\sigma}$ (mag)
(1) | (2) | (3) | (4) | (5) | (6) |  | (7) | (8)
NB387 (CHORUS) | 3863 | 55 | $2.178\pm 0.023$ | 1.561 | 25.67 |  | — | —
NB527 (CHORUS) | 5260 | 79 | $3.327\pm 0.032$ | 1.613 | 26.39 |  | — | —
NB816 (HSC-SSP) | 8177 | 113 | $5.726\pm 0.046$ | 2.261 | 25.75 |  | 2.278 | 25.61
NB921 (HSC-SSP) | 9215 | 135 | $6.580\pm 0.056$ | 2.278 | 25.48 |  | 2.278 | 25.31

Note. — Columns: (1) NB filter. (2)-(3) Central wavelength ($\lambda_{\rm c}$)
and FWHM of the NB filter transmission curve. (4) Redshift range of the LAEs
whose Ly$\alpha$ emission enters the NB filter. (5)-(6) Effective area and
$5\sigma$ limiting NB magnitude ($m_{{\rm NB},5\sigma}$) in the UD-COSMOS
field. The $m_{{\rm NB},5\sigma}$ value is measured with a $2\arcsec$-diameter
aperture in each patch, and then averaged over the field. See Inoue et al.
(2020) and Hayashi et al. (2020) for the spatial variance of $m_{{\rm
NB},5\sigma}$. The $m_{{\rm NB},5\sigma}$ value in NB387 is corrected for the
systematic zero-point offset by 0.45 mag (see Section 2.2). (7)-(8) Same as
Columns (5)-(6), but for the UD-SXDS field. The values in this Table are cited
from Ono et al. (2021).

Bright sources in the NB and BB images must be masked since they contaminate
diffuse emission. We thus mask pixels flagged with either DETECT or
BRIGHT_OBJECT using the masks provided by the HSC pipeline (termed original
masks). A pixel is flagged with DETECT or BRIGHT_OBJECT when the pixel is
covered by a detected ($\geq 5\sigma$) object or is affected by nearby bright
sources, respectively. However, because a part of bright sources were missed
in the original masks due to bad photometry, the HSC-SSP team offered new
masks that mitigated this problem (hereafter termed revised
masks).111https://hsc-release.mtk.nao.ac.jp/doc/index.php/bright-star-masks-2/
We adopt the revised masks in addition to the original masks to flag
BRIGHT_OBJECT. We use the revised $g$-, $r2$-, $z$-, and $y$-band masks for
NB387, NB527, NB816, and NB921 images, respectively, because the revised masks
are offered only in the BB filters. For the BB images, we use the revised
masks defined for each BB filter. We visually confirm that these criteria
successfully cover bright sources and contaminants in the images.

###  2.2 LAE Samples

We use the LAE catalog constructed by Ono et al. (2021) as a part of the
Systematic Identification of LAEs for Visible Exploration and Reionization
Research Using Subaru HSC (SILVERRUSH) project (Ouchi et al. 2018, see also
Shibuya et al. 2018a, b; Konno et al. 2018; Harikane et al. 2018; Inoue et al.
2018; Higuchi et al. 2019; Harikane et al. 2019; Kakuma et al. 2021). Ono et
al. (2021) selected LAE candidates based on color and removed contaminants by
a convolutional neural network (CNN) and visual inspection. Their final
catalog includes (542, 959, 395, 150) LAEs at $z=(2.2$, 3.3, 5.7, 6.6) in the
UD-COSMOS field, and (560, 75) LAEs at $z=(5.7$, 6.6) in the UD-SXDS field.

The NB images of the UD-COSMOS and UD-SXDS fields are deepest at the center
and become shallower toward the edges (Hayashi et al., 2020; Inoue et al.,
2020). We thus exclude LAEs outside of the boundaries that are shown with the
black dashed circles in Figures 2 and 3.

We estimate the Ly$\alpha$ line luminosities ($L_{\rm Ly\alpha}$) following
Shibuya et al. (2018a, see also ). First, we measure the NB (BB) magnitudes
$m_{\rm NB}$ ($m_{\rm BB}$) of the LAEs at $z=2.2$, 3.3, 5.7, an 6.6 in the
NB387 ($g$-band), NB527 ($r2$-band), NB816 ($z$-band), and NB921 ($y$-band)
filters, respectively. The magnitudes are measured with a $2\arcsec$-diameter
aperture because it efficiently covers the point spread function (PSF), whose
FWHM is $0\farcs 8-1\farcs 1$ (Ono et al., 2021). The NB387 magnitudes are
corrected for the systematic zero-point offset by 0.45 mag, following the
recommendation by the HSC-SSP team.222https://hsc-
release.mtk.nao.ac.jp/doc/index.php/known-problems-2/#hsc-link-10 Next, we
follow Shibuya et al. (2018a) in our derivation of the Ly$\alpha$ line fluxes
($f_{\rm Ly\alpha}$) from $m_{\rm NB}$ and $m_{\rm BB}$, assuming a flat UV
continuum and the IGM attenuation model taken from Inoue et al. (2014).
Lastly, the values of $L_{\rm Ly\alpha}$ are derived via $L_{\rm
Ly\alpha}=4\pi d_{\rm L}(z_{\rm LAE})^{2}f_{\rm Ly\alpha}$, where $d_{\rm
L}(z_{\rm LAE})$ denotes the luminosity distance to the LAE at redshift
$z_{\rm LAE}$.

Although the completeness of the LAEs is as high as $\gtrsim 90$ % at $m_{\rm
NB}\lesssim 24.5$ in the CNN of Ono et al. (2021), faint LAEs may be missed in
the observations and selection. To ensure completeness, we use only LAEs whose
$L_{\rm Ly\alpha}$ values are larger than the modes (peaks) of the $L_{\rm
Ly\alpha}$ histograms, which are represented as $L_{\rm Ly\alpha}^{\rm min}$
in Table 2. This sample, termed the all sample, consists of (289, 762, 210,
56) LAEs at $z=(2.2$, 3.3, 5.7, 6.6) and (393, 24) LAEs at $z=(5.7$, 6.6) in
the UD-COSMOS and UD-SXDS fields, respectively. To accurately compare the LAEs
of different redshifts at similar $L_{\rm Ly\alpha}$ values, we further
exclude faint LAEs from the all sample such that the mean $L_{\rm Ly\alpha}$
values are equal to $10^{42.9}$ erg s-1 at each redshift. This selection
results in (37, 123, 125) LAEs at $z=(2.2$, 3.3, 5.7) and 313 LAEs at $z=5.7$
in the UD-COSMOS and UD-SXDS fields, respectively, which we hereafter refer to
as the bright subsample. At $z=6.6$, since the mean $L_{\rm Ly\alpha}$ values
of the all sample are $10^{43.0}$ erg s-1 in both the UD-COSMOS and UD-SXDS
fields, we also use the all sample as the bright subsample. In summary, we use
a total of 1781 and 717 LAEs in the UD-COSMOS$+$UD-SXDS fields as the all
sample and the bright subsample, respectively. The $L_{\rm Ly\alpha}$ values
and sample sizes are summarized in Table 2. The sky distributions of the LAEs
in the UD-COSMOS and UD-SXDS fields are presented in Figures 2 and 3,
respectively.

A major update of our data compared to Kakuma et al. (2021) is that we add new
1051 LAEs at $z=2.2$ and 3.3 using the CHORUS data. The catalog of $z=5.7$ and
6.6 LAEs are also updated in that we use the latest catalog constructed by Ono
et al. (2021) based on HSC-SSP S18A images, while Kakuma et al. (2021) used
the S16A catalog taken from Shibuya et al. (2018a). Although the number of the
LAEs increased owing to the improved limiting magnitudes, we excluded faint
LAEs from those included in Ono et al. (2021), which resulting in comparable
numbers of the LAEs between Kakuma et al. (2021)’s and our all samples. The
sky and $L_{\rm Ly\alpha}$ distributions are also similar between these two
samples. We use the same images as Kakuma et al. (2021, i.e., those taken from
the HSC-SSP S18A release) in NB816 and NB921, while our masking prescription
may be slightly different.

Table 2: Summary of the sample.

|  all sample |  |  bright subsample
---|---|---|---
| UD-COSMOS |  | UD-SXDS |  | UD-COSMOS |  | UD-SXDS
$z_{\rm LAE}$ | $L_{\rm Ly\alpha}^{\rm min}$ | $L_{\rm Ly\alpha}^{\rm mean}$ | $N_{\rm LAE}$ |  | $L_{\rm Ly\alpha}^{\rm min}$ | $L_{\rm Ly\alpha}^{\rm mean}$ | $N_{\rm LAE}$ |  | $L_{\rm Ly\alpha}^{\rm min}$ | $L_{\rm Ly\alpha}^{\rm mean}$ | $N_{\rm LAE}$ |  | $L_{\rm Ly\alpha}^{\rm min}$ | $L_{\rm Ly\alpha}^{\rm mean}$ | $N_{\rm LAE}$
(1) | (2) | (3) | (4) |  | (5) | (6) | (7) |  | (8) | (9) | (10) |  | (11) | (12) | (13)
2.2 | 42.2 | 42.5 | 289 |  | — | — | — |  | 42.6 | 42.9 | 37 |  | — | — | —
3.3 | 42.1 | 42.5 | 762 |  | — | — | — |  | 42.6 | 42.9 | 123 |  | — | — | —
5.7 | 42.6 | 42.8 | 210 |  | 42.5 | 42.8 | 393 |  | 42.7 | 42.9 | 125 |  | 42.6 | 42.9 | 313
6.6 | 42.8 | 43.0 | 56 |  | 42.9 | 43.0 | 24 |  | 42.8∗ | 43.0∗ | 56∗ |  | 42.9∗ | 43.0∗ | 24∗
Total |  |  | 1317 |  |  |  | 464 |  |  |  | 341 |  |  |  | 376

Note. — Columns: (1) Redshift. (2)-(3) Minimum and mean Ly$\alpha$ luminosity
of the all sample LAEs in the UD-COSMOS field, measured with a
$2\arcsec$-diameter aperture and shown in units of log erg s-1. (4) Number of
the LAEs. The total number of the LAEs over $z=2.2-6.6$ are shown in the
bottom row. (5)-(7) Same as Columns (2)-(4), but for the all subsample in the
UD-SXDS field. (8)-(13) Same as Columns (2)-(7), but for the bright subsample.

$*$ At $z=6.6$, we treat the all sample also as the bright subsample in each
field (see Section 2.2).

Figure 2: Sky distributions of the LAEs in the UD-COSMOS field at $z=2.2$,
3.3, 5.7, and 6.6, from left to right. The open red circles show the positions
of the LAEs included in the all sample but not in the bright subsample, while
the filled red circles indicate those of the bright subsample LAEs (i.e., all
sample LAEs are represented by the open$+$filled red circles). We use the LAEs
inside the black dashed circles. The background white shades show the NB387,
NB527, NB816, and NB921 images, from left to right. The gray shades show the
regions where the pixel is masked or an UD image is not offered. Note that the
distribution of the bright subsample at $z=6.6$ is not displayed, because we
treat the all sample also as the bright subsample at $z=6.6$. Figure 3: Same
as Figure 2, but for the UD-SXDS field at $z=5.7$ (left) and 6.6 (right).

###  2.3 NonLAE Samples

Although the intensity mapping technique can remove spurious signals from low-
redshift interlopers, other systematics, such as the sky background and PSF,
may still contaminate LAE signals. To estimate the contribution from these
systematics, we use foreground sources, which we hereafter term “NonLAEs.”
Because NonLAEs should correlate only with the systematics, but not with
Ly$\alpha$ emission from LAEs, we estimate these systematics by applying the
intensity mapping technique to NonLAEs.

We construct NonLAE samples as follows. First, we detect sources in the NB
images using SExtractor (Bertin & Arnouts, 1996). Second, we select only
sources that are sufficiently bright ($\lesssim 26$ mag) in the $g$-, $r2$-,
and $i2$-bands, to remove spurious sources and artifacts. Third, we randomly
select the sources such that they have the same sky, FWHM, and $m_{\rm NB}$
distributions as those of the LAEs in each field at each redshift. In this
way, $\sim 10^{3}$ sources are selected, which we define as NonLAEs (see
Figure 4 for the FWHM-$m_{\rm NB}$ distributions of the $z=3.3$ LAEs and
corresponding NonLAEs). We note that only $<1$ % of the NonLAEs meet the color
selection criteria of LAEs defined by Ono et al. (2021, see their Section 2).

Figure 4: Center: NB FWHM (${\rm FWHM_{NB}}$) and magnitude ($m_{\rm NB}$)
distributions of the all sample at $z=3.3$ (red circles), and the
corresponding NonLAEs (black dots). Top: $m_{\rm NB}$ histogram of the LAEs
and NonLAEs (red and black bars, respectively). Right: Same as the top panel,
but for ${\rm FWHM_{NB}}$.

##  3 Intensity Mapping Analysis

In this section, we derive the SB radial profiles of Ly$\alpha$ emission
around the LAEs using the intensity mapping technique.

###  3.1 Cross-correlation Functions

We compute SB as a cross-correlation function between given band (XB) emission
intensities and given objects (OBJs), ${\rm SB_{XB\times OBJ}}$333The right-
hand side of Equation 1 represents a cross-correlation function, which has
usually been referred to as $\xi$ in previous work (e.g. Croft et al., 2016,
2018; Bielby et al., 2017; Momose et al., 2021b, a). However, we refer to this
function as ‘SB,’ since the cross-correlation function is equivalent to
surface brightness in our analyses., via

${\rm SB_{XB\times OBJ,\nu}}(r)=\frac{1}{N_{r,{\rm OBJ}}}\sum_{i=1}^{N_{r,{\rm
OBJ}}}\mu^{(\rm XB)}_{\nu,i}$ (1)

and

${\rm SB_{XB\times OBJ}}(r)={\rm SB_{XB\times OBJ,\nu}}(r)\times{\rm
FWHM_{XB}}.$ (2)

The pixel of the $i$-th pixel-OBJ pair has a pixel value of $\mu^{(\rm
XB)}_{\nu,i}$ in the XB image in units of erg s-1 cm-2 Hz-1 arcsec-2. In
NB387, we multiply $\mu^{(\rm NB387)}_{\nu,i}$ by 1.5 to correct fot the zero-
point offset of 0.45 mag (see Section 2.2). The summation runs over the
$N_{r,{\rm OBJ}}$ pixel-OBJ pairs that are separated by a spatial distance
$r$. We use pixels at distances of between $1\farcs 5$ and $40\arcsec$ from
each OBJ (corresponding to the outer part of the CGM and outside), which were
then divided into six radial bins. A total of $\sim(1-2)\times 10^{9}$ pixels
were used for the calculation at each redshift. ${\rm FWHM_{XB}}$ represents
the XB filter width (in units of Hz) corrected for IGM attenuation, derived
via

${\rm FWHM_{XB}}=\frac{\int_{0}^{\infty}e^{-\tau_{\rm eff}(\nu)}T_{\rm
XB}(\nu){\rm d}\nu/\nu}{T_{\rm XB}(\nu_{\alpha})/\nu_{\alpha}},$ (3)

where $T_{\rm XB}(\nu)$ denotes the transmittance of the XB filter, and
$\nu_{\alpha}$ is the the observed Ly$\alpha$ line frequency. We adopt the IGM
optical depth $\tau_{\rm eff}(\nu)$ from Inoue et al. (2014).

The statistical uncertainty of ${\rm SB_{XB\times OBJ}}$ is estimated by the
bootstrap method. We randomly resample objects while keeping the sample size
and calculated ${\rm SB_{XB\times OBJ}}$. We then repeat the resampling
$10^{4}$ times, adopting the $1\sigma$ standard deviation of the ${\rm
SB_{XB\times OBJ}}$ values as the $1\sigma$ statistical uncertainty of the
original ${\rm SB_{XB\times OBJ}}$.

###  3.2 Ly$\alpha$ Surface Brightness

We estimate the SB of Ly$\alpha$ emission (${\rm SB_{Ly\alpha}}$) around the
LAEs as follows. First, we subtract the systematics (${\rm SB_{NB\times
NonLAE}}$) from the emission from the LAEs (${\rm SB_{NB\times LAE}}$) via

${\rm SB_{NB}=SB_{NB\times LAE}-SB_{NB\times NonLAE}}.$ (4)

Uncertainties in ${\rm SB_{NB}}$ propagate from those of ${\rm SB_{NB\times
LAE}}$ and ${\rm SB_{NB\times NonLAE}}$. We present the radial profiles of the
${\rm SB_{NB\times LAE}}$, ${\rm SB_{NB\times NonLAE}}$, and ${\rm SB_{NB}}$
of the all sample and bright subsample in Figures 5 and 6, respectively.

We use Fisher’s method (Fisher, 1970) to estimate the S/N ratios of ${\rm
SB_{NB}}$ over all the radial bins, following Kakuma et al. (2021). Generally,
a $p$-value is expressed as

$p=\int_{{\rm S/N}}^{\infty}\mathcal{N}(x;\mu=0,\sigma=1){\rm d}x,$ (5)

where $\mathcal{N}(x;\mu=0,\sigma=1)$ is a Gaussian distribution with an
expected value $\mu=0$ and a variance $\sigma^{2}=1$. We thus use this
equation to convert the S/N ratio in the $i$-th radial bin (S/Ni) into the
$p$-value in that bin ($p_{i}$). The $\chi^{2}$ value over all the radial bins
($1\leq i\leq N$), $\hat{\chi}^{2}$, is then calculated as
$\hat{\chi}^{2}=-2\sum_{i=1}^{N}\ln(p_{i})$. Since $\hat{\chi}^{2}$ follows a
$\chi^{2}_{2N}$ distribution with $2N$ degrees of freedom, $\chi^{2}(x;{\rm
dof}=2N)$, the $p$-value over all the radial bins, $\hat{p}$, is derived as

$\hat{p}=\int_{\hat{\chi}^{2}}^{\infty}\chi^{2}(x;{\rm dof}=2N){\rm d}x.$ (6)

We convert this to the S/N ratio over all the radial bins by solving Equation
(5) for S/N. We use the radial bins at $<1$ cMpc.

The black vertical dashed lines in Figures 5 and 6 indicate $R_{\rm vir}$.
Following the observational results by Ouchi et al. (2010) and Kusakabe et al.
(2018), we assume that the DMHs hosting LAEs have halo masses ($M_{\rm halo}$)
of $10^{11}\ M_{\odot}$ at all the redshifts, which corresponds to a $R_{\rm
vir}$ value of $\sim 150$ ckpc.

Figure 5: SB radial profiles of the all sample. Left: ${\rm SB_{NB}}$ in the
UD-COSMOS field at $z=2.2$, 3.3, 5.7, and 6.6, from top to bottom. The red
triangles and circles show SB before and after systematics subtraction (${\rm
SB_{NB\times LAE}}$ and ${\rm SB_{NB}}$), respectively. The black solid lines
represent the systematics that are estimated with the NonLAEs (${\rm
SB_{NB\times NonLAE}}$). The red error bars and gray shades are the $1\sigma$
uncertainties estimated by the bootstrap method. The vertical black dashed
line represents $R_{\rm vir}$ of a DMH with $M_{\rm halo}=10^{11}$
$M_{\odot}$, while the horizontal black line represents ${\rm SB}=0$. The
upper, middle, and lower horizontal axes in each panel show radii in units of
arcsec, pkpc, and ckpc, respectively. The data points of ${\rm SB_{NB\times
LAE}}$ are slightly shifted along the horizontal axis for clarity. Right: Same
as the left column, but in the UD-SXDS field at $z=5.7$ and 6.6 in the top and
bottom panels, respectively. The S/N ratios of ${\rm SB_{NB}}$ are 4.1, 1.6,
and 2.2 at $z=3.3$, 5.7, and 6.6, respectively (they are based on ${\rm
SB_{NB}}$ averaged over the UD-COSMOS and UD-SXDS fields, at $z=5.7$ and 6.6).
There is no clear detection at $z=2.2$. Figure 6: Same as Figure 5, but for
the bright subsample. The S/N ratios of ${\rm SB_{NB}}$ are 1.5 and 1.4 at
$z=2.2$ and 5.7, respectively, while there is no clear detection at $z=3.3$.
The results at $z=6.6$ (bottom panels) are the same as those in Figure 5 since
we treat the all sample also as the bright subsample at $z=6.6$.

As presented in Figures 5 and 6, at $z=2.2$, we see no clear detection for the
all sample, although ${\rm SB_{NB\times LAE}}$ slightly exceeds ${\rm
SB_{NB\times NonLAE}}$ for the bright subsample with ${\rm S/N}=1.5$. At
$z=3.3$, ${\rm SB_{NB}}$ of the all sample is significantly positive over wide
scales from $\sim 100$ ckpc to 1 cMpc with ${\rm S/N}=4.1$. At $z=5.7$, ${\rm
SB_{NB}}$ is significantly positive at $\sim 80-10^{3}$ ckpc for the all
sample and bright subsample in the UD-COSMOS field with ${\rm S/N}=4.6$ and
4.1, respectively. Averaging ${\rm SB_{NB}}$ of the all sample in both fields
with weights of the number of the LAEs in each field, we identify a
tentatively positive signal with ${\rm S/N}=1.6$. At $z=6.6$, ${\rm SB_{NB}}$
are positive at $\sim 80-200$ ckpc in both fields. Averaging ${\rm SB_{NB}}$
over the two fields, we tentatively identify a positive signal with ${\rm
S/N}=2.3$. We find that the emission SB at $z=3.3-6.6$ is as diffuse as $\sim
10^{-20}-10^{-19}$ erg s-1 cm-2 arcsec-2. We note that ${\rm SB_{NB\times
LAE}}$ in the UD-SXDS field are systematically lower than ${\rm SB_{NB\times
NonLAE}}$, as reported in Kakuma et al. (2021).

UV continuum emission might contribute to the ${\rm SB_{NB}}$ in addition to
the Ly$\alpha$ line emission. We thus estimate SB of the UV continuum
emission, ${\rm SB_{cont,\nu}}$, using

${\rm SB_{cont,\nu}<SB_{BB,\nu}\equiv SB_{BB\times LAE,\nu}-SB_{BB\times
random,\nu}},$ (7)

where ${\rm SB_{BB\times LAE,\nu}}$ (${\rm SB_{BB\times random,\nu}}$)
represents the SB value that is derived from the cross-correlation function
between the BB images and LAEs (random sources) in units of erg s-1 cm-2 Hz-1
arcsec-2. To estimate the sky background, we use random sources, not NonLAEs,
since it is difficult to match the ${\rm FWHM_{NB}}$-$m_{\rm NB}$
distributions of NonLAEs with those of the LAEs due to faintness of the LAEs
in the BB. Since ${\rm SB_{BB\times random,\nu}}$ neglects neglects signals
from the PSF, ${\rm SB_{BB,\nu}}$ should be treated as the upper limit of
${\rm SB_{cont,\nu}}$. As shown in Figure 7 as an example, the values of ${\rm
SB_{cont,\nu}}$ at $z=2.2-6.6$ are consistent with null detection within
$\sim(1-2)\sigma$ uncertainties. Additionally, we confirm that the UV
continuum emission contributing to ${\rm SB_{NB}}$, i.e., ${\rm
SB_{cont,\nu}}\times{\rm FWHM_{NB}}$, is negligible compared to ${\rm
SB_{NB}}$. Therefore, we hereafter assume that ${\rm SB_{NB}}$ is equivalent
to the Ly$\alpha$ SB (${\rm SB_{Ly\alpha}}$).

Figure 7: Radial profiles of the UV continuum SB of the all sample (open) and
bright subsample (filled) at $z=3.3$. The profiles are averages weighted with
the number of the LAEs in the UD-COSMOS and UD-SXDS fields. The data points
for bright subsample are slightly shifted along the horizontal axis for
clarity. The data points should be treated as upper limits, since ${\rm
SB_{cont,\nu}}$ neglects signals from the PSF. ${\rm SB_{cont,\nu}}$ agree
with null detection within $1\sigma$ uncertainties for both sample, which
disfavors the UV continuum contribution to ${\rm SB_{NB}}$.

In summary, we identify very diffuse ($\sim 10^{-20}-10^{-19}$ erg s-1 cm-2
arcsec-2) Ly$\alpha$ signals beyond $R_{\rm vir}$ around the all sample LAEs
at $z=3.3$ at the $4.1\sigma$ level. We also potentially detect positive
signals around the all sample LAEs at $z=5.7$ and 6.6 and bright subsample
LAEs with ${\rm S/N}\sim 2$. These results imply the potential existence of
very diffuse and extended Ly$\alpha$ emission around $z=2.2-6.6$ LAEs.

Again, a major update compared to Kakuma et al. (2021) is that we apply the
intensity mapping technique newly to the CHORUS data to investigate extended
Ly$\alpha$ emission at $z=2.2-3.3$. In particular, we identify extended
Ly$\alpha$ emission around the $z=3.3$ LAEs with the S/N levels comparable or
higher than those at $z=5.7$ and 6.6, probably thanks to the large number of
LAEs at $z=3.3$ ($N_{\rm LAE}=762$). Our results at $z=5.7$ and 6.6 are
consistent with those obtained in Kakuma et al. (2021), which is as expected
from similar properties of the LAEs of Kakuma et al. (2021) and our all
sample, such as $N_{\rm LAE}$, the $L_{\rm Lya}$ ranges, and the sky
distributions (Section 2.2). In the next section, we compare our results with
previous work including Kakuma et al. (2021) in detail.

###  3.3 Comparison with Previous Work

We compare our ${\rm SB_{Ly\alpha}}$ radial profiles with those of previous
studies at $z\sim 2.2$, 3.3, 5.7, and 6.6 (Figure 8). We compile the data
taken from Momose et al. (2014), Momose et al. (2016), Leclercq et al. (2017),
Wisotzki et al. (2018), Wu et al. (2020), and Kakuma et al. (2021), which are
summarized in Table 3. Because SB is affected by the cosmological dimming
effect, all the ${\rm SB_{Ly\alpha}}$ profiles, including ours, are scaled by
$(1+z)^{-4}$ to $z=2.2$, 3.3, 5.7, and 6.6 in each panel. We also shift the
radii of the ${\rm SB_{Ly\alpha}}$ profiles in units of cpkc by $(1+z)$, while
fixing the radii in pkpc. For our samples at $z=5.7$ and 6.6, we hereafter
present ${\rm SB_{Ly\alpha}}$ averaged at each redshift over the UD-COSMOS and
UD-SXDS fields weighting by the number of LAEs in each field, unless otherwise
stated.

Although the ${\rm SB_{Ly\alpha}}$ profiles are measured under different
seeing sizes, the typical image PSF FWHMs are as small as $\lesssim 1\farcs
5$, corresponding to $\lesssim 40-60$ ckpc at $z=2-7$ (e.g. Momose et al.,
2014; Ono et al., 2021). Since we focus on ${\rm SB_{Ly\alpha}}$ profiles at
larger scales of $\gtrsim 100$ ckpc, PSF differences are unlikely to affect
the following discussion.

Momose et al. (2016) found that ${\rm SB_{Ly\alpha}}$ profiles depend on
$L_{\rm Ly\alpha}$ of the galaxy. To avoid this dependency, we take the data
from Momose et al. (2014), Momose et al. (2016), and Kakuma et al. (2021),
because their LAEs have $L_{\rm Ly\alpha}$ values similar to those of our all
samples in the same $2\arcsec$-diameter aperture size. We additionally take
the data from Leclercq et al. (2017), Wisotzki et al. (2018), and Wu et al.
(2020), but these samples have different $L_{\rm Ly\alpha}$ values measured in
different aperture sizes. Therefore, for precise comparisons, we normalize the
${\rm SB_{Ly\alpha}}$ profiles of these samples such that ${\rm
SB_{Ly\alpha}}$ integrated over a central $2\arcsec$-diameter aperture ($=4\pi
d_{\rm L}^{2}\int_{0^{\prime\prime}}^{1^{\prime\prime}}{\rm
SB_{Ly\alpha}}(r)\cdot 2\pi r{\rm d}r$) becomes equal to $L_{\rm Ly\alpha}$ of
our all sample at each redshift.

Table 3: Summary of the observational studies used for comparison.

Reference | $z$ | $L_{\rm Ly\alpha}$ | Sample | Method
---|---|---|---|---
$z\sim 2.2$ |  |  |  |
This Work (all sample) | 2.2 | 42.5 | 289 LAEs | intensity mapping (mean)
Momose et al. (2016) | 2.2 | 42.6 | 710 LAEs ($L_{\rm Ly\alpha}\geq 10^{42.4}$ erg s-1) | stacking (mean)
$z\sim 3.3$ |  |  |  |
This Work (all sample) | 3.3 | 42.5 | 762 LAEs | intensity mapping (mean)
Momose et al. (2014) | 3.1 | 42.7 | 316 LAEs | stacking (mean)
Leclercq et al. (2017) | 3.28 | [42.5] | LAE MUSE#106∗ | individual detection
Wisotzki et al. (2018) | $3-4$ | [42.5] | 18 LAEs ($L_{\rm Ly\alpha}>10^{42}$ erg s-1) | stacking (median)
Matsuda et al. (2012) | 3.1 | — | 894 LAEs ($26<BV<27$)† | stacking (median)
$z\sim 5.7$ |  |  |  |
This Work (all sample) | 5.7 | 42.8 | 650 LAEs | intensity mapping (mean)
Momose et al. (2014) | 5.7 | 42.7 | 397 LAEs | stacking (mean)
Leclercq et al. (2017) | 5.98 | [42.8] | LAE MUSE#547∗ | individual detection
Wisotzki et al. (2018) | $5-6$ | [42.8] | 6 LAEs ($L_{\rm Ly\alpha}>10^{42}$ erg s-1) | stacking (median)
Wu et al. (2020) | 5.7 | [42.8] | 310 LAEs | stacking (median)
Kakuma et al. (2021) | 5.7 | 42.9 | 425 LAEs | intensity mapping (mean)
$z\sim 6.6$ |  |  |  |
This Work (all sample) | 6.6 | 43.0 | 80 LAEs | intensity mapping (mean)
Momose et al. (2014) | 6.6 | 42.7 | 119 LAEs | stacking (mean)
Kakuma et al. (2021) | 6.6 | 42.8 | 396 LAEs | intensity mapping (mean)

Note. — Columns: (1) Reference. (2) Redshift. (3) Mean or median Ly$\alpha$
luminosity of the sample within a $2\arcsec$-diameter aperture in units of log
erg s-1. We normalize the ${\rm SB_{Ly\alpha}}$ profiles of Leclercq et al.
(2017), Wisotzki et al. (2018), and Wu et al. (2020) such that ${\rm
SB_{Ly\alpha}}$ integrated over a central $2\arcsec$-diameter aperture becomes
equal to $L_{\rm Ly\alpha}$ of our all sample at each redshift, which are
indicated by the brackets. (4) Sample used for comparison. The parentheses
indicate specific subsamples. (5) Method for deriving the ${\rm
SB_{Ly\alpha}}$ profiles.

∗ ID of the individual LAE. See Section 3.3 for the variance among the
individual LAEs.

† $BV\equiv(2B+V)/3$, where $B$ and $V$ are $B$\- and $V$-band magnitudes.

Figure 8: Comparison of the ${\rm SB_{Ly\alpha}}$ radial profiles at $z=2.2$
(top left), 3.3 (top right), 5.7 (bottom left), and 6.6 (bottom right). All
the ${\rm SB_{Ly\alpha}}$ profiles are corrected by $(1+z)^{-4}$ to match to
each redshift. The gray shades illustrate the regions whose radius is smaller
than $R_{\rm vir}$ of a DMH with $M_{\rm halo}=10^{11}\ M_{\odot}$. The filled
red circles show the ${\rm SB_{Ly\alpha}}$ profiles of our all sample. The
profiles at $z=5.7$ and 6.6 are averaged over the UD-COSMOS and UD-SXDS
fields. The other red symbols represent the profiles taken from previous
studies: Matsuda et al. (2012, Ma12: crosses), Momose et al. (2014, Mo14:
squares), Momose et al. (2016, Mo16: squares), Leclercq et al. (2017, Le17:
diamonds), Wisotzki et al. (2018, Wi18: pentagons), Wu et al. (2020, Wu20:
hexagons), and Kakuma et al. (2021, Ka21: stars). See Table 3 for details of
these samples. We omit the data points below the detection limits defined in
the literature. The ${\rm SB_{Ly\alpha}}$ profiles taken from Leclercq et al.
(2017), Wisotzki et al. (2018), and Wu et al. (2020) are normalized such that
$L_{\rm Ly\alpha}$ values measured in a central $2\arcsec$ diameter equal to
those of our all sample at each redshift. Some data points are slightly
shifted along the horizontal axes for clarity. The ${\rm SB_{Ly\alpha}}$
profiles of our all sample are in good agreements with those taken from the
literature.

The top left panel in Figure 8 presents ${\rm SB_{Ly\alpha}}$ radial profiles
at $z=2.2$. We compare the results of Momose et al. (2016, $L_{\rm
Ly\alpha}=10^{42.6}$ erg s-1 subsample) against our all sample. We found that
the ${\rm SB_{Ly\alpha}}$ profile of Momose et al. (2016) is in good agreement
with that of our all sample at $r\sim 100$ ckpc.

In the top right panel of Figure 8, we show the ${\rm SB_{Ly\alpha}}$ profiles
at $z=3.3$. We compare the results of Momose et al. (2014, $z=3.1$ LAEs),
Leclercq et al. (2017, an individual LAE MUSE#106), and Wisotzki et al. (2018,
$L_{\rm Ly\alpha}>10^{42}$ erg s-1 subsample at $z=3-4$), and our all sample.
The ${\rm SB_{Ly\alpha}}$ profiles from the literature approximately agree
with that of our all sample even at $r\sim R_{\rm vir}$. We additionally show
the results from a subsample of Matsuda et al. (2012) with continuum
magnitudes $BV\equiv(2B+V)/3$ of $26<BV<27$ (typical value range for LAEs),
where $B$ and $V$ are $B$\- and $V$-band magnitudes measured with the
Subaru/SC. The ${\rm SB_{Ly\alpha}}$ profile of Matsuda et al. (2012) also
agrees with that of our all sample around $r\sim R_{\rm vir}$. We note that,
although the results of Leclercq et al. (2017) are represented by their
individual LAE MUSE#6905, their LAEs have similar ${\rm SB_{Ly\alpha}}$
profiles when the amplitudes are normalized to match the $L_{\rm Ly\alpha}$
values at $r\leq 1\arcsec$.

The ${\rm SB_{Ly\alpha}}$ profiles at $z=5.7$ are displayed in the bottom left
panel of Figure 8. We compare the results of Momose et al. (2014, $z=5.7$
LAEs), Leclercq et al. (2017, an individual LAE MUSE#547), Wisotzki et al.
(2018, $L_{\rm Ly\alpha}>10^{42}$ erg s-1 subsample at $z=5-6$), Wu et al.
(2020, $z=5.7$ LAEs), Kakuma et al. (2021, $z=5.7$ LAEs), and our all sample.
The ${\rm SB_{Ly\alpha}}$ profile of our all sample agrees well with those of
Momose et al. (2014), Leclercq et al. (2017), Wisotzki et al. (2018) and Wu et
al. (2020) at $r\sim 80-200$ ckpc, and with that of Kakuma et al. (2021) up to
$r\sim 1$ cMpc.

The bottom right panel of Figure 8 shows ${\rm SB_{Ly\alpha}}$ profiles at
$z=6.6$ taken from Momose et al. (2014, $z=6.6$ LAEs), Kakuma et al. (2021,
$z=6.6$ LAEs), and our all sample. The ${\rm SB_{Ly\alpha}}$ profile of our
all sample is consistent with those of Momose et al. (2014) and Kakuma et al.
(2021) up to the scales of $r\sim 100$ ckpc and 1 cMpc, respectively.

In summary, our ${\rm SB_{Ly\alpha}}$ profiles are in good agreement with
those of the previous studies at each redshift, provided that the LAEs have
similar $L_{\rm Ly\alpha}$ values at $r\leq 1\arcsec$. Our ${\rm
SB_{Ly\alpha}}$ profiles ($r\gtrsim 80$ ckpc) are smoothly connected with the
inner ($r\lesssim 100$ ckpc) profiles taken from the literature at $\sim 100$
ckpc, and extend to larger scales.

In Figure 9, we compare the ${\rm SB_{cont,\nu}}$ radial profiles between
Momose et al. (2014, 2016) and our all sample at $z=2.2$, 3.3, and 5.7 (the
data of Momose et al. 2014, 2016 are of the same LAEs as used in Figure 8). We
find that the ${\rm SB_{cont,\nu}}$ profiles roughly agree at $r\lesssim 400$
ckpc, although the uncertainties are large. The ${\rm SB_{cont,\nu}}$ profiles
are much less extended than ${\rm SB_{Ly\alpha}}$ profiles, which was also
suggested by Momose et al. (2014), Momose et al. (2016), and Wu et al. (2020).
We note the profiles at $z=6.6$ are not compared here because our sample is
$\sim 0.3$ dex brighter than that of Momose et al. (2014).

Figure 9: Comparison of the ${\rm SB_{cont,\nu}}$ radial profiles at $z=2.2$
(top), 3.3 (middle), and 5.7 (bottom). The red circles and cyan diamonds show
the ${\rm SB_{cont,\nu}}$ profiles of our all sample and Momose et al. (2014,
2016), respectively (we use the data of Momose et al. 2016 at $z=2.2$, and of
$z=3.1$ and 5.7 LAEs of Momose et al. 2014 at $z=3.3$ and 5.7, respectively).
The gray shades illustrate the regions inside $R_{\rm vir}$ of a DMH with
$M_{\rm halo}=10^{11}\ M_{\odot}$. The ${\rm SB_{cont,\nu}}$ profiles of our
all sample are in good agreements with those of Momose et al. (2014, 2016).

##  4 Discussion

###  4.1 Redshift Evolution of Extended Ly$\alpha$ Emission Profiles

In this section, we investigate the redshift evolution of ${\rm
SB_{Ly\alpha}}$ profiles of extended Ly$\alpha$ emission. In Figure 10, we
compare the ${\rm SB_{Ly\alpha}}$ profiles of our bright subsamples at
$z=2.2-6.6$ as a function of radius in units of ckpc. We also present the
results taken from Momose et al. (2014, $z=5.7$ and 6.6 LAEs) and Wisotzki et
al. (2018, $L_{\rm Ly\alpha}>10^{42}$ erg s-1 subsample at $z=3-4$). Because
our bright subsamples have uniform $L_{\rm Ly\alpha}$ values ($\sim
10^{42.9}-10^{43.0}$ erg s-1) over $z=2.2-6.6$, the $L_{\rm Ly\alpha}$
differences between the redshifts are unlikely to influence the following
discussion.

Figure 10: Comparison between ${\rm SB_{Ly\alpha}}$ radial profiles at
$z=2.2-6.6$. Top: Observed ${\rm SB_{Ly\alpha}}$ profiles. The blue hexagons,
green pentagons, orange squares, and red circles represent the profiles of our
bright subsamples at $z=2.2$, 3.3, 5.7, and 6.6, respectively. The green,
orange, and red crosses denote profiles taken from Wisotzki et al. (2018),
Momose et al. (2014), and Momose et al. (2014), which are normalized to
$z=3.3$, 5.7, and 6.6, respectively. The ${\rm SB_{Ly\alpha}}$ profile of
Wisotzki et al. (2018) is normalized such that $L_{\rm Lya}$ in $r\leq
1\arcsec$ matches to that of our all sample (see Section 3.3). We again note
that our profiles here are based on bright subsamples and thus differ from
those shown in Figure 8. Bottom: Same as the top panel, but showing the
profiles corrected for the cosmological dimming effect, i.e., the intrinsic
${\rm SB_{Ly\alpha}}$ profiles. We tentatively identify an increasing trend
roughly by $(1+z)^{3}$ toward high redshifts, albeit with the large
uncertainties.

The top panel of Figure 10 shows the observed ${\rm SB_{Ly\alpha}}$ profiles.
We identify no significant difference among the ${\rm SB_{Ly\alpha}}$ profiles
beyond the $1\sigma$ uncertainties at $r\sim 100-1000$ ckpc over $z=2.2-6.6$,
while the uncertainties are large. This finding is consistent with that of
Kakuma et al. (2021) at $z=5.5-6.6$. There is no significant difference also
in the profiles at $r<100$ ckpc, which was also suggested by MUSE observations
(Leclercq et al., 2017, see also Figure 11 of Byrohl et al. 2020) at $3<z<6$.

Observed ${\rm SB_{Ly\alpha}}$ profiles are affected by the cosmological
dimming effect. To correct this effect, we shift the observed ${\rm
SB_{Ly\alpha}}$ profiles vertically by $(1+z)^{4}/(1+3.3)^{4}$ and
horizontally by $(1+3.3)/(1+z)$, which are hereafter termed as the intrinsic
profiles (we match the profiles to $z=3.3$ just for visibility). The intrinsic
profiles are presented in the bottom panel of Figure 10. There is a tentative
increasing trend in the intrinsic profiles toward $z=6.6$, although those at
$z=2.2$ and 3.3 remain comparable due to the large uncertainties.

To quantitatively investigate the evolution, we derive the ${\rm
SB_{Ly\alpha}}$ intrinsic profiles amplitudes at $r=200$ ckpc (${\rm
SB_{Ly\alpha,intr}}$) as a function of redshift. We fit the relation with
${\rm SB_{Ly\alpha,intr}}\propto(1+z)^{b}$, where $b$ is a constant, weighting
${\rm SB_{Ly\alpha,intr}}$ with $\sqrt{N_{\rm LAE}}$ at each redshift. We find
that the best-fit value of $b$ is 3.1, which implies that the intrinsic ${\rm
SB_{Ly\alpha}}$ profile amplitudes increase toward high redshifts roughly by
$(1+z)^{3}$ at a given radius in units of ckpc. This trend might correspond to
increasing density of hydrogen gas toward high redshifts due to the evolution
of the cosmic volume. Nevertheless, it is still difficult to draw a conclusion
due to large uncertainties. We cannot rule out other or additional
possibilities, such as higher Ly$\alpha$ escape fractions toward high
redshifts (e.g., Hayes et al., 2011; Konno et al., 2016). It is also necessary
to investigate potential impact of the cosmic reionization on the neutral
hydrogen density at $r\sim 100-1000$ ckpc around LAEs. Deeper observations in
the future will help to elucidate the evolution and its physical
interpretation.

###  4.2 Physical Origins of Extended Ly$\alpha$ Emission

In this section, we compare the observational results obtained in Section 3
with theoretical models taken from the literature to investigate the mechanism
of extended Ly$\alpha$ emission production. We investigate extended Ly$\alpha$
emission focusing on 1) where Ly$\alpha$ photons originate; 2) which processes
produce Ly$\alpha$ photons; and 3) how these photons transfer in the
surrounding materials.

1) First, we distinguish Ly$\alpha$ emission according to where it originates
from:

  1. 1.

The ISM of the targeted galaxy (central galaxy).

  2. 2.

The CGM surrounding the central galaxy.

  3. 3.

Satellite galaxies.

  4. 4.

Other halos, which refer to halos distinct from that hosting the central
galaxy.

2) We consider that Ly$\alpha$ photons are produced in the processes of
recombination and/or collisional excitation (cooling radiation), as stated in
the Introduction.

3) Lastly, Ly$\alpha$ photons subsequently transfer through surrounding
hydrogen gas while being resonantly scattered, which affects the observed SB
profiles. Hence, we should be conscious of whether models takes into account
the scattering process or not (see Byrohl et al. 2020 for the impact of
scattering on ${\rm SB_{Ly\alpha}}$ profiles).

Table 4: Summary of the theoretical studies used for comparison.

Reference | $z_{\rm model}$ | $z_{\rm plot}$ | Origin | Process | Scattering | $M_{\rm halo}\ (M_{\odot})$
---|---|---|---|---|---|---
Mas-Ribas et al. (2017b) | 5.7, 6.6 | 2.2, 3.3, 5.7, 6.6 | CGM/Sat. | Rec. | n | $10^{11\ *}$
Kakiichi & Dijkstra (2018) | $2-3$ | 2.2 | Cen. | Rec. | y | $\sim 10^{12}$
Byrohl et al. (2020) | $2-5$ | 3.3 | Cen./Sat./Other. | Rec./Cool. | y | $\sim 10^{11\ \\#}$
Dijkstra & Kramer (2012) | $\sim 2.65$ | 3.3 | Cen. | Rec. | y | —
Lake et al. (2015) | 3.1 | 3.3 | Cen./Other. | Rec./Cool. | y | $10^{11.5}$
Zheng et al. (2011) | 5.7 | 5.7 | Cen./Other. | Rec. | y | $10^{11.2}$

Note. — Columns: (1) Reference. (2) Redshift assumed in the model. (3)
Redshift where we plot the model in Figure 11. (4) Origins of Ly$\alpha$
photons (where Ly$\alpha$ photons are produced): the central galaxy (Cen.),
CGM, satellite galaxies (Sat.), and other halos (Other.) (5) Physical
processes of Ly$\alpha$ emission production (how Ly$\alpha$ emission is
produced): recombination (Rec.), and cooling radiation (Cool.). (6) Whether
the model takes into account Ly$\alpha$ RT (i.e., resonant scattering) or not.
(7) Assumed halo mass in units of $M_{\odot}$.

∗ We fix $M_{\rm halo}$ to $10^{11}\ M_{\odot}$, while $M_{\rm halo}$ is a
free parameter in the model.

# This value was converted from the assumed stellar mass $M_{\star}$ of
$10^{8.5}-10^{9.5}\ M_{\odot}$ with the $M_{\star}/M_{\rm halo}$ ratio
obtained in Behroozi et al. (2019).

Figure 11: Comparison of the ${\rm SB_{Ly\alpha}}$ radial profiles between the
observational and theoretical studies at $z=2.2$ (top left), 3.3 (top right),
5.7 (bottom left), and 6.6 (bottom right). The filled red circles depict the
observational results taken from our all sample ($z=2.2-6.6$), while the open
red circles show the results from previous observational work (Momose et al.
2016 at $z=2.2$, Wisotzki et al. 2018 at $z=3.3$, and Momose et al. 2014 at
$z=5.7$ and 6.6). The ${\rm SB_{Ly\alpha}}$ amplitude of Wisotzki et al.
(2018) is normalized such that $L_{\rm Lya}$ in $r\leq 1\arcsec$ matches to
that of our all sample (see Section 3.3). We represent the theoretical work
with the lines, whose colors indicate the Ly$\alpha$ emission originations:
the central galaxy (black), CGM (green), satellite galaxies (blue), and other
halos (orange). The models are taken from Mas-Ribas et al. (2017b, green and
blue dashed lines at $z=2.2$, 3.3, 5.7, 6.6), Kakiichi & Dijkstra (2018, black
dotted line at $z=2.2$), Byrohl et al. (2020, black, blue, and orange solid
lines at $z=3.3$), Lake et al. (2015, black and orange dash-dotted lines at
$z=3.3$), Dijkstra & Kramer (2012, black dotted line at $z=3.3$), and Zheng et
al. (2011, black and orange solid lines at $z=5.7$). See Table 4 for details
of these models. We normalize the ${\rm SB_{Ly\alpha}}$ profiles from the
model predictions such that the central galaxy ($r\leq 1\arcsec$) has the
$L_{\rm Ly\alpha}$ value similar to that observed at each redshift. The gray
shades illustrate the regions inside $R_{\rm vir}$ of a DMH with $M_{\rm
halo}=10^{11}\ M_{\odot}$.

In Figure 11, we compare the ${\rm SB_{Ly\alpha}}$ profiles that are observed
(Section 3.2) and those predicted by theoretical work. The observational
results are taken from Momose et al. (2016, $L_{\rm Ly\alpha}=10^{42.6}$ erg
s-1 subsample at $z=2.2$), Wisotzki et al. (2018, $L_{\rm Ly\alpha}>10^{42}$
erg s-1 subsample at $z=3-4$), Momose et al. (2014, $z=5.7$ and 6.6 LAEs), and
our all samples ($z=2.2-6.6$). Predicted profiles are taken from Zheng et al.
(2011), Dijkstra & Kramer (2012), Lake et al. (2015), Mas-Ribas et al.
(2017b), Kakiichi & Dijkstra (2018), and Byrohl et al. (2020), which are
summarized in Table 4. We plot these observational and theoretical results at
the nearest redshifts among $z=2.2$, 3.3, 5.7, and 6.6, except that the model
of Mas-Ribas et al. (2017b) is presented at all the redshifts, with the
correction for the cosmological dimming effect. We normalize all the profiles
such that they match at $r=1\arcsec$ in amplitude for precise comparison under
the same Ly$\alpha$ luminosity of the central galaxy (the models of Zheng et
al. 2011, Lake et al. 2015 and Byrohl et al. 2020 are normalized in each total
profile summing up different origins). We compare these observational and
theoretical results in the context of the Ly$\alpha$ photon origins in the
following subsections.

####  4.2.1 Central Galaxy

To discuss the contribution from the central galaxy, we compare the models of
Zheng et al. (2011), Dijkstra & Kramer (2012), Lake et al. (2015), Kakiichi &
Dijkstra (2018), and Byrohl et al. (2020). They applied Ly$\alpha$ radiative
transfer (RT) modeling to hydrodynamic cosmological galaxy formation
simulations at $z=5.7$, 3.1, $2-3$, and $2-5$, respectively, to investigate
Ly$\alpha$ photons produced by SF in the central galaxy and resonantly
scattered into the CGM. We take the ‘one-halo’ term model from Zheng et al.
(2011), where Ly$\alpha$ photons are scattered not only in the CGM but also
the IGM. The model of Lake et al. (2015) takes into account the contribution
from cooling radiation in addition to that from SF. The model of Kakiichi &
Dijkstra (2018) is taken in an approximated form of their Equation (28). The
model of Byrohl et al. (2020) considers recombination caused by ionizing
photons from SF and UVB as well as cooling via collisional de-excitation.

The black lines in Figure 11 represent the contribution from the central
galaxy: Kakiichi & Dijkstra (2018, dotted, $z=2.2$), Byrohl et al. (2020,
solid, $z=3.3$), Lake et al. (2015, dash-dotted, $z=3.3$), Dijkstra & Kramer
(2012, dotted, $z=3.3$), and Zheng et al. (2011, dotted, $z=5.7$). These
models successfully reproduce the observed ${\rm SB_{Ly\alpha}}$ profiles
inside the CGM ($r<R_{\rm vir}$) at $z=3.3$ and 5.7, implying resonant
scattering as a major source powering Ly$\alpha$ emission. This finding is
also supported by previous studies in other aspects, such as halo properties
(e.g., a halo luminosity-mass relation; Kusakabe et al. 2019) and kinematics
(e.g., a correlation between the peak velocity shift and the width of a
Ly$\alpha$ line, Leclercq et al. 2020; red peak dominated Ly$\alpha$ spectra,
Chen et al. 2021).

The one-halo term of Zheng et al. (2011) reproduces the observed ${\rm
SB_{Ly\alpha}}$ profile also from $r\sim R_{\rm vir}$ up to $\sim 1$ cMpc at
$z=5.7$. They considered resonant scattering in the IGM in addition to the
CGM. This effect leads to a plateau-like feature in a ${\rm SB_{Ly\alpha}}$
profile at $r\sim 0.3-1$ cMpc (see also Jeeson-Daniel et al. 2012), which is
similar to those at $z=3.3$ and 5.7. Therefore, we may interpret that the
observed extended Ly$\alpha$ emission outside the CGM is produced by resonant
scattering in the IGM, although the uncertainties are large. At $z=3.3$,
however, the models adopted here (Lake et al. 2015, Dijkstra & Kramer 2012,
and Byrohl et al. 2020) produce values lying far below the observed ${\rm
SB_{Ly\alpha}}$ profile beyond $300$ ckpc and do not reproduce a plateau-like
shape, unlike the one-halo term of Zheng et al. (2011). One possibility is
that the contribution from the central galaxy decreases from $z=5.7$ to 3.3,
but there is no evidence to confirm this. Alternatively, we suppose that this
discrepancy can be attributed to the different assumptions and incorporated
physics in the models, such as treatment of stellar radiation, dusts, and
scattering, especially beyond $R_{\rm vir}$. Nevertheless, since the one-halo
term of Zheng et al. (2011) reproduces the observed ${\rm SB_{Ly\alpha}}$
profile at $z=5.7$, we cannot rule out the possibility that scattered
Ly$\alpha$ photons originating from the central galaxy contribute to extended
Ly$\alpha$ emission beyond $R_{\rm vir}$. We need additional inputs on the
neutral hydrogen gas distribution outside $R_{\rm vir}$ to further determine
the contribution of resonant scattering.

The assumptions on $M_{\rm halo}$ values are unlikely to affect our discussion
here, because the models above use roughly similar $M_{\rm halo}$ values:
$10^{11.2}$ and $10^{11.5}\ M_{\odot}$ for Zheng et al. (2011) and Lake et al.
(2015), respectively. Byrohl et al. (2020) assume a stellar mass $M_{\star}$
range of $10^{8.5}-10^{9.5}\ M_{\odot}$ at $z=3$, which corresponds to $M_{\rm
halo}\sim 10^{11}\ M_{\odot}$ given the $M_{\star}/M_{\rm halo}$ ratio
obtained in Behroozi et al. (2019, see also ). These values are similar to
those obtained in the previous observations (e.g., Ouchi et al., 2010;
Kusakabe et al., 2018).

####  4.2.2 CGM

We next use the model from Mas-Ribas et al. (2017b) to investigate the
Ly$\alpha$ emission produced in the CGM. They constructed an analytical model
of fluorescent emission in the CGM caused by ionizing radiation from SF in the
central galaxy at $z=5.7$ and 6.6 (we apply this model also to $z=3.3$ and
2.2; see also Mas-Ribas & Dijkstra 2016). The model of Mas-Ribas et al.
(2017b) includes three free parameters: 1) the CGM structure, 2) SF rate
(SFR), and 3) radius $R_{\rm max}$. 1) We adopt the simplified clumpy outflow
model of Steidel et al. (2010) as the CGM structure (the choice here has only
a small impact on ${\rm SB_{Ly\alpha}}$ profiles; see Figure 2 of Mas-Ribas &
Dijkstra 2016). 2) We normalize the model with ${\rm SFR}=1$, 10, 10, and 20
$M_{\odot}$ yr-1 at $z=2.2$, 3.3, 5.7, 6.6, respectively. 3) ${\rm
SB_{Ly\alpha}}(b)$ is derived as the integration of Ly$\alpha$ emissivity at
radius $r$ over $b\leq r\leq R_{\rm max}$, where $b$ denotes the impact
parameter (see Equations 4 and 2 of Mas-Ribas & Dijkstra 2016 and Mas-Ribas et
al. 2017b, respectively). Here $R_{\rm max}$ represents the extent to which
Ly$\alpha$ emission contributes to ${\rm SB_{Ly\alpha}}$. We assume that
$R_{\rm max}$ is equivalent to $R_{\rm vir}$ for a DMH of $M_{\rm
halo}=10^{11}\ M_{\odot}$ ($=46$, 34, 22, and 20 pkpc at $z=2.2$, 3.3, 5.7,
and 6.6, respectively). Their model ignores the effect of resonant scattering.

The green dashed lines in Figure 11 represent the model of Mas-Ribas et al.
(2017b). The ${\rm SB_{Ly\alpha}}$ profiles predicted by their model are in
good agreement with those observed inside $R_{\rm vir}$ at all the redshifts.
However, because ${\rm SB_{Ly\alpha}}$ profiles drops sharply at $r=R_{\rm
max}$ according to this model, the emission originating from the CGM cannot
contribute to extended emission beyond $R_{\rm vir}$. Mas-Ribas & Dijkstra
(2016) and Mas-Ribas et al. (2017b) arbitrarily adopted much larger $R_{\rm
max}$ values ($>100$ pkpc) to reproduce the profiles obtained in Momose et al.
(2014), bust such large $R_{\rm max}$ values correspond to $M_{\rm
halo}>10^{12}\ M_{\odot}$, which is much larger than those observed (e.g.,
Ouchi et al., 2010; Kusakabe et al., 2018). If $R_{\rm max}$ is larger than
$R_{\rm vir}$, materials should exist outside the CGM ($\gtrsim 3$ times
larger scales than $R_{\rm vir}$) and produce fluorescent emission
contributing to extended Ly$\alpha$ emission. In either case, the observed
${\rm SB_{Ly\alpha}}$ profiles beyond 300 cpkc at $z=3.3$ cannot be reproduced
even with larger values of $R_{\rm max}$.

Overall, fluorescence in the CGM can power Ly$\alpha$ emission inside $R_{\rm
vir}$ according to the model of Mas-Ribas et al. (2017b), while it plays only
a marginal role beyond $R_{\rm vir}$. This behavior was also suggested by MUSE
UDF data (Gallego et al., 2018; Bacon et al., 2021). We note that the model of
Mas-Ribas et al. (2017b) ignores the scattering effect, which leads to a sharp
drop of ${\rm SB_{Ly\alpha}}$ at $r=R_{\rm max}$. Hence, it is necessary to
incorporate resonant scattering to extend Ly$\alpha$ emission when we rely on
the CGM fluorescence scenario.

The contribution from cooling radiation in the CGM remains unclear. Byrohl et
al. (2020) argued that cooling radiation dominates $\sim 30$ % of the total
Ly$\alpha$ emission at $r\gtrsim 20$ pkpc. On the other hand, Rosdahl &
Blaizot (2012) found that Ly$\alpha$ emission from cooling radiation is
centrally ($r<10$ pkpc) concentrated for a DMH with $M_{\rm halo}=10^{11}\
M_{\odot}$. We need additional models at larger scales to further discuss
whether cooling radiation contributes to extended Ly$\alpha$ emission or not.

####  4.2.3 Satellite Galaxies

We adopt another model from Mas-Ribas et al. (2017b), which predicts the
contribution from SF in satellite galaxies (see also Mas-Ribas et al. 2017a).
There are three free parameters in their model: 1) a clustering description,
2) Ly$\alpha$ escape fraction ($f_{\rm Lya}^{\rm esc}$), and 3) $R_{\rm max}$.
1) We assume that clustering follows a power-law two-point cross-correlation
function $\xi(r)$ of $\xi(r)=(r/r_{0})^{-\alpha}$ with the scale length
$r_{0}=4$ cMpc and index $\alpha=1.8$ (e.g., Ouchi et al., 2010; Harikane et
al., 2016; Bielby et al., 2017). 2) We fix $f_{\rm Lya}^{\rm esc}$ to 0.4,
while ${\rm SB_{Ly\alpha}}$ linearly depends on $f_{\rm Lya}^{\rm esc}$. 3)
Satellite galaxies are assumed to exist from $r=10$ pkpc to $R_{\rm max}$,
which we assume equal to $R_{\rm vir}$ in the same way as for the fluorescence
model (Section 4.2.2). We additionally take the ‘outer halo’ model from Byrohl
et al. (2020).

The blue dashed and solid lines represent the models of Mas-Ribas et al.
(2017b) and Byrohl et al. (2020), respectively. We find that the contribution
from satellite galaxies are negligible compared to the other contributions,
except at $r\sim 100$ ckpc at $z=2.2$. While we choose a power-low correlation
function to describe the clustering, other choices, such as the Navarro-Frenk-
White (NFW) profile (Navarro et al., 1997), reduce the ${\rm SB_{Ly\alpha}}$
values at $r>40$ ckpc (see the left panel of Figure 2 of Mas-Ribas et al.
2017a). When a power-law correlation function is assumed, larger $r_{0}$
values increase the overall profiles. However, unrealistically large values of
$r_{0}$ and $R_{\rm max}$ are necessary to reproduce extended Ly$\alpha$
emission beyond $R_{\rm vir}$ with satellite galaxies alone. Additionally, the
model overpredicts the observed ${\rm SB_{cont,\nu}}$ values when the model is
tuned to reproduce the observed ${\rm SB_{Ly\alpha}}$ profiles (Mas-Ribas et
al., 2017a).

For these reasons, we conclude that satellite galaxies are unlikely to
contribute to extended Ly$\alpha$ emission beyond $R_{\rm vir}$. This
conclusion is supported by the fact that emission is more extended in
Ly$\alpha$ than in UV continuum (Section 3.2; see also Momose et al. 2014,
2016; Wu et al. 2020), because the ${\rm SB_{cont,\nu}}$ profiles should be
extended similarly as the ${\rm SB_{Ly\alpha}}$ profiles if satellite galaxies
contribute to ${\rm SB_{Ly\alpha}}$.

####  4.2.4 Other Halos

Lastly, we compare the models for other halos taken from Zheng et al. (2011),
Lake et al. (2015), and Byrohl et al. (2020). From Lake et al. (2015) we
specifically adopt the model in which Ly$\alpha$ emission originates from
‘knots,’ the regions with high Ly$\alpha$ emissivity around the central
galaxy. The ‘two-halo’ term model is taken from Zheng et al. (2011).

We show these models with the orange lines in Figure 11: Byrohl et al. (2020,
solid, $z=3.3$), Lake et al. (2015, dash-dotted, $z=3.3$), and Zheng et al.
(2011, solid, $z=5.7$). At $z=3.3$, the contribution from other halos
predicted by Byrohl et al. (2020) agrees with the observed ${\rm
SB_{Ly\alpha}}$ profiles within the $1\sigma$ uncertainties. Although the
knots model of Lake et al. (2015) is limited to $r\lesssim 300$ cpkc, it
roughly reproduces the observed profiles at $r\sim(200-300)$ ckpc. On the
other hand, the two-halo term of Zheng et al. (2011) at $z=5.7$ significantly
overestimates the ${\rm SB_{Ly\alpha}}$ values beyond 100 ckpc. The amplitudes
of the models of Zheng et al. (2011) and Byrohl et al. (2020) differ by $\sim
1$ dex, similarly as we found in Section 4.2.1.

Kakuma et al. (2021) argued that the difference is caused because they masked
out bright objects. However, this interpretation is not necessarily
appropriate, since the profile of Byrohl et al. (2020) has an amplitude
similar to (or rather slightly higher than) profiles observed at $z=3.3$. The
two-halo term of Zheng et al. (2011) also overpredicts the ${\rm
SB_{Ly\alpha}}$ values of Momose et al. (2014). Nevertheless, we cannot rule
out the possibility that other halos contribute to extended Ly$\alpha$
emission beyond $R_{\rm vir}$, since the model of Byrohl et al. (2020) agrees
with the observed ${\rm SB_{Ly\alpha}}$ profiles. This suggestion is
consistent with Bacon et al. (2021), who identified very extended ($>300$
arcsec2 or $>2\times 10^{4}$ pkpc2) Ly$\alpha$ emission at $z\sim 3$ using
MUSE data; they found that 70 % of the total Ly$\alpha$ luminosity originates
from filamentary structures beyond the CGM. They argued that the extended
Ly$\alpha$ emission can be reproduced by a population of extremely faint
($<10^{40}$ erg s-1) galaxies under certain conditions, which correspond to
other halos considered in this subsection.

We note that the ${\rm SB_{Ly\alpha}}$ profile of our all sample at $z=5.7$
increases when we use only LAEs in the UD-COSMOS field, resulting in a smaller
amplitude gap between the two-halo term of Zheng et al. (2011) and our
profile. However, our profile based solely on the UD-COSMOS field disagrees
with those from the previous observational studies. We thus suppose that the
averaged profile is more appropriate for comparison against the models.

####  4.2.5 Overall Interpretation

In summary, ${\rm SB_{Ly\alpha}}$ profiles inside the CGM ($<R_{\rm vir}$) are
possibly explained either by scattered Ly$\alpha$ emission originating from
the central galaxy and/or fluorescent emission in the CGM. Meanwhile, extended
Ly$\alpha$ emission beyond $R_{\rm vir}$ is possibly powered either by
resonant scattering at large scales and/or contributed from other halos.
Fluorescence in the CGM and satellite galaxies are not sufficient to reproduce
the observed ${\rm SB_{Ly\alpha}}$ profiles beyond $R_{\rm vir}$.

We emphasize that the processes and origins of Ly$\alpha$ emission may differ
among LAEs. They may also vary according to radius and redshift even when we
focus on averaged profiles around different LAEs. However, our systematic
investigation of extended Ly$\alpha$ emission at $z=2-7$ is advantageous for a
comprehensive understanding of the processes and origins of extended
Ly$\alpha$ emission. More simulations focusing on large scales will help to
further distinguish the processes and origins.

##  5 Summary

In this paper, we investigated very extended Ly$\alpha$ emission around LAEs
at $z=2.2-6.6$ by applying the intensity mapping technique to the Subaru/HSC-
SSP and CHORUS data. Our major findings are summarized below:

  1. 1.

We calculated cross-correlation functions between 1781 LAEs at $z=2.2$, 3.3,
5.7 and 6.6 with Ly$\alpha$ emission traced by the NB387, NB527, NB816, and
NB921 images. A total of $\sim 1-2$ billion pixels were used to derive the
correlation function for each redshift. The deep ($m_{\rm NB,5\sigma}\sim 26$
mag) and wide ($\sim 4$ deg2 over the UD-COSMOS and UD-SXDS fields) images of
the HSC enabled us to detect very diffuse Ly$\alpha$ emission. We utilized
foreground objects (NonLAEs) to carefully estimate the systematics, including
the sky background and PSF. Subtracting these systematics, we identified the
Ly$\alpha$ emission of $\sim 10^{-20}-10^{-19}$ erg s-1 cm-2 arcsec-2 with
${\rm S/N}=4.1$ around the $z=3.3$ LAEs. Ly$\alpha$ emission extends beyond
the radial scale of the $R_{\rm vir}$ of a DMH with $10^{11}\ M_{\odot}$
($\sim 100$ ckpc), and up to $\sim 1$ cMpc. We also tentatively detected
Ly$\alpha$ emission beyond the $R_{\rm vir}$ scales at $z=5.7$, and 6.6 with
${\rm S/N}\sim 2$. Extended Ly$\alpha$ emission was also tentatively detected
around the $z=2.2$ LAEs when faint LAEs were excluded from the sample.

  2. 2.

We confirmed that the Ly$\alpha$ surface brightness (${\rm SB_{Ly\alpha}}$)
radial profiles around our LAEs agree well with those obtained in the previous
studies, when the LAEs have similar Ly$\alpha$ luminosity values.

  3. 3.

We compared the observed ${\rm SB_{Ly\alpha}}$ profiles across $z=2.2-6.6$,
finding no significant difference among the redshifts beyond the
uncertainties. Meanwhile, there is a potential increasing trend toward high
redshifts in the intrinsic ${\rm SB_{Ly\alpha}}$ profiles, which are corrected
for the cosmological dimming effect. The increasing trend roughly follows
$(1+z)^{3}$, which might be explained by the increasing density of the neutral
hydrogen gas due to the evolution of the cosmic volume.

  4. 4.

We compared the ${\rm SB_{Ly\alpha}}$ profiles obtained from the observational
and theoretical studies. We found that the observed ${\rm SB_{Ly\alpha}}$
profiles inside the CGM can be reproduced by the models in which Ly$\alpha$
photons originating from the central galaxy subsequently transfer into the CGM
via resonant scattering, or in which Ly$\alpha$ emission is produced in the
CGM via fluorescence due to ionizing photons. Extended Ly$\alpha$ emission
beyond $R_{\rm vir}$ may be reproduced by resonant scattering at large scales,
and/or emission originating from clustered halos around the targeted galaxy.
The CGM and satellite galaxies are unlikely to contribute to extended
Ly$\alpha$ emission beyond $R_{\rm vir}$.

This work, in conjunction with the previous observational studies, might
suggest that very extended diffuse Ly$\alpha$ emission beyond $R_{\rm vir}$
ubiquitously exist around LAEs at $z\sim 2-7$, not only around massive
galaxies. Deeper images obtained by larger-area surveys in the future should
enable further investigation of very extended Ly$\alpha$ emission at more
diffuse levels. Applying the intensity mapping technique to the emission of
multiple lines, such as H$\alpha$ and [Oiii], will help to further distinguish
physical processes and origins of extended Ly$\alpha$ emission, because they
trace different components (see Figures 6 and 12 of Mas-Ribas et al. 2017a and
Fujimoto et al. 2019). These emission lines will be observed with next-
generation facilities, such as the James-Webb Space Telescope (JWST), the
Nancy Grace Roman Space Telescope (NGRST), and Spectro-Photometer for the
History of the Universe, Epoch of Reionization and Ices Explorer (SPHEREx).

## Acknowledgements

We are grateful to Chris Byrohl for helpful comments and discussions
especially on simulations and physical interpretations.

The 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 the 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.

The Pan-STARRS1 Surveys (PS1) and the PS1 public science archive have been
made possible through contributions by 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, the 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 grant No.
AST-1238877, the University of Maryland, Eotvos Lorand University (ELTE), the
Los Alamos National Laboratory, and the Gordon and Betty Moore Foundation.

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.lsst.org.

This paper is based on data collected at the Subaru Telescope and retrieved
from the HSC data archive system, which is operated by the Subaru Telescope
and Astronomy Data Center (ADC) at National Astronomical Observatory of Japan.
Data analysis was in part carried out with the cooperation of Center for
Computational Astrophysics (CfCA), National Astronomical Observatory of Japan.
The Subaru Telescope is honored and grateful for the opportunity of observing
the Universe from Maunakea, which has the cultural, historical and natural
significance in Hawaii.

The NB387 filter was supported by KAKENHI (23244022) Grant-in-Aid for
Scientific Research (A) through the JSPS. The NB527 filter was supported by
KAKENHI (24244018) Grant-in-Aid for Scientific Research (A) through the JSPS.
The NB816 filter was supported by Ehime University. The NB921 filter was
supported by KAKENHI (23244025) Grant-in-Aid for Scientific Research (A)
through the JSPS.

This work is supported by the World Premier International Research Center
Initiative (WPI Initiative), MEXT, Japan, as well as KAKENHI Grant-in-Aid for
Scientific Research (A) (20H00180, and 21H04467) through the JSPS. This work
was supported by JSPS KAKENHI Grant Numbers 17H01114, 19J01222, 20J11993, and
21K13953.

This work was supported by the joint research program of the Institute for
Cosmic Ray Research (ICRR), University of Tokyo.

## References

  * Aihara et al. (2019) Aihara, H., AlSayyad, Y., Ando, M., et al. 2019, PASJ, 71, 114, doi: 10.1093/pasj/psz103
  * Arrigoni Battaia et al. (2019) Arrigoni Battaia, F., Hennawi, J. F., Prochaska, J. X., et al. 2019, MNRAS, 482, 3162, doi: 10.1093/mnras/sty2827
  * Bacon et al. (2021) Bacon, R., Mary, D., Garel, T., et al. 2021, A&A, 647, A107, doi: 10.1051/0004-6361/202039887
  * Behroozi et al. (2019) Behroozi, P., Wechsler, R. H., Hearin, A. P., & Conroy, C. 2019, MNRAS, 488, 3143, doi: 10.1093/mnras/stz1182
  * Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393, doi: 10.1051/aas:1996164
  * Bielby et al. (2017) Bielby, R. M., Shanks, T., Crighton, N. H. M., et al. 2017, MNRAS, 471, 2174, doi: 10.1093/mnras/stx1772
  * Bond et al. (2010) Bond, N. A., Feldmeier, J. J., Matković, A., et al. 2010, ApJ, 716, L200, doi: 10.1088/2041-8205/716/2/L200
  * Borisova et al. (2016) Borisova, E., Cantalupo, S., Lilly, S. J., et al. 2016, ApJ, 831, 39, doi: 10.3847/0004-637X/831/1/39
  * Bosch et al. (2018) Bosch, J., Armstrong, R., Bickerton, S., et al. 2018, PASJ, 70, S5, doi: 10.1093/pasj/psx080
  * Byrohl et al. (2020) Byrohl, C., Nelson, D., Behrens, C., et al. 2020, arXiv e-prints, arXiv:2009.07283.  https://arxiv.org/abs/2009.07283
  * Cantalupo et al. (2014) Cantalupo, S., Arrigoni-Battaia, F., Prochaska, J. X., Hennawi, J. F., & Madau, P. 2014, Nature, 506, 63, doi: 10.1038/nature12898
  * Cantalupo et al. (2005) Cantalupo, S., Porciani, C., Lilly, S. J., & Miniati, F. 2005, ApJ, 628, 61, doi: 10.1086/430758
  * Carilli (2011) Carilli, C. L. 2011, ApJ, 730, L30, doi: 10.1088/2041-8205/730/2/L30
  * Chen et al. (2021) Chen, Y., Steidel, C. C., Erb, D. K., et al. 2021, arXiv e-prints, arXiv:2104.10173.  https://arxiv.org/abs/2104.10173
  * Claeyssens et al. (2019) Claeyssens, A., Richard, J., Blaizot, J., et al. 2019, MNRAS, 489, 5022, doi: 10.1093/mnras/stz2492
  * Comaschi & Ferrara (2016a) Comaschi, P., & Ferrara, A. 2016a, MNRAS, 455, 725, doi: 10.1093/mnras/stv2339
  * Comaschi & Ferrara (2016b) —. 2016b, MNRAS, 463, 3078, doi: 10.1093/mnras/stw2199
  * Croft et al. (2018) Croft, R. A. C., Miralda-Escudé, J., Zheng, Z., Blomqvist, M., & Pieri, M. 2018, MNRAS, 481, 1320, doi: 10.1093/mnras/sty2302
  * Croft et al. (2016) Croft, R. A. C., Miralda-Escudé, J., Zheng, Z., et al. 2016, MNRAS, 457, 3541, doi: 10.1093/mnras/stw204
  * Dawson et al. (2013) Dawson, K. S., Schlegel, D. J., Ahn, C. P., et al. 2013, AJ, 145, 10, doi: 10.1088/0004-6256/145/1/10
  * Dijkstra & Kramer (2012) Dijkstra, M., & Kramer, R. 2012, MNRAS, 424, 1672, doi: 10.1111/j.1365-2966.2012.21131.x
  * Eisenstein et al. (2011) Eisenstein, D. J., Weinberg, D. H., Agol, E., et al. 2011, AJ, 142, 72, doi: 10.1088/0004-6256/142/3/72
  * Fardal et al. (2001) Fardal, M. A., Katz, N., Gardner, J. P., et al. 2001, ApJ, 562, 605, doi: 10.1086/323519
  * Faucher-Giguère et al. (2010) Faucher-Giguère, C.-A., Kereš, D., Dijkstra, M., Hernquist, L., & Zaldarriaga, M. 2010, ApJ, 725, 633, doi: 10.1088/0004-637X/725/1/633
  * Feldmeier et al. (2013) Feldmeier, J. J., Hagen, A., Ciardullo, R., et al. 2013, ApJ, 776, 75, doi: 10.1088/0004-637X/776/2/75
  * Fisher (1970) Fisher, Ronald Aylmer, S. 1970, Statistical methods for research workers, 14th edn. (Edinburgh : Oliver and Boyd)
  * Fonseca et al. (2017) Fonseca, J., Silva, M. B., Santos, M. G., & Cooray, A. 2017, MNRAS, 464, 1948, doi: 10.1093/mnras/stw2470
  * Fujimoto et al. (2019) Fujimoto, S., Ouchi, M., Ferrara, A., et al. 2019, ApJ, 887, 107, doi: 10.3847/1538-4357/ab480f
  * Furlanetto et al. (2005) Furlanetto, S. R., Schaye, J., Springel, V., & Hernquist, L. 2005, ApJ, 622, 7, doi: 10.1086/426808
  * Gallego et al. (2018) Gallego, S. G., Cantalupo, S., Lilly, S., et al. 2018, MNRAS, 475, 3854, doi: 10.1093/mnras/sty037
  * Garel et al. (2021) Garel, T., Blaizot, J., Rosdahl, J., et al. 2021, MNRAS, 504, 1902, doi: 10.1093/mnras/stab990
  * Goerdt et al. (2010) Goerdt, T., Dekel, A., Sternberg, A., et al. 2010, MNRAS, 407, 613, doi: 10.1111/j.1365-2966.2010.16941.x
  * Gong et al. (2011) Gong, Y., Cooray, A., Silva, M. B., Santos, M. G., & Lubin, P. 2011, ApJ, 728, L46, doi: 10.1088/2041-8205/728/2/L46
  * Goto et al. (2009) Goto, T., Utsumi, Y., Furusawa, H., Miyazaki, S., & Komiyama, Y. 2009, MNRAS, 400, 843, doi: 10.1111/j.1365-2966.2009.15486.x
  * Haiman et al. (2000) Haiman, Z., Spaans, M., & Quataert, E. 2000, ApJ, 537, L5, doi: 10.1086/312754
  * Harikane et al. (2016) Harikane, Y., Ouchi, M., Ono, Y., et al. 2016, ApJ, 821, 123, doi: 10.3847/0004-637X/821/2/123
  * Harikane et al. (2018) Harikane, Y., Ouchi, M., Shibuya, T., et al. 2018, ApJ, 859, 84, doi: 10.3847/1538-4357/aabd80
  * Harikane et al. (2019) Harikane, Y., Ouchi, M., Ono, Y., et al. 2019, ApJ, 883, 142, doi: 10.3847/1538-4357/ab2cd5
  * Hayashi et al. (2020) Hayashi, M., Shimakawa, R., Tanaka, M., et al. 2020, PASJ, 72, 86, doi: 10.1093/pasj/psaa076
  * Hayashino et al. (2004) Hayashino, T., Matsuda, Y., Tamura, H., et al. 2004, AJ, 128, 2073, doi: 10.1086/424935
  * Hayes et al. (2011) Hayes, M., Schaerer, D., Östlin, G., et al. 2011, ApJ, 730, 8, doi: 10.1088/0004-637X/730/1/8
  * Hayes et al. (2013) Hayes, M., Östlin, G., Schaerer, D., et al. 2013, ApJ, 765, L27, doi: 10.1088/2041-8205/765/2/L27
  * Hayes et al. (2014) Hayes, M., Östlin, G., Duval, F., et al. 2014, ApJ, 782, 6, doi: 10.1088/0004-637X/782/1/6
  * Higuchi et al. (2019) Higuchi, R., Ouchi, M., Ono, Y., et al. 2019, ApJ, 879, 28, doi: 10.3847/1538-4357/ab2192
  * Inoue et al. (2014) Inoue, A. K., Shimizu, I., Iwata, I., & Tanaka, M. 2014, MNRAS, 442, 1805, doi: 10.1093/mnras/stu936
  * Inoue et al. (2018) Inoue, A. K., Hasegawa, K., Ishiyama, T., et al. 2018, PASJ, 70, 55, doi: 10.1093/pasj/psy048
  * Inoue et al. (2020) Inoue, A. K., Yamanaka, S., Ouchi, M., et al. 2020, PASJ, 72, 101, doi: 10.1093/pasj/psaa100
  * Itoh et al. (2018) Itoh, R., Ouchi, M., Zhang, H., et al. 2018, ApJ, 867, 46, doi: 10.3847/1538-4357/aadfe4
  * Jeeson-Daniel et al. (2012) Jeeson-Daniel, A., Ciardi, B., Maio, U., et al. 2012, MNRAS, 424, 2193, doi: 10.1111/j.1365-2966.2012.21378.x
  * Jiang et al. (2013) Jiang, L., Egami, E., Fan, X., et al. 2013, ApJ, 773, 153, doi: 10.1088/0004-637X/773/2/153
  * Kakiichi & Dijkstra (2018) Kakiichi, K., & Dijkstra, M. 2018, MNRAS, 480, 5140, doi: 10.1093/mnras/sty2214
  * Kakuma et al. (2021) Kakuma, R., Ouchi, M., Harikane, Y., et al. 2021, arXiv e-prints, arXiv:1906.00173.  https://arxiv.org/abs/1906.00173
  * Kikuta et al. (2019) Kikuta, S., Matsuda, Y., Cen, R., et al. 2019, PASJ, 71, L2, doi: 10.1093/pasj/psz055
  * Kollmeier et al. (2010) Kollmeier, J. A., Zheng, Z., Davé, R., et al. 2010, ApJ, 708, 1048, doi: 10.1088/0004-637X/708/2/1048
  * Konno et al. (2016) Konno, A., Ouchi, M., Nakajima, K., et al. 2016, ApJ, 823, 20, doi: 10.3847/0004-637X/823/1/20
  * Konno et al. (2018) Konno, A., Ouchi, M., Shibuya, T., et al. 2018, PASJ, 70, S16, doi: 10.1093/pasj/psx131
  * Kovetz et al. (2017) Kovetz, E. D., Viero, M. P., Lidz, A., et al. 2017, arXiv e-prints, arXiv:1709.09066.  https://arxiv.org/abs/1709.09066
  * Kusakabe et al. (2018) Kusakabe, H., Shimasaku, K., Ouchi, M., et al. 2018, PASJ, 70, 4, doi: 10.1093/pasj/psx148
  * Kusakabe et al. (2019) Kusakabe, H., Shimasaku, K., Momose, R., et al. 2019, PASJ, 71, 55, doi: 10.1093/pasj/psz029
  * Lake et al. (2015) Lake, E., Zheng, Z., Cen, R., et al. 2015, ApJ, 806, 46, doi: 10.1088/0004-637X/806/1/46
  * Laursen & Sommer-Larsen (2007) Laursen, P., & Sommer-Larsen, J. 2007, ApJ, 657, L69, doi: 10.1086/513191
  * Laursen et al. (2011) Laursen, P., Sommer-Larsen, J., & Razoumov, A. O. 2011, ApJ, 728, 52, doi: 10.1088/0004-637X/728/1/52
  * Leclercq et al. (2017) Leclercq, F., Bacon, R., Wisotzki, L., et al. 2017, A&A, 608, A8, doi: 10.1051/0004-6361/201731480
  * Leclercq et al. (2020) Leclercq, F., Bacon, R., Verhamme, A., et al. 2020, A&A, 635, A82, doi: 10.1051/0004-6361/201937339
  * Li et al. (2016) Li, T. Y., Wechsler, R. H., Devaraj, K., & Church, S. E. 2016, ApJ, 817, 169, doi: 10.3847/0004-637X/817/2/169
  * Martin et al. (2014) Martin, D. C., Chang, D., Matuszewski, M., et al. 2014, ApJ, 786, 107, doi: 10.1088/0004-637X/786/2/107
  * Mas-Ribas & Dijkstra (2016) Mas-Ribas, L., & Dijkstra, M. 2016, ApJ, 822, 84, doi: 10.3847/0004-637X/822/2/84
  * Mas-Ribas et al. (2017a) Mas-Ribas, L., Dijkstra, M., Hennawi, J. F., et al. 2017a, ApJ, 841, 19, doi: 10.3847/1538-4357/aa704e
  * Mas-Ribas et al. (2017b) Mas-Ribas, L., Hennawi, J. F., Dijkstra, M., et al. 2017b, ApJ, 846, 11, doi: 10.3847/1538-4357/aa8328
  * Matsuda et al. (2012) Matsuda, Y., Yamada, T., Hayashino, T., et al. 2012, MNRAS, 425, 878, doi: 10.1111/j.1365-2966.2012.21143.x
  * Mitchell et al. (2021) Mitchell, P. D., Blaizot, J., Cadiou, C., et al. 2021, MNRAS, 501, 5757, doi: 10.1093/mnras/stab035
  * Momose et al. (2021a) Momose, R., Shimasaku, K., Nagamine, K., et al. 2021a, ApJ, 912, L24, doi: 10.3847/2041-8213/abf04c
  * Momose et al. (2021b) Momose, R., Shimizu, I., Nagamine, K., et al. 2021b, ApJ, 911, 98, doi: 10.3847/1538-4357/abe1b9
  * Momose et al. (2014) Momose, R., Ouchi, M., Nakajima, K., et al. 2014, MNRAS, 442, 110, doi: 10.1093/mnras/stu825
  * Momose et al. (2016) —. 2016, MNRAS, 457, 2318, doi: 10.1093/mnras/stw021
  * Navarro et al. (1997) Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493, doi: 10.1086/304888
  * Oke & Gunn (1983) Oke, J. B., & Gunn, J. E. 1983, ApJ, 266, 713, doi: 10.1086/160817
  * Ono et al. (2021) Ono, Y., Itoh, R., Shibuya, T., et al. 2021, ApJ, 911, 78, doi: 10.3847/1538-4357/abea15
  * Östlin et al. (2009) Östlin, G., Hayes, M., Kunth, D., et al. 2009, AJ, 138, 923, doi: 10.1088/0004-6256/138/3/923
  * Ouchi et al. (2020) Ouchi, M., Ono, Y., & Shibuya, T. 2020, ARA&A, 58, 617, doi: 10.1146/annurev-astro-032620-021859
  * Ouchi et al. (2010) Ouchi, M., Shimasaku, K., Furusawa, H., et al. 2010, ApJ, 723, 869, doi: 10.1088/0004-637X/723/1/869
  * Ouchi et al. (2018) Ouchi, M., Harikane, Y., Shibuya, T., et al. 2018, PASJ, 70, S13, doi: 10.1093/pasj/psx074
  * Patrício et al. (2016) Patrício, V., Richard, J., Verhamme, A., et al. 2016, MNRAS, 456, 4191, doi: 10.1093/mnras/stv2859
  * Péroux & Howk (2020) Péroux, C., & Howk, J. C. 2020, ARA&A, 58, 363, doi: 10.1146/annurev-astro-021820-120014
  * Pullen et al. (2014) Pullen, A. R., Doré, O., & Bock, J. 2014, ApJ, 786, 111, doi: 10.1088/0004-637X/786/2/111
  * Rauch et al. (2008) Rauch, M., Haehnelt, M., Bunker, A., et al. 2008, ApJ, 681, 856, doi: 10.1086/525846
  * Rosdahl & Blaizot (2012) Rosdahl, J., & Blaizot, J. 2012, MNRAS, 423, 344, doi: 10.1111/j.1365-2966.2012.20883.x
  * Scoville et al. (2007) Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1, doi: 10.1086/516585
  * Sekiguchi et al. (2005) Sekiguchi, K., Akiyama, M., Furusawa, H., et al. 2005, in Multiwavelength Mapping of Galaxy Formation and Evolution, ed. A. Renzini & R. Bender, 82, doi: 10.1007/10995020_12
  * Shibuya et al. (2018a) Shibuya, T., Ouchi, M., Konno, A., et al. 2018a, PASJ, 70, S14, doi: 10.1093/pasj/psx122
  * Shibuya et al. (2018b) Shibuya, T., Ouchi, M., Harikane, Y., et al. 2018b, PASJ, 70, S15, doi: 10.1093/pasj/psx107
  * Silva et al. (2013) Silva, M. B., Santos, M. G., Gong, Y., Cooray, A., & Bock, J. 2013, ApJ, 763, 132, doi: 10.1088/0004-637X/763/2/132
  * Smit et al. (2017) Smit, R., Swinbank, A. M., Massey, R., et al. 2017, MNRAS, 467, 3306, doi: 10.1093/mnras/stx245
  * Smith et al. (2019) Smith, A., Ma, X., Bromm, V., et al. 2019, MNRAS, 484, 39, doi: 10.1093/mnras/sty3483
  * Smith et al. (2018) Smith, A., Tsang, B. T. H., Bromm, V., & Milosavljević, M. 2018, MNRAS, 479, 2065, doi: 10.1093/mnras/sty1509
  * Steidel et al. (2011) Steidel, C. C., Bogosavljević, M., Shapley, A. E., et al. 2011, ApJ, 736, 160, doi: 10.1088/0004-637X/736/2/160
  * Steidel et al. (2010) Steidel, C. C., Erb, D. K., Shapley, A. E., et al. 2010, ApJ, 717, 289, doi: 10.1088/0004-637X/717/1/289
  * Swinbank et al. (2007) Swinbank, A. M., Bower, R. G., Smith, G. P., et al. 2007, MNRAS, 376, 479, doi: 10.1111/j.1365-2966.2007.11454.x
  * Tumlinson et al. (2017) Tumlinson, J., Peeples, M. S., & Werk, J. K. 2017, ARA&A, 55, 389, doi: 10.1146/annurev-astro-091916-055240
  * Verhamme et al. (2012) Verhamme, A., Dubois, Y., Blaizot, J., et al. 2012, A&A, 546, A111, doi: 10.1051/0004-6361/201218783
  * Wisotzki et al. (2016) Wisotzki, L., Bacon, R., Blaizot, J., et al. 2016, A&A, 587, A98, doi: 10.1051/0004-6361/201527384
  * Wisotzki et al. (2018) Wisotzki, L., Bacon, R., Brinchmann, J., et al. 2018, Nature, 562, 229, doi: 10.1038/s41586-018-0564-6
  * Wu et al. (2020) Wu, J., Jiang, L., & Ning, Y. 2020, ApJ, 891, 105, doi: 10.3847/1538-4357/ab7333
  * Xue et al. (2017) Xue, R., Lee, K.-S., Dey, A., et al. 2017, ApJ, 837, 172, doi: 10.3847/1538-4357/837/2/172
  * Zhang et al. (2020) Zhang, H., Ouchi, M., Itoh, R., et al. 2020, ApJ, 891, 177, doi: 10.3847/1538-4357/ab7917
  * Zheng et al. (2011) Zheng, Z., Cen, R., Weinberg, D., Trac, H., & Miralda-Escudé, J. 2011, ApJ, 739, 62, doi: 10.1088/0004-637X/739/2/62