# UV & Ly$\alpha$ halos of Ly$\alpha$ emitters across environments at
$z=2.84$111This research is based on data collected at the Subaru Telescope,
which is operated by the National Astronomical Observatory of Japan.

Satoshi Kikuta National Astronomical Observatory of Japan, 2-21-1, Osawa,
Mitaka, Tokyo, 181-8588, Japan  Center for Computational Sciences, University
of Tsukuba, 1-1-1, Ten-nodai, Tsukuba, Ibaraki, 305-8577, Japan  Yuichi
Matsuda Graduate University for Advanced Sciences (SOKENDAI), 2-21-1, Osawa,
Mitaka, Tokyo, 181-8588, Japan  National Astronomical Observatory of Japan,
2-21-1, Osawa, Mitaka, Tokyo, 181-8588, Japan  Shigeki Inoue Department of
Cosmosciences, Graduates School of Science, Hokkaido University, Sapporo,
Hokkaido 060-0810, Japan  Charles C. Steidel California Institute of
Technology, MS 249-17, Pasadena, CA 91125, USA  Renyue Cen School of Physics,
Zhejiang University, Hangzhou, 310058, China  Zheng Zheng Department of
Physics and Astronomy, University of Utah, 115 South 1400 East, Salt Lake
City, UT 84112, USA  Hidenobu Yajima  Center for Computational Sciences,
University of Tsukuba, 1-1-1, Ten-nodai, Tsukuba, Ibaraki, 305-8577, Japan
Rieko Momose Carnegie Observatories, 813 Santa Barbara Street, Pasadena, CA
91101, USA  Department of Astronomy, School of Science, The University of
Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo, 113-0033, Japan  Masatoshi Imanishi
Graduate University for Advanced Sciences (SOKENDAI), 2-21-1, Osawa, Mitaka,
Tokyo, 181-8588, Japan  National Astronomical Observatory of Japan, 2-21-1,
Osawa, Mitaka, Tokyo, 181-8588, Japan  Yutaka Komiyama Department of Advanced
Sciences, Faculty of Science and Engineering, Hosei University, 3-7-2 Kajino-
cho, Koganei-shi, Tokyo 184-8584, Japan

(Received 2022 May 20; Accepted 2023 February 24)

###### Abstract

We present UV & Ly$\alpha$ radial surface brightness (SB) profiles of
Ly$\alpha$ emitters (LAEs) at $z=2.84$ detected with the Hyper Suprime-Cam
(HSC) on the Subaru Telescope. The depth of our data, together with the wide
field coverage including a protocluster, enable us to study the dependence of
Ly$\alpha$ halos (LAHs) on various galaxy properties, including Mpc-scale
environments. UV and Ly$\alpha$ images of 3490 LAEs are extracted, and
stacking the images yields SB sensitivity of $\sim 1\times
10^{-20}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ in Ly$\alpha$,
reaching the expected level of optically thick gas illuminated by the UV
background at $z\sim 3$. Fitting of the two-component exponential function
gives the scale-lengths of $1.56\pm 0.01$ and $10.4\pm 0.3$ pkpc. Dividing the
sample according to their photometric properties, we find that while the
dependence of halo scale-length on environment outside of the protocluster
core is not clear, LAEs in the central regions of protoclusters appear to have
very large LAHs which could be caused by combined effects of source
overlapping and diffuse Ly$\alpha$ emission from cool intergalactic gas
permeating the forming protocluster core irradiated by active members. For the
first time, we identify “UV halos” around bright LAEs which are probably due
to a few lower-mass satellite galaxies. Through comparison with recent
numerical simulations, we conclude that, while scattered Ly$\alpha$ photons
from the host galaxies are dominant, star formation in satellites evidently
contributes to LAHs, and that fluorescent Ly$\alpha$ emission may be boosted
within protocluster cores at cosmic noon and/or near bright QSOs.

galaxies: high-redshift, galaxies: formation, intergalactic medium

††journal: ApJ††facilities: Subaru (HSC)††software:  SExtractor (Bertin &
Arnouts, 1996), IRAF

##  1 Introduction

Gas inflow from the cosmic web fuels star formation and active galactic
nucleus (AGN) activity in galaxies. Following these activities, outflows play
a significant role in expelling or heating a large amount of gas, thereby
reducing further star formation and supermassive black hole (SMBH) growth. The
circumgalactic medium (CGM) contains vital information on these flow
components and thus holds an important key to revealing galaxy evolution (see
Tumlinson et al., 2017, for a recent review). The CGM of $z\gtrsim 2$ star-
forming galaxies is now routinely detected as diffuse Ly$\alpha$ nebulae, or
Ly$\alpha$ halos (LAHs), around star-forming galaxies such as Ly$\alpha$
emitters (LAEs) and Lyman break galaxies (LBGs) at high-redshift both
individually (Rauch et al., 2008; Wisotzki et al., 2016; Leclercq et al.,
2017; Erb et al., 2018; Bacon et al., 2021; Kusakabe et al., 2022) and through
a stacking technique (Hayashino et al., 2004; Steidel et al., 2011; Matsuda et
al., 2012; Feldmeier et al., 2013; Momose et al., 2014, 2016; Xue et al.,
2017; Lujan Niemeyer et al., 2022a, b). Together with a technique based on
absorption lines in the spectra of neighboring background sources (Adelberger
et al., 2005; Steidel et al., 2010; Rudie et al., 2012; Chen et al., 2020;
Muzahid et al., 2021), LAHs have provided a crucial empirical window into the
CGM of distant galaxies.

To extract useful information on the CGM from LAHs, the physical origins of
the Ly$\alpha$ emission should be identified. Ly$\alpha$ surface brightness
(SB) profiles of LAHs hold the key since they are determined by the
distribution and kinematics of gas and the relative importance of various
Ly$\alpha$ production mechanisms such as scattering of Ly$\alpha$ photons from
host galaxies, star formation in neighboring galaxies, collisional excitation
of inflow gas powered by gravitational energy (sometimes called gravitational
cooling radiation), and recombination following photoionization by external
sources, often referred to as “fluorescence.” Theoretical studies have
attempted to reproduce and predict observed LAHs by considering these
mechanisms (Haiman & Rees, 2001; Dijkstra & Loeb, 2009; Goerdt et al., 2010;
Kollmeier et al., 2010; Faucher-Giguère et al., 2010; Zheng et al., 2011;
Dijkstra & Kramer, 2012; Rosdahl & Blaizot, 2012; Yajima et al., 2013; Cen &
Zheng, 2013; Cantalupo et al., 2014; Lake et al., 2015; Mas-Ribas & Dijkstra,
2016; Gronke & Bird, 2017; Mitchell et al., 2021; Byrohl et al., 2021).
Powerful outflows (so-called “superwind”, Taniguchi & Shioya, 2000; Mori et
al., 2004) have also been proposed to excite gas, but often for more
energetic/massive counterparts such as Ly$\alpha$ blobs (LABs) and nebulae
around QSOs and radio galaxies222There is no clear demarcation, but
conventionally LABs refer to extended Ly$\alpha$ nebulae that are particularly
bright ($L_{\mathrm{Ly\alpha}}>10^{43}\mathrm{~{}erg~{}s^{-1}}$) but without
obvious AGN activity at optical wavelengths..

Dependence of LAH shapes on e.g., their hosts’ halo mass and large-scale
overdensity is naturally expected because both gas and sources of Ly$\alpha$
and ionizing photons are more abundant in massive halos and/or denser
environments (Zheng et al., 2011; Mas-Ribas & Dijkstra, 2016; Kakiichi &
Dijkstra, 2018). Current simulations cannot treat all relevant physics with
sufficient accuracy and it is only very recently that such predictions are
reported in the literature with a statistical number of simulated galaxies
(e.g., Byrohl et al., 2021). Observations of LAHs can help theorists pin down
which Ly$\alpha$ production processes are at work by revealing the dependence
of LAH SB profile shapes on various properties such as UV and Ly$\alpha$
luminosity, Ly$\alpha$ equivalent width (EWLyα,0, Momose et al. 2014, 2016;
Wisotzki et al. 2016, 2018; Leclercq et al. 2017), and the large-scale
environment (Matsuda et al., 2012; Xue et al., 2017). The results in the
literature are, however, far from converging (see e.g., Figure 12 of Leclercq
et al. (2017)). It is often parametrized by an exponential function
$\propto\exp(-r/r_{\mathrm{h}})$ with a scale-length $r_{\mathrm{h}}$, which
is fit to observed profiles. The reported scale-lengths of individual LAEs as
a function of UV magnitude show a large scatter (from $<1$ physical kpc (pkpc
hereafter) to $\sim 10$ pkpc) and the relation for stacked LAEs shows large
differences as well. In the case of large-scale environment, relevant
observations of LAHs of LAEs are still scarce. First, Steidel et al. (2011)
found very large LAHs with a scale-length of 25 pkpc around LBGs in three
protoclusters at $z=2.3$-3.1. Following this result, Matsuda et al. (2012)
suggested that the scale-length of LAHs of LAEs are proportional to galaxy
overdensity squared $\delta^{2}$. On the other hand, Xue et al. (2017) found
no such dependence with LAEs in two overdense regions at $z=2.66$ and
$z=3.78$.

A major problem with some previous work is poor sensitivity. While only a few
studies investigated LAHs with deep images of a fairly large sample of LAEs
($N_{\mathrm{LAE}}>2000$; Matsuda et al. 2012; Momose et al. 2014), others
used the insufficient number of LAEs ($N_{\mathrm{LAE}}\sim$ a few$\times
100$–1000) and/or images taken with 4m telescopes (e.g., Feldmeier et al.,
2013; Xue et al., 2017). Because LAHs beyond the virial radii of LAEs have
extremely low SB
($<10^{-19}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$), interpretation
of relations derived from shallow data would not be straightforward. Even with
sufficiently deep data, the extent of LAHs is difficult to measure.
Disagreements among different studies can be attributed in part to different
fitting methods, fitting range, radial bin size of SB profiles, and
sensitivity of observational data among each study. To alter this situation
and to provide a firm observational basis for theorists, a well-controlled
statistical sample of LAEs drawn from a wide dynamic range of environments,
with sufficiently deep images, is required.

In this paper, we present a new LAH study with deep narrow-band (NB468) data
obtained using the Hyper Suprime-Cam (HSC; Miyazaki et al., 2012) on the
Subaru Telescope toward the HS1549 protocluster at $z=2.84$ (Trainor &
Steidel, 2012; Mostardi et al., 2013; Kikuta et al., 2019) to probe what
shapes LAHs. Thanks to the HSC’s large field of view ($\phi\sim 1.5$ deg,
corresponding to 160 comoving Mpc at $z=2.84$), we can construct a large LAE
sample across environments from a protocluster to surrounding lower density
fields at the same time. The sample size of our study of $N=3490$ is one of
the largest to date, giving robust UV and Ly$\alpha$ SB profiles to be
compared with simulations. As a result, we for the first time detect “UV
halos” which directly prove the contribution of star formation in satellite
galaxies. Moreover, we detect very extended LAHs for the protocluster LAEs
which suggest an important role of locally enhanced ionizing radiation fields
for LAHs. This paper is structured as follows. In Section 2, we describe our
LAE sample, followed by how we divide them for the stacking analyses described
in Section 3. The results of the analyses are shown in Section 4. Based on
these, we present discussion in Section 5 and summarize the work in Section 6.
Throughout this paper, we use the AB magnitude system and assume a cosmology
with $\Omega_{\mathrm{m}}=0.3$, $\Omega_{\Lambda}=0.7$, and
$H_{0}=70\mathrm{~{}km~{}s^{-1}~{}Mpc^{-1}}$, unless otherwise noted.

##  2 Data and sample

We used the LAE sample described by Kikuta et al. (2019). Here, we highlight
key properties of the data and refer readers to the aforementioned paper for
details. The target protocluster contains a hyperluminous QSO, namely
HS1549+1919 ($L_{\mathrm{1450\AA}}=1.5\times 10^{14}~{}L_{\odot}$, Trainor &
Steidel 2012; Mostardi et al. 2013) at its center. The field was observed in
the g-band (central wavelength $\lambda_{\mathrm{c}}=4712$Å,
$\mathrm{FWHM}=1479$Å) and NB468 ($\lambda_{\mathrm{c}}=4683$Å,
$\mathrm{FWHM}=88$Å) narrow-band filters. The global sky subtraction
method333https://hsc.mtk.nao.ac.jp/pipedoc_e/e_tips/skysub.html#global-sky was
used to estimate and subtract the sky on scales larger than that of individual
CCDs in the mosaic, with a grid size of 6000 pixels (17′) not to subtract
diffuse emission. The FWHMs of stellar sources in the final images are
0$.\\!\\!\arcsec$77 (0$.\\!\\!\arcsec$65) for g-band (NB468). The NB468 image
was smoothed with a circular Gaussian function to match the FWHM of stellar
sources in the g-band image. A 5$\sigma$ limiting magnitude measured with
1$.\\!\\!\arcsec$5 diameter aperture is 27.4 (26.6) mag for the g-band (NB468)
image. Our criterion for NB excess, $\mathrm{g}-\mathrm{NB468}>0.5$,
corresponds to $\mathrm{EW_{obs}}>45$ Å after considering the 0.1 mag offset.
SExtractor (Bertin & Arnouts, 1996) was used to perform 1$.\\!\\!\arcsec$5
aperture photometry with double-image mode, using the NB468 image as the
detection band and a background mesh size of 64 pixels ($=11$″) for local sky
estimation. This small mesh size is only used for LAE detection, since this is
optimal for detecting compact sources such as distant galaxies. As a result,
we detected 3490 LAEs within $r<36\arcmin$ from the QSO position. Their sky
distribution is shown in Figure 1.

To study the dependence of LAHs on various photometric galaxy properties, we
divided the sample into several groups according to the following five
quantities: UV magnitude, Ly$\alpha$ luminosity, rest-frame Ly$\alpha$
equivalent width (EW0), environment, and distance from the HLQSO, as
summarized in Table 1. The first three quantities, which are derived with
1$.\\!\\!\arcsec$5 aperture ($=9$ pixels $=2\times\mathrm{[PSF~{}FWHM]}$), are
obviously not independent of one another. UV continuum and Ly$\alpha$
luminosity are both good proxies for star formation rate (SFR), but the
latter’s resonant nature and a different impact of dust extinction thereby
make interpretation difficult (e.g., Scarlata et al., 2009). Moreover, UV
slope or hardness of UV emission is a strong function of age and metallicity,
and thus the Ly$\alpha$ equivalent width changes accordingly (Schaerer, 2003;
Hashimoto et al., 2017a). There is a known observational relation between the
EW distribution and UV limiting magnitude, known as “the Ando effect” (a
fainter UV threshold tends to include more high-EW LAEs; Ando et al.,
2006)444Note that this effect could be purely attributed to the selection bias
and not intrinsic (Nilsson et al., 2009a; Zheng et al., 2010).. Since
Ly$\alpha$ emission can be powered by mechanisms other than star formation in
galaxies, UV magnitude is the most robust to use here as a tracer for SFR.
Binning in UV magnitude, Ly$\alpha$ luminosity, and Ly$\alpha$ equivalent
width (as well as distance from the HLQSO) were done so that each subsample
has approximately the same number of LAEs ($N_{\mathrm{LAE}}\sim 700$).

The projected distance from the HLQSO is used to test whether QSO radiation
affects the LAHs of surrounding LAEs. HS1549 is so luminous that the entire
field covered by the HSC could experience a higher ionizing radiation field
than the cosmic average at $z\sim 3$ (Haardt & Madau, 2012) if QSO radiation
has had time to propagate, and additional ionization induced by the QSO could
increase Ly$\alpha$ luminosity of LAEs in the field (see Section 5.2.1 for
discussion). Here, we use a projected distance from the QSO. Note, however,
that the NB468 filter has an FWHM of $\Delta\lambda=88$Å or $\Delta z=0.075$,
corresponding to 19 pMpc width when centered at $z=2.84$. This brings
uncertainty in a line-of-sight distance and therefore also in real (3D)
distance. The boundaries defining the LAE subsamples are indicated by
concentric circles in Figure 1. Lastly, grouping based on environment is done
using projected (surface) LAE overdensity $\delta\equiv(n-\bar{n})/\bar{n}$
measured locally with an aperture radius of 1$.\\!\\!\arcmin$8 ($=0.83$ pMpc)
to be consistent with the measurement of Matsuda et al. (2012). Here, $n$ and
$\bar{n}$ are, respectively, the number of LAEs within a circle with a
1$.\\!\\!\arcmin$8 radius centered at the position of interest, and its
average over the entire field. This division is visually illustrated in Figure
1 by gray contours (see also Figures 1 and 2 of Kikuta et al. 2019). The
boundary of $\delta=2.5$ is set to only include protocluster members in the
densest subgroup by visual inspection. The next boundary $\delta=1$ is set by
a trade-off between tracing sufficiently dense regions and including adequate
numbers to allow sufficient S/N in stacks. The remainder are set so as to
roughly equalize the number of LAEs in each bin.

Figure 2 shows cumulative distributions of UV and Ly$\alpha$ luminosity, rest-
frame Ly$\alpha$ equivalent width, distance from the HLQSO, and overdensity
$\delta$ for each division, illustrating correlations between these
quantities. We note that the odd behavior of the thin blue and purple curves
in the second panel from left in the third row of Figure 2 is likely to be
artificial; EWs of LAEs not detected in g band (above $2\sigma$) are just
lower limits. The samples of faint low-EW LAEs are also incomplete due to the
lack of dynamic range in the measurement of $\mathrm{g-NB468}$, distorting the
distribution. The protocluster subsample (LAEs with $\delta>2.5$; thick red
curves in panels shown with $\delta$) evidently stands out among others, while
the projected distance from the HLQSO seems not to make a significant
difference except for $\delta$. These differences should be kept in mind when
interpreting the result in Section 4. We visualize these quantities also in
Figure 1 which roughly includes all of above information.

Figure 1:  Sky distribution of our LAEs (colored dots). The HLQSO is located
at the origin (0,0). North is up and east to the left. The colors indicate
their rest-frame Ly$\alpha$ equivalent width and the sizes indicate their UV
absolute magnitude (the larger, the brighter). Gray contours indicate their
overdensity $\delta$, with each level showing a different group defined in
Table 1. Concentric black circles indicate radii with which the distance
subsamples are defined.  Figure 2:  From top to bottom, cumulative
distributions of UV magnitude, Ly$\alpha$ luminosity, rest-frame Ly$\alpha$
EW, distance from the HLQSO, and overdensity of all subsamples described in
Table 1 and the top-left panel and panels in the second row. In each panel,
the quantity used to divide the sample is labeled. Thicker curves present [UV
brighter, Ly$\alpha$ brighter, lower EW, farther from the HLQSO, denser]
subsamples. Note that in the second row, last column, the $2.5<\delta$
subsample contains 4 sources with
$L_{\mathrm{Ly\alpha}}>10^{44}\mathrm{~{}erg~{}s^{-1}}$ (including the HLQSO)
and thus the red curve does not reach to unity.  Table 1: Subsample used in
our stacking analysis and the result of fitting performed in Section 4.2.

quantity | criteria | median | $N$ | $C_{\mathrm{1,Ly\alpha}}$ | $r_{\mathrm{1,Ly\alpha}}$ | $C_{\mathrm{2,Ly\alpha}}$ | $r_{\mathrm{2,Ly\alpha}}$ | $C_{\mathrm{3,Ly\alpha}}$ | $\alpha$
---|---|---|---|---|---|---|---|---|---
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10)
– | All | – | 3490 (430) | $2.742\pm 0.011$ | $1.564\pm 0.013$ | $1.062\pm 0.043$ | $10.384\pm 0.257$ | $1.575\pm 0.009$ | $2.285\pm 0.005$
| $M_{\mathrm{UV}}<-19.2$ | -19.62 | 690 (0) | $3.782\pm 0.018$ | $2.106\pm 0.019$ | $1.369\pm 0.085$ | $12.795\pm 0.478$ | $2.058\pm 0.013$ | $2.057\pm 0.005$
| $-19.2<M_{\mathrm{UV}}<-18.6$ | -18.88 | 696 (0) | $2.689\pm 0.020$ | $1.893\pm 0.031$ | $0.949\pm 0.128$ | $10.144\pm 0.758$ | $1.575\pm 0.017$ | $2.200\pm 0.009$
UV magnitude | $-18.6<M_{\mathrm{UV}}<-18.0$ | -18.31 | 773 (0) | $2.852\pm 0.024$ | $1.494\pm 0.025$ | $0.661\pm 0.078$ | $11.403\pm 0.850$ | $1.785\pm 0.025$ | $2.439\pm 0.013$
| $-18.0<M_{\mathrm{UV}}<-17.4$ | -17.73 | 648 (0) | $2.754\pm 0.031$ | $1.284\pm 0.028$ | $0.467\pm 0.059$ | $13.307\pm 1.166$ | $1.850\pm 0.037$ | $2.622\pm 0.020$
| $-17.4<M_{\mathrm{UV}}$ | -16.92 | 683 (430) | $2.884\pm 0.034$ | $1.088\pm 0.027$ | $0.398\pm 0.048$ | $14.282\pm 1.259$ | $2.155\pm 0.053$ | $2.843\pm 0.026$
| $42.25<\log L_{\mathrm{Ly\alpha}}\mathrm{[erg~{}s^{-1}]}$ | 42.40 | 647 (0) | $6.253\pm 0.025$ | $1.566\pm 0.014$ | $2.609\pm 0.130$ | $9.446\pm 0.271$ | $3.630\pm 0.021$ | $2.292\pm 0.005$
| $42.05<\log L_{\mathrm{Ly\alpha}}<42.25$ | 42.14 | 833 (5) | $3.382\pm 0.024$ | $1.484\pm 0.026$ | $1.423\pm 0.155$ | $8.092\pm 0.479$ | $2.092\pm 0.022$ | $2.392\pm 0.010$
Ly$\alpha$ luminosity | $41.95<\log L_{\mathrm{Ly\alpha}}<42.05$ | 41.99 | 610 (26) | $2.726\pm 0.027$ | $1.495\pm 0.026$ | $0.558\pm 0.054$ | $14.929\pm 1.011$ | $1.610\pm 0.025$ | $2.383\pm 0.014$
| $41.85<\log L_{\mathrm{Ly\alpha}}<41.95$ | 41.90 | 645 (80) | $2.484\pm 0.030$ | $1.317\pm 0.034$ | $0.642\pm 0.096$ | $10.271\pm 0.966$ | $1.607\pm 0.032$ | $2.544\pm 0.020$
| $\log L_{\mathrm{Ly\alpha}}<41.85$ | 41.79 | 755 (319) | $2.386\pm 0.032$ | $1.124\pm 0.030$ | $0.403\pm 0.044$ | $14.702\pm 1.182$ | $1.600\pm 0.039$ | $2.696\pm 0.025$
|  $12<\mathrm{EW_{0,Ly\alpha}}<30$ Å | 21.1 Å | 644 (0) | $2.045\pm 0.016$ | $2.404\pm 0.036$ | $0.821\pm 0.085$ | $13.848\pm 0.848$ | $1.081\pm 0.011$ | $1.942\pm 0.008$
|  $30<\mathrm{EW_{0,Ly\alpha}}<55$ Å | 42.4 Å | 735 (0) | $2.439\pm 0.020$ | $1.883\pm 0.030$ | $1.094\pm 0.090$ | $11.721\pm 0.576$ | $1.298\pm 0.013$ | $2.086\pm 0.008$
Ly$\alpha$ EW0 |  $55<\mathrm{EW_{0,Ly\alpha}}<90$ Å | 70.5 Å | 698 (0) | $2.898\pm 0.026$ | $1.480\pm 0.031$ | $1.137\pm 0.137$ | $8.862\pm 0.609$ | $1.753\pm 0.023$ | $2.375\pm 0.012$
|  $90<\mathrm{EW_{0,Ly\alpha}}<160$ Å | 121 Å | 727 (100) | $3.017\pm 0.027$ | $1.381\pm 0.025$ | $0.680\pm 0.071$ | $11.880\pm 0.811$ | $1.904\pm 0.028$ | $2.497\pm 0.014$
| 160 Å$<\mathrm{EW_{0,Ly\alpha}}$ | 216 Å | 686 (330) | $3.433\pm 0.032$ | $1.178\pm 0.024$ | $0.646\pm 0.078$ | $10.951\pm 0.869$ | $2.450\pm 0.043$ | $2.734\pm 0.018$
|  $d_{\mathrm{Q}}<6.2$ pMpc | 4.05 pMpc | 679 (81) | $2.709\pm 0.023$ | $1.659\pm 0.026$ | $0.991\pm 0.068$ | $13.178\pm 0.586$ | $1.418\pm 0.016$ | $2.161\pm 0.010$
|  $6.2<d_{\mathrm{Q}}<9.5$ pMpc | 7.95 pMpc | 739 (81) | $2.710\pm 0.024$ | $1.540\pm 0.034$ | $1.353\pm 0.158$ | $8.356\pm 0.533$ | $1.600\pm 0.020$ | $2.304\pm 0.011$
Projected distance |  $9.5<d_{\mathrm{Q}}<12.0$ pMpc | 10.7 pMpc | 633 (70) | $2.687\pm 0.029$ | $1.442\pm 0.043$ | $1.639\pm 0.248$ | $6.983\pm 0.537$ | $1.641\pm 0.023$ | $2.358\pm 0.013$
from the HLQSO |  $12.0<d_{\mathrm{Q}}<14.8$ pMpc | 13.5 pMpc | 778 (111) | $2.838\pm 0.023$ | $1.563\pm 0.026$ | $0.876\pm 0.091$ | $10.690\pm 0.676$ | $1.688\pm 0.021$ | $2.341\pm 0.011$
|  $14.8<d_{\mathrm{Q}}<16.9$ pMpc | 15.9 pMpc | 661 (87) | $2.766\pm 0.025$ | $1.567\pm 0.026$ | $0.856\pm 0.063$ | $13.593\pm 0.665$ | $1.490\pm 0.018$ | $2.236\pm 0.011$
| $2.5<\delta$ | 3.78 | 55 (2) | $2.295\pm 0.058$ | $2.440\pm 0.079$ | $0.615\pm 0.069$ | $43.237\pm 5.108$ | $0.924\pm 0.026$ | $1.715\pm 0.019$
| $1.0<\delta<2.5$ | 1.33 | 433 (57) | $2.686\pm 0.031$ | $1.574\pm 0.039$ | $1.122\pm 0.140$ | $10.020\pm 0.737$ | $1.524\pm 0.024$ | $2.264\pm 0.014$
Environment | $0.3<\delta<1.0$ | 0.63 | 944 (116) | $2.781\pm 0.022$ | $1.543\pm 0.024$ | $1.005\pm 0.070$ | $11.743\pm 0.525$ | $1.533\pm 0.017$ | $2.257\pm 0.010$
| $-0.15<\delta<0.3$ | 0.05 | 1076 (146) | $2.693\pm 0.021$ | $1.555\pm 0.026$ | $1.181\pm 0.102$ | $9.482\pm 0.483$ | $1.545\pm 0.016$ | $2.279\pm 0.009$
| $-1.0<\delta<-0.15$ | -0.30 | 982 (109) | $2.812\pm 0.022$ | $1.513\pm 0.029$ | $1.249\pm 0.143$ | $8.284\pm 0.525$ | $1.703\pm 0.019$ | $2.352\pm 0.010$

Note. — Column (1) and (2): quantity and criteria used to define subsamples,
Column (3): median value of each quantity of each subsample, Column (4): the
number of LAEs which satisfy the criteria described in the column (2). The
numbers in parenthesis are those of g-band undetected sources. Column (5):
$C_{\mathrm{1,Ly\alpha}}$ in units of
$10^{-17}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$, Column (6):
$r_{\mathrm{1,Ly\alpha}}$ in units of physical kpc, Column (7):
$C_{\mathrm{2,Ly\alpha}}$ in units of
$10^{-18}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$, Column (8):
$r_{\mathrm{2,Ly\alpha}}$ in units of physical kpc, Column (9):
$C_{\mathrm{3,Ly\alpha}}$ in units of
$10^{-16}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$, Column (10):
power-law index $\alpha$. Note that UV magnitude, Ly$\alpha$ luminosity, and
its EW are derived with 1$.\\!\\!\arcsec$5 aperture. The uncertainties of
fitting parameters include fitting errors only.

##  3 Analyses

###  3.1 Image Stacking

Before stacking, additional sky subtraction needs to be made, as the global
sky subtraction (see Section 2) alone is generally insufficient to eliminate
the influence of artifacts such as halos of bright stars. Since we are
interested in diffuse and extended components, we use a background mesh size
of 176 pixels ($=30$ arcsec) for additional background evaluation of g-band
and NB images using SExtractor and then subtract the sky. Then a Ly$\alpha$
(continuum) image was created by subtracting the g-band (NB468) image from the
NB468 (g-band) image after scaling by their relative zero points and
considering the difference in the filters’ transmission curves assuming flat
continuum (see Appendix B of Mawatari et al., 2012, for details). A
segmentation image of the continuum image was used for masking. The
segmentation image is an output of SExtractor and it specifies which pixel is
detected as a source. If detected, that pixel has a non-zero integer value
that corresponds to the ID number in the output catalog. We set DETECT_MINAREA
and DETECT_THRESH parameters to 5 and 2.0. We confirmed that the overall
results do not depend significantly on the choice of the threshold value.

We created cutout Ly$\alpha$ and UV continuum images centered on each LAE. The
centers of the LAEs are identified as their centroids in the NB image. The
mask is applied to each cutout image. If the LAE at the center of the cutout
image is detected in the continuum image, masking for the object is turned off
so as not to underestimate Ly$\alpha$ and UV continuum emission near the
center. Stacking is executed using the IRAF task “imcombine” with
median/average and with/without sigma clipping to further eliminate unrelated
signals. Ly$\alpha$ SB profiles of stacked images are measured in a series of
annuli with a width of 2 pixels.

###  3.2 Uncertainties and Limitations

To estimate the noise level of the stacked image, we created cutouts of the
Ly$\alpha$ and continuum image centered on randomly selected points in the
field. After applying the continuum source mask, these “sky cutouts” are
stacked to make stacked sky images and their radial SB profiles are measured
in the same manner as stacked LAE images. This time the turning off of the
masking of some sources is not employed. The 1$\sigma$ noise level of each
annulus is estimated by repeating this 1000 times and deriving the standard
deviation of the distribution of total count in the annulus555When estimating
the noise level, using the same aperture shape (in this case, annuli with a
width of 2 pixels) is important because spatial correlations between
neighboring pixels affect the noise level differently with different aperture
shapes.. As shown in Figure 3 (Top), we confirmed that the noise level
decreases almost as $\propto N^{-1/2}$. The result also demonstrates that the
choice of stacking method does not make a significant difference in the noise
level except for the case of average stacking without sigma clipping (in this
case, too many artifacts remain; green points for “ave no clip” except for
$N=3000$ are far above the graph’s upper bound). Thus, we decided to
extrapolate this relation between the noise and $N$ (the number of images for
stacking) to estimate the noise level of stacked images with any $N$ rather
than to iterate 1000 times for every possible $N$ in the column (4) of Table
1. The same method is used to estimate the noise level of the continuum image,
which behaves almost $\propto N^{-1/2}$ as well. In figure 3 (Bottom), the
noise level of the sky stack as a function of radius is shown for the
$N_{\mathrm{stack}}=700$ case. Due to the pixel spatial correlation, they do
not perfectly decrease as (the number of pixels in each annulus)-0.5. Stacking
methods again do not change the result.

Figure 3: (Top) Estimating the sky noise level of the Ly$\alpha$ image. The
x-axis shows how many sky cutouts were stacked, and the y-axis shows the
estimated noise level evaluated in 50th annulus from the center (thus $r=$ 100
pixel or 17 arcsec) containing 1270 pixels. Different points indicate the
different stacking methods (average VS. median, with 3$\sigma$ or 5$\sigma$
clipping VS. no sigma clipping), which are well converged except for the case
with average without sigma clipping. The red curve is a fitting function with
a form $a\times N^{b}$ to the blue circles, which is consistent with inverse
square root proportionality ($b=-0.493\pm 0.0074$). Symbols for ”ave no clip”
are mostly out of the upper boundary. (Bottom) The estimated noise level of
each annulus for the $N_{\mathrm{stack}}=700$ case. Their behavior with
respect to the stacking methods is the same as the above figure.

At the same time, the average values of the sum of the counts in each annulus
were measured. Due to systematic errors and sky residuals, the average sky
counts are not exactly equal to zero. To correct this effect, we subtract the
average sky value when we derive the radial profile. Typical sky value of the
Ly$\alpha$ and continuum images are $\sim-5\times
10^{-21}~{}\mathrm{erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ and $\sim 1.2\times
10^{-32}~{}\mathrm{erg~{}s^{-1}~{}cm^{-2}~{}Hz^{-1}~{}arcsec^{-2}}$.

Since the Ly$\alpha$ image was created by subtracting the g-band image from
the NB image, any difference between the PSFs of the images could produce
spurious patterns around sources in the Ly$\alpha$ image. Even if the simple
Gaussian smoothing done in Section 2 can match the FWHM of stellar sources, it
cannot exactly match the shape of the PSFs of the two images. Moreover, the
shape of the PSF at a large radius may introduce additional errors. To examine
the detailed shapes of the PSF in the two images, we first select bright
unsaturated sources from a source catalog using the SExtractor output
CLASS_STAR, which is a parameter characterizing the stellarity of sources.
CLASS_STAR is 1 if an object is a point source and drops to 0 if extended.
Here we use following criteria: CLASS_STAR$>0.95$ and $18<\mathrm{g}<22$. In
total, 3980 sources are stacked to determine the central part of the PSFs in
the images666Initially we divide this sample into two; one for sources
distributed in the inner part of the field and the other for the outer part.
The profiles of the stacked image of the two subsamples are almost identical.
Thus we conclude that variation of the PSF within the field is minor and
ignore the effect in the following analyses.. To determine the much fainter
outer part of the PSFs, we extracted stars with magnitude
$13<\mathrm{g_{SDSS}}<15$ from the SDSS DR14 catalog (Abolfathi et al., 2018).
After excluding stars with bright nearby companions and/or obviously extended
sources when seen with our deep images, images of 113 bright stars are
stacked. Since point sources in this magnitude range start to saturate, the
PSF of the brighter sources is connected at $r=20$ pixels or 3.4 arcsec with
that of fainter sources, following a method described in Infante-Sainz et al.
(2019). Derived PSFs from 0.17 arcsec to 40 arcsec are shown in Figure 4. The
PSF of the NB image is slightly smaller than that of the g band at $r=$ a few
arcsec. The PSFs of both bands beyond several arcsec agree very well. They are
not Gaussian-like and have power-law tails with a slope of $\sim-2.8$.

Figure 4: Normalized PSFs of the NB (red circles) and g (blue squares) band
images. The inner and outer PSFs are connected at $r=3.4$ arcsec (20 pixels).
Gray curves are the extrapolated inner and outer PSFs beyond the junction
radius.

To check whether the slight difference between the broadband and NB PSFs
affects our surface brightness measurement, we created a stacked “non-LAE”
image, following a method described in Momose et al. (2014). “Non-LAE” sources
are defined as objects not selected as LAEs which have almost the same
distribution in the FWHMNB468 vs. NB468 magnitude plane as the real LAEs
(Figure 5). Since the majority of LAEs are distributed in the range
$0.75\mathrm{~{}arcsec}<\mathrm{FWHM}<3.25$ arcsec and
$24<\mathrm{NB468}<26.5$, we select non-LAE sources from this range for
stacking. Any signal detected in the stacked Ly$\alpha$ image of non-LAEs can
be used to estimate the effect not only of the PSF difference but also of
other unknown systematics such as errors associated with flat-fielding and sky
subtraction as discussed in Feldmeier et al. (2013).

Figure 5:  FWHM VS NB468 plot of LAE (red) and non-LAE (blue). Top and right
panels show histograms of NB468 and FWHM, respectively. Their distribution is
almost identical.

##  4 Results

###  4.1 Stacked profiles and effects of systematics

Figure 6 shows the median-stacked Ly$\alpha$ and continuum images of all LAEs
and non-LAEs without sigma clipping. The SB profiles of them are shown in
Figure 7. We confirmed that the profiles do not depend on the stacking methods
(except for average stacking without sigma clipping; see Figure 3); they show
$>1\sigma$ deviation only at very low-S/N regime near $r\sim 100$ pkpc.
Hereafter, we present the results for median stacking without sigma clipping.
The non-LAE has a negative ring-like structure around the center. This
probably arises from the slight differences in their colors and PSFs of g-band
and NB images. Still, the absolute value of the Ly$\alpha$ SB profile of the
non-LAE is about an order of magnitude smaller than that of LAEs in Figure 7.
Beyond 2 arcsec, the SB profile of the non-LAE stack is almost consistent with
the sky value and thus we conclude the effect of the PSF difference is
negligible, in particular at large radii of $r>2$ arcsec of most interest to
the present work. The PSF, which is shown with the gray curve in Figure 7,
drops much more rapidly than the Ly$\alpha$ profile of LAEs. From the above
arguments, we conclude that LAHs around LAEs at $z=2.84$ are robustly detected
down to $\sim 1\times
10^{-20}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ and the effects of
systematic errors cannot have a significant impact on the derived Ly$\alpha$
SB profiles out to $\sim 100$ pkpc.

On the other hand, the UV SB profile of LAEs seems to be negative beyond 20
pkpc. Similar patterns can be also seen in some previous studies performing
stacking analyses (Matsuda et al., 2012; Momose et al., 2016) but the exact
reasons were never identified in the literature. When estimating the sky value
with SExtractor, pixels with counts above some threshold are masked. However,
since our LAEs are selected with the NB image as a detection band, some LAEs
in our sample are too faint in the UV continuum image to be masked. In
addition, when creating the continuum image, the Ly$\alpha$ contribution is
subtracted even if the source is not detected or significantly affected by the
sky noise in the g-band image for UV-faint sources. These effects cause
oversubtraction in the continuum image of UV-faint LAEs, affecting the UV SB
profiles of subsamples that contain UV-faint LAEs. In our sample, there are
430 (1391) LAEs with $<2\sigma$ ($<5\sigma$) detection in the g-band image,
respectively. The numbers of LAEs with $<2\sigma$ detection in all subsamples
are also shown in Table 1. Thus, UV SB profiles of subsamples that contain
many UV non-detected LAEs should be interpreted with caution 777It is also
possible that a color term difference within the subsamples affects our SB
measurement, although a correlation between $\beta$ and UV absolute magnitude
$\mathrm{M_{UV}}$ of LAEs is, though still debated, weak
($\mathrm{d}\beta/\mathrm{dM_{UV}}\sim-0.1$–0.0 Hathi et al., 2016; Hashimoto
et al., 2017b).

Figure 6:  Stacked Ly$\alpha$ (top) and continuum (bottom) image of all LAEs
(left) and non-LAEs (right). The size of each image is $\sim 25\arcsec\times
25\arcsec$ or $\sim 200\times 200$ pkpc at $z=2.84$. Contours correspond to
$3\times 10^{-18},1\times 10^{-18},1\times
10^{-19}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ in the Ly$\alpha$
images and $3\times 10^{-31},1\times 10^{-31},3\times
10^{-32}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}Hz^{-1}~{}arcsec^{-2}}$ in
contiuum images.  Figure 7:  Radial SB profile of LAEs (solid) and non-LAEs
(dashed) in the Ly$\alpha$ (orange) and continuum (blue curve) images are
shown. Gray dashed curve shows 1$\sigma$ noise level. On the bottom side, SB
profiles below $1\times
10^{-21}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ level are shown in
linear scale. The normalized PSF of the NB image is shown with the gray solid
curve.

To check whether or not the mesh size for sky estimation matters, we derived a
stacked Ly$\alpha$ profile of all LAEs with sky mesh sizes different from 30
arcsec. A larger mesh size enables us to probe possible large-scale emission
around LAEs, at the same time increasing errors due to residual non-
astrophysical signals (which stem from e.g. halos around bright stars). A
smaller mesh size may lead to oversubtraction of the real signal while
reducing the errors described above. To find a better compromise, we tested
sky mesh sizes of 1 arcmin, 2 arcmin, and 11 arcsec (64 pixels, the default
mesh size used in LAE selection in Section 2). The 1$\sigma$ errors and
residual sky emission to be subtracted were derived in the same way as
described in Section 3.2. In Figure 8, we showed the results of this test.
Except for the case of 11 arcsec (blue curve), the derived Ly$\alpha$ SB
profiles are all consistent with each other within uncertainty. Also, no
systematic trend is evident with increasing mesh size in the slight offset
seen in the outer part. This suggests the effect of oversubtraction of diffuse
Ly$\alpha$ emission is minor at this sensitivity, on this scale. However, as
the larger sky mesh sizes are utilized, residual artificial emission around
bright stars in the Ly$\alpha$ images becomes clearer as well. On the other
hand, a mesh size of 64 pixels $=11$ arcsec $=85$ pkpc clearly oversubtract
halo emission. Considering that LAHs are detected out to $\sim 100$ pkpc, a
mesh size of 85 pkpc, which is comparable to the extent of LAHs, should not be
used. The same trend is seen in the continuum image as well. Thus we decided
to use 0.5 arcmin $=$ 30 arcsec sky mesh.

Figure 8:  Radial Ly$\alpha$ SB profiles of all LAEs with different sky mesh
size. Red, yellow, green, and blue curves respectively indicate profiles
derived with sky mesh sizes of 2 arcmin, 1 arcmin, 0.5 arcmin, and 11 arcsec.

To quantify the extent of SB profiles, we performed fits to both UV and
Ly$\alpha$ SB profiles using the following exponential function(s) and power-
law function:

$\displaystyle\mathrm{PSF}*\left[C_{1}\exp{\left(-\frac{r}{r_{1}}\right)}+C_{2}\exp{\left(-\frac{r}{r_{2}}\right)}\right]$
(1) $\displaystyle\mathrm{PSF}*C_{3}r^{-\alpha},$ (2)

where “$\mathrm{PSF}*$” means convolution with the measured PSF of NB468.
While exponential functions have commonly been used in previous observational
work, a power-law function is motivated by an analytical model by Kakiichi &
Dijkstra (2018). $C_{2}$ is set to zero for 1-component exponential fitting
and let $r_{1}<r_{2}$ if otherwise, thus $r_{1}$ is scale-length for the core
component and $r_{2}$ is for the halo component. Unlike most previous work, we
do not assume that the scale-lengths for the core component of Ly$\alpha$ and
UV SB profile are the same, and they were fitted separately. The result shown
in Figure 9 clearly demonstrates the need for non-zero $C_{2}$ or a power-law
function to fit the Ly$\alpha$ SB profile, while a single exponential function
will do for the UV SB profile fitting. The 2-component exponential fit to the
UV SB profile does not converge, and the power-law fit clearly deviates from
the observed UV profile.

Figure 9:  Radial Ly$\alpha$ (left) and UV (right) SB profiles of all LAEs are
shown with solid orange/blue curves with errorbars. Dashed, dot-dashed, and
dotted curves respectively show the result of fitting with two- and one-
component exponential functions and a power-law function. The normalized PSF
of the NB image is shown with the thin gray curve. Downward triangles show
1$\sigma$ error levels after residual sky subtraction.

###  4.2 Subsamples

The stacked Ly$\alpha$ and UV continuum images of all subsamples are shown in
Figure 10; their corresponding Ly$\alpha$ and UV continuum SB profiles are
shown in Figures 11 and 12, respectively with fitting curves. Rest-frame
Ly$\alpha$ equivalent widths calculated in each annulus are also plotted with
orange dots in Figure 12. The resulting fit parameters for the Ly$\alpha$ SB
profiles are given in Table 1, and those of UV SB profiles in Table 2. In the
right panel of Figure 10, the effect of oversubtraction discussed in Section
4.1 is clearly manifested by the black ring-like structures around the central
emission in the stacked UV continuum images of the UV/Ly$\alpha$ faintest and
highest-EW subsamples (rightmost three panels from top to middle). Again, the
UV SB (and the EW) profiles of UV-faint subsamples should be interpreted with
caution.

Figure 10:  Stacked Ly$\alpha$ (left, orange) and continuum (right, blue)
images of different subsamples in the asinh color stretch. From top to bottom,
we show stacked images of UV magnitude, $L_{\mathrm{Ly\alpha}}$, EW0,Lyα,
distance from the HLQSO, and environment subsamples. From left to right,
median [UV magnitude, $L_{\mathrm{Ly\alpha}}$, EW0,Lyα, distance from the
HLQSO, $\delta$] of each subsample respectively become [fainter, fainter,
larger, larger, smaller]. The size of each image is $\sim 200\times 200$ pkpc.
Contours correspond to $3\times 10^{-18},1\times 10^{-18},1\times
10^{-19}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ in the Ly$\alpha$
images and $3\times 10^{-31},1\times 10^{-31},3\times
10^{-32}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}Hz^{-1}~{}arcsec^{-2}}$ in
continuum images. .  Figure 11:  Radial Ly$\alpha$ SB profiles of all
subsamples (specified in Table 1, solid orange curve with errorbars) with
fitting curves. Green dashed, red dot-dashed, and gray dotted curves show the
result of fitting with two- and one-component exponential functions and a
power-law function. Downward triangles show 1$\sigma$ error levels after
residual sky subtraction. Thin gray curve shows the normalized PSF shape.
Figure 12:  Radial UV SB profiles of all subsamples (specified in Table 1,
solid blue curve with errorbars) with fitting curves. Green dashed and red
dot-dashed curves show the result of fitting with two- and one-component
exponential functions. The gray dashed curve is Ly$\alpha$ SB profiles
converted from $F_{\mathrm{Ly\alpha}}$ units to $f_{\nu}$ units by simply
dividing by the FWHM of the NB filter. Orange dot indicates rest-frame
Ly$\alpha$ equivalent width in Åcalculated in each annulus, with its value on
the right axis. Downward and upward triangles show 1$\sigma$ limits of UV
emission and equivalent width (some datapoints in the rightmost panels are
above the upper boundary) after residual sky subtraction. Thin gray curve
shows the normalized PSF shape. Note the difference in the range of the x-axis
from that of Figure 11.

We detected Ly$\alpha$ emission more extended than UV (stellar) emission in
all subsamples (Figure 11), while most UV SB profiles can be fitted with the
one-component exponential function (Figure 12). However, UV SB profiles of the
UV/Ly$\alpha$ brightest and lowest-EW subsamples clearly require two-component
exponential functions and indeed can be fitted remarkably well. To our
knowledge, this is the first detection of “UV halos” in high-redshift LAEs. We
confirmed that this detection is robust against the choice of used stacking
methods, the sky mesh size (as long as it is not too small), and the masking
threshold. In addition, these subsamples are securely detected in the g-band
and thus the effect of oversubtraction (Section 4.1) should be minor. For
subsamples for which two-component exponential fitting is well converged, we
show the resultant fitting parameters in Table 3. Figures 19 and 20 in
Appendix A show results in a different manner for a clearer comparison within
each photometric property. In Figure 19, clear systematic differences can be
seen in bins of UV, $L_{\mathrm{Ly\alpha}}$, and EW0,Lyα in a way that
UV/Ly$\alpha$-bright LAEs and low EW0,Lyα LAEs have larger LAHs. We see hints
of UV halos also in the second and third UV brightest subsamples in the upper
right panel of Figure 19. On the other hand, the difference of profiles for
the projected distance and local environment subsamples is not obvious except
for the protocluster subsample (those with $\delta>2.5$). There is no
significant difference in UV SB profiles in both cases except for the
protocluster subsample, which contains more UV bright galaxies as seen in
Figure 2.

The scale-lengths and power-law index of fitting functions can be used for
more quantitative discussion. In Figure 13, we show the results of fitting for
each photometric property. $r_{\mathrm{1,UV}}$ plotted in Figures 13 and 14
are those obtained from a one-component fit (including subsamples with UV-
halo), while $r_{\mathrm{1,Ly\alpha}}$ are from a two-component fit. In the
left three columns of the top row, both $r_{\mathrm{1,Ly\alpha}}$ and
$r_{\mathrm{1,UV}}$ show nearly monotonic behavior. While the power-law index
$\alpha$ also shows consistent behavior, the power-law functions often deviate
from the real data beyond a few tens of pkpc in Figure 11. On the other hand,
$r_{\mathrm{2,Ly\alpha}}$ behaves not that simply. This would result from both
astrophysical and observational reasons as we discuss in Section 5.2.

Figure 14 compares scale-lengths of UV/Ly$\alpha$ core/halo components. In the
left panel, the scale-lengths for the core components are compared. The
centroids of UV continuum emission of LAEs are known to show some offset from
those of Ly$\alpha$ emission which are defined as the image centers in this
work. The vast majority have offset lower than $0.2$ arcsec (Shibuya et al.,
2014; Leclercq et al., 2017) which is comparable to HSC’s pixel scale of 0.17
arcsec. This should not affect the measurement of $r_{\mathrm{2,UV}}$ but
leads to an overestimation of $r_{\mathrm{1,UV}}$. In addition, some
subsamples contain g-band nondetected sources (see Table 1) and could be
affected by oversubtraction. Although one should be aware of these potential
issues, $r_{\mathrm{1,Ly\alpha}}$ seems to be almost always larger than
$r_{\mathrm{1,UV}}$ and they correlate with Spearman’s rank correlation
coefficient of 0.74 and p-value of $8.2\times 10^{-5}$. On the other hand,
similarly to the trend seen in the panels in the middle row of Figure 13,
$r_{\mathrm{2,Ly\alpha}}$ do not have a clear correlation with
$r_{\mathrm{1,UV}}$ nor $r_{\mathrm{1,Ly\alpha}}$ ($p$-value 0.61 and 0.65,
respectively).

Table 2: The result of the fitting of one-component exponential functions to
UV profiles.

criteria | $C_{\mathrm{1,UV}}$ | $r_{\mathrm{1,UV}}$
---|---|---
(1) | (2) | (3)
All | $8.146\pm 0.072$ | $1.486\pm 0.019$
$M_{\mathrm{UV}}<-19.2$ | $26.830\pm 0.144$ | $1.632\pm 0.013$
$-19.2<M_{\mathrm{UV}}<-18.6$ | $13.896\pm 0.178$ | $1.285\pm 0.027$
$-18.6<M_{\mathrm{UV}}<-18.0$ | $8.049\pm 0.179$ | $1.203\pm 0.046$
$-18.0<M_{\mathrm{UV}}<-17.4$ | $5.552\pm 0.309$ | $0.495\pm 0.214$
$-17.4<M_{\mathrm{UV}}$ | – | –
$42.25<\log L_{\mathrm{Ly\alpha}}\mathrm{[erg~{}s^{-1}]}$ | $25.078\pm 0.157$ | $1.538\pm 0.014$
$42.05<\log L_{\mathrm{Ly\alpha}}<42.25$ | $13.285\pm 0.164$ | $1.278\pm 0.026$
$41.95<\log L_{\mathrm{Ly\alpha}}<42.05$ | $8.606\pm 0.210$ | $1.147\pm 0.050$
$41.85<\log L_{\mathrm{Ly\alpha}}<41.95$ | $6.118\pm 0.220$ | $1.042\pm 0.073$
$\log L_{\mathrm{Ly\alpha}}<41.85$ | – | –
$12<\mathrm{EW_{0,Ly\alpha}}<30$ Å | $21.208\pm 0.142$ | $1.710\pm 0.016$
$30<\mathrm{EW_{0,Ly\alpha}}<55$ Å | $11.772\pm 0.141$ | $1.614\pm 0.028$
$55<\mathrm{EW_{0,Ly\alpha}}<90$ Å | $8.182\pm 0.170$ | $1.353\pm 0.045$
$90<\mathrm{EW_{0,Ly\alpha}}<160$ Å | $6.173\pm 0.264$ | $0.690\pm 0.102$
160 Å$<\mathrm{EW_{0,Ly\alpha}}$ | – | –
$d_{\mathrm{Q}}<6.2$ pMpc | $8.451\pm 0.168$ | $1.395\pm 0.043$
$6.2<d_{\mathrm{Q}}<9.5$ pMpc | $8.716\pm 0.154$ | $1.458\pm 0.039$
$9.5<d_{\mathrm{Q}}<12.0$ pMpc | $7.798\pm 0.163$ | $1.495\pm 0.047$
$12.0<d_{\mathrm{Q}}<14.8$ pMpc | $8.229\pm 0.149$ | $1.475\pm 0.040$
$14.8pMpc<d_{\mathrm{Q}}<16.9$ pMpc | $7.512\pm 0.148$ | $1.625\pm 0.046$
$2.5<\delta$ | $9.423\pm 0.471$ | $1.815\pm 0.123$
$1.0<\delta<2.5$ | $8.319\pm 0.212$ | $1.391\pm 0.054$
$0.3<\delta<1.0$ | $8.050\pm 0.143$ | $1.445\pm 0.038$
$-0.15<\delta<0.3$ | $7.828\pm 0.129$ | $1.507\pm 0.036$
$-1.0<\delta<-0.15$ | $8.431\pm 0.133$ | $1.530\pm 0.034$

Note. — Column (1): criteria used to define subsamples, Column (2):
$C_{\mathrm{1,UV}}$ in units of
$10^{-31}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}Hz^{-1}~{}arcsec^{-2}}$, Column
(3): $r_{\mathrm{1,UV}}$ in units of physical kpc. The uncertainties of
fitting parameters include fitting errors only.

Table 3: The result of the fitting of two-component exponential functions to
UV profiles.

criteria |  $C_{\mathrm{1,UV}}$, $C_{\mathrm{2,UV}}$ |  $r_{\mathrm{1,UV}}$, $r_{\mathrm{2,UV}}$
---|---|---
(1) | (2) | (3)
$M_{\mathrm{UV}}<-19.2$ | $26.730\pm 0.389$ | $1.280\pm 0.048$
| $1.969\pm 0.495$ | $4.845\pm 0.532$
$42.25<\log L_{\mathrm{Ly\alpha}}\mathrm{[erg~{}s^{-1}]}$ | $22.466\pm 1.150$ | $0.874\pm 0.110$
| $5.483\pm 0.495$ | $3.072\pm 0.289$
$12<\mathrm{EW_{0,Ly\alpha}}<30$ Å | $20.804\pm 0.531$ | $1.343\pm 0.071$
| $1.818\pm 0.660$ | $4.570\pm 0.679$
$30<\mathrm{EW_{0,Ly\alpha}}<55$ Å | $12.102\pm 0.169$ | $1.407\pm 0.050$
| $0.250\pm 0.080$ | $10.586\pm 2.201$

Note. — Column (1): criteria used to define subsamples, Column (2): First and
second rows respectively show $C_{\mathrm{1,UV}}$ and $C_{\mathrm{2,UV}}$ in
units of $10^{-31}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}Hz^{-1}~{}arcsec^{-2}}$,
Column (3): First and second rows respectively show $r_{\mathrm{1,UV}}$ and
$r_{\mathrm{2,UV}}$ in units of physical kpc. The uncertainties of fitting
parameters include fitting errors only.

Figure 13:  Resulting scale-lengths and power-law index vs. median of each
subsample. From left to right, we show the results for UV magnitude,
Ly$\alpha$ luminosity, Ly$\alpha$ equivalent width, distance from the HLQSO,
and environment subsamples. $r_{\mathrm{1,UV}}$ are those obtained from a one-
component fit, while $r_{\mathrm{1,Ly\alpha}}$ are from a two-component fit.
Figure 14:  Comparison between scale-length for the UV first component
$r_{\mathrm{1,UV}}$ and the Ly$\alpha$ first and second component
$r_{\mathrm{1,Ly\alpha}},r_{\mathrm{2,Ly\alpha}}$. $r_{\mathrm{1,UV}}$ are
those obtained from a one-component fit, while $r_{\mathrm{1,Ly\alpha}}$ are
from a two-component fit. In each panel, we list Spearman’s rank correlation
coefficient $\rho$ and p-value $p$. The diagonal line in the left panel shows
the 1:1 relation.

##  5 Discussion

###  5.1 Sources of Differences from Previous Observational Studies

While we clearly detect LAHs more extended than UV continuum for all subsets
of LAEs, some previous works report non-detections of such components (e.g.,
Bond et al., 2010; Feldmeier et al., 2013). We describe a number of possible
reasons for such discrepancies between this study and others. The Ly$\alpha$
morphology cannot be properly captured by just comparing simple quantities
such as their FWHMs or half-light radii without enough sensitivity or taking a
very large aperture for total luminosity estimate (Nilsson et al., 2009b).
Detailed analyses on SB profiles are desirable but with enough sensitivity of
(at least) $\sim 10^{-19}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$,
and even higher sensitivity is required for safer arguments beyond a mere
detection. The scale-lengths of an exponential function(s) are most widely
used in the literature for such analyses (e.g., Steidel et al., 2011; Matsuda
et al., 2012; Momose et al., 2014, 2016; Wisotzki et al., 2016; Leclercq et
al., 2017; Xue et al., 2017). However, resulting scale-lengths vary a lot
depending not only on the data quality but also on the details of analyses:
results depend on sample selection criteria, sky subtraction, masking, and
stacking methods, the range of radius used for fitting, the radial binning
size, fitting functions, whether to assume
$r_{\mathrm{1,UV}}=r_{\mathrm{1,Ly\alpha}}$ or not, etc. First, we
investigated the impact of the sensitivity using randomly chosen LAEs with a
smaller sample size. We randomly took 100 LAEs and obtained
$r_{\mathrm{2,Ly\alpha}}$ of stacked Ly$\alpha$ images, and repeated this
process 1000 times. While the medians of obtained distribution of
$r_{\mathrm{2,Ly\alpha}}$ did not show a systematic trend, the distribution
had a 3.4 times larger standard deviation compared to that obtained for the
700 LAE case (see Section 5.2.1), although this number may just represent
diversity in our sample. Results from lower-sensitivity data could be thus
more uncertain. Secondly, most of the observed Ly$\alpha$ SB profiles are
downwardly convex (Figure 11) due to the flattening at $\sim 15$ pkpc and thus
the scale-length becomes smaller when the outer (inner) boundary of the
fitting range is smaller (larger). Indeed, if we limit our fitting range up to
$<30$ pkpc, obtained $r_{\mathrm{2,Ly\alpha}}$ is underestimated by 35%. The
outer boundary is also affected by the sensitivity; without deep data, one has
no choice other than to set it to smaller values where the signal is detected.
Thirdly, as we showed in Figure 8, an insufficiently small mesh size for local
sky background estimate leads to underestimation of the scale-length of LAHs,
but in many cases the mesh size (or even a brief summary about sky
subtraction) has not been presented in the text. As for fitting functions, a
two-component fit can more robustly capture the shape of the Ly$\alpha$ SB
profiles (see Appendix C of Xue et al., 2017). However, there are few cases
where such analyses have been present with enough sensitivity at $z\sim 3$.

For example, Xue et al. (2017) reported a halo scale-length of LAEs of $\sim
5$–9 pkpc and found no evidence for environmental dependence based on NB
surveys of two protoclusters, at $z=3.78$ and $z=2.66$. These observations are
an order of magnitude shallower than the present study; Xue et al. (2017) used
the Mayall 4m telescope for the $z=3.78$ data, and the Subaru telescope for
the $z=2.66$ data but used an intermediate band filter (IA445 on Suprime-Cam,
$\Delta\lambda=201$Å) which has lower line sensitivity. The number of LAEs
used to examine environmental dependence was at most $139$ (for an
intermediate density sample). The criteria used to select LAEs picks up
relatively high EW LAEs, with EW${}_{\mathrm{0,Ly\alpha}}>50$Å, in the
$z=2.66$ protocluster field. Our study suggests that this could have biased
the results toward smaller LAHs (Figure 19). Lower ionizing radiation field
strength and/or lower abundance of cool gas in their protoclusters may also
produce smaller LAHs (see Section 5.2.1 and 5.2.2). As for observations with
integral field spectrographs, Wisotzki et al. (2016); Leclercq et al. (2017)
probed $r<30$ pkpc of LAEs at $z>3$ and obtained $r_{\mathrm{2,Ly\alpha}}$ of
$\lesssim 15$ pkpc with the majority with $r_{\mathrm{2,Ly\alpha}}<5$ pkpc.
Analyses on an individual basis would suffer from greater noise and sample
variance than ours. On the other hand, Chen et al. (2021) probed 59 star-
forming galaxies (including non-LAEs) at $z=2$–3 with sufficiently deep
Keck/KCWI observational data and conducted stacking. They got
$r_{\mathrm{1,Ly\alpha}}=3.71^{+0.06}_{-0.04}$ pkpc and
$r_{\mathrm{2,Ly\alpha}}=15.6^{+0.5}_{-0.4}$ pkpc for their stacked Ly$\alpha$
SB profiles. Their sample is typically about an order of magnitude more
massive and star-forming than our sample and this could explain the size
difference. To summarize, comparing results obtained using inhomogeneous
analyses in the literature without consistent reanalysis is difficult, and
thus we will not attempt a detailed comparison here.

###  5.2 Dependence of Scale-length on Galaxy Properties

In Figure 13, UV/Ly$\alpha$-bright or low-EW LAEs tend to have larger
$r_{1,Ly\alpha}$ and $r_{1,UV}$ than UV/Ly$\alpha$-faint or high-EW LAEs. This
is qualitatively consistent with the HST-based results of Leclercq et al.
(2017) although our results from a ground-based telescope tend to show larger
values (0.1–1 pkpc vs $>1$ pkpc). Considering that UV/Ly$\alpha$ luminous
and/or low-EW LAEs tend to be more massive, the trend is also consistent with
the known trend between $M_{\mathrm{UV}}$ and effective radius in the UV
(e.g., Shibuya et al., 2019), or more generally the so-called size-luminosity
or size-mass relation. A larger scale-length in massive LAEs also results from
suppression of UV/Ly$\alpha$ light due to more abundant dust especially in the
central region (Laursen et al., 2009) than in less massive LAEs. Ly$\alpha$
photons are then further affected by resonant scattering and differential
extinction caused thereby, leading to
$r_{\mathrm{1,Ly\alpha}}>r_{\mathrm{1,UV}}$.

On the other hand, $r_{\mathrm{2,Ly\alpha}}$ do not show a simple trend with
respect to any photometric properties. This is again almost consistent with
Leclercq et al. (2017). This fact can be attributed to both observational and
astrophysical reasons. First, we may simply lack the sensitivity to reveal a
real trend. Even with our deep images, the Ly$\alpha$ SB profiles in Figure 11
at $r>50$ pkpc are not well constrained. In addition, the astrophysics
involved in diffuse emission in LAHs is notoriously complicated as we discuss
in Section 5.4. Dominant mechanisms for Ly$\alpha$ production might be
different over different mass/luminosity ranges, making a simple trend (if
any) difficult to discern. For example, as compiled in Kusakabe et al. (2019,
their Figure 7), the total Ly$\alpha$ luminosity of the halo component may
depend in a different way on halo mass with respect to production mechanisms
(e.g., collisional excitation in cold streams vs. scattering). Lastly, both
observations (Wisotzki et al., 2016; Leclercq et al., 2017) and numerical
studies (Lake et al., 2015; Byrohl et al., 2021) have shown that the
Ly$\alpha$ SB profiles of individual LAEs are very diverse, even within
galaxies with similar integrated properties. To conclude whether there is any
trend, even larger samples and/or deeper observations are needed. With better
data, we could more easily select whether two-component exponential functions
or the power-law functions are preferred.

####  5.2.1 Curious Behavior of Distance Subsamples: QSO Radiative History
Imprinted?

Subsamples based on the projected distance from the HLQSO ($d_{\mathrm{Q}}$)
show a significant variation with a minimum at $d_{\mathrm{Q}}\sim 10$ pMpc
(Figure 13, the second panel from right in the second row). This could just be
due to the stochasticity discussed above, but if real, it could be related to
the QSO’s radiative history. $r_{\mathrm{2,Ly\alpha}}$ of the two largest
($r_{\mathrm{2,Ly\alpha}}>13$ pkpc, $d_{\mathrm{Q}}<6.2$ pMpc and 14.8
$<d_{\mathrm{Q}}<16.9$ pMpc subsample) and the smallest
($r_{\mathrm{2,Ly\alpha}}<7$ pkpc, $9.5<d_{\mathrm{Q}}<12$ pMpc subsample)
differ significantly; when we repeatedly select 700 LAEs at random from the
whole sample, stacking their Ly$\alpha$ images, and measuring
$r_{\mathrm{2,Ly\alpha}}$ 1000 times, we get $r_{\mathrm{2,Ly\alpha}}<7$ pkpc
0.8% of the time and $r_{\mathrm{2,Ly\alpha}}>13$ pkpc 12% of the time (with
the median value of $r_{\mathrm{2,Ly\alpha}}=10.5$ pkpc) 888If we exclude the
55 protocluster LAEs from the $d_{\mathrm{Q}}<6.2$ pMpc subsample,
$r_{\mathrm{2,Ly\alpha}}$ becomes 11.3 pkpc (the $\sim 30$ percentile).. If
the HLQSO was active $\sim 50$ Myrs ago, followed by $\sim 30$ Myrs of
inactivity, and was re-ignited $\lesssim 20$ Myrs ago, the ionizing photons
emitted by the QSO would have traveled a distance of $>15$ pMpc and $<6$ pMpc
from the QSO 999These time estimates are lower limits calculated with
projected distance and the speed of light. Propagation of ionization fronts
could be delayed in some situations (Shapiro et al., 2006).. These photons can
ionize the envelopes of the LAEs and boost their Ly$\alpha$ luminosity,
explaining the observed behavior 101010see also Trainor & Steidel (2013);
Borisova et al. (2016) where the authors used QSOs associated with
spectroscopic high-EW($>240$ Å) LAEs to place limits on QSO lifetime. We see
consistent behavior also in the fraction of high EW LAEs as a function of
distance from the HLQSO (Kida et al., 2019). Those EW are derived with
1$.\\!\\!\arcsec$5 aperture. However, boosting EW of such a central part of
LAEs with $L_{\mathrm{Ly\alpha}}>10^{41}\mathrm{erg~{}s^{-1}}$ at 16 pMpc
distance would be energetically not feasible with the current HLQSO
luminosity.. Assuming the HLQSO 50 Myrs ago had the same luminosity as it has
today ( luminosity near a rest-frame wavelength 1450 Å, $\nu
L_{\mathrm{\nu,1450}}=5.7\times 10^{47}\mathrm{~{}erg~{}s^{-1}}$, Trainor &
Steidel 2012) and isotropic radiation with escape fraction of unity, the
ionizing radiation at 16 pMpc from the QSO can still dominate over the cosmic
average UV background at $z\sim 3$,
$\Gamma_{\mathrm{bkg}}^{\mathrm{z=3}}=1.0\times 10^{-12}\mathrm{~{}s^{-1}}$
(Becker & Bolton, 2013) by a factor of a few. Cantalupo et al. (2005)
calculated the fluorescent Ly$\alpha$ emission due to QSOs in addition to the
cosmic background and gave a fitting formula for an effective boost factor
$b_{\mathrm{eff}}$ (their Equation 14–16) which can be used to estimate SB of
illuminated gas clouds. In our case, the resulting SB would be
$\mathrm{SB}=(0.74+0.50(11.5(r/16\mathrm{~{}pMpc})^{-2}(1.6^{\alpha}/\alpha))^{0.89})\mathrm{SB_{HM}}$,
where $r$ is the distance from the HLQSO in pMpc, $\alpha$ is the QSO’s
spectral slope ($L_{\nu}\propto\nu^{-\alpha}$), and
$\mathrm{SB_{HM}}=3.67\times
10^{-20}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ is the expected SB
without QSO boost. Assuming $\alpha\sim 1$ gives $\mathrm{SB}=2.7\times
10^{-19}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$. Thus, QSO-induced
fluorescence is energetically possible to have caused variation in
$r_{\mathrm{2,Ly\alpha}}$ in the projected distance subsamples, at least in
the optimistic case. In reality, our narrow-band selection picks up LAEs with
line-of-sight distance uncertainty of $\sim 19$ pMpc, which is the same level
as the radius of the FoV of our images, and this randomizes light-travel time
from the QSO to each LAE. Such effects further complicate the situation, but
we have shown that under some circumstances with appropriate QSO light curve
and line of sight distribution of LAEs, the observed trend of
$r_{\mathrm{2,Ly\alpha}}$ might be explained. Upcoming instruments such as
Prime Focus Spectrograph (PFS; a wide-field multi-fiber spectrograph, Takada
et al., 2014) on the Subaru Telescope can test this fluorescence scenario by
obtaining systemic redshifts of LAEs and thus reducing uncertainties on their
3D distances from the QSO.

####  5.2.2 Origin of the Large LAH of Protocluster LAEs

We showed that the dependence of LAHs on environment is not large where
$\delta<2.5$, but in the protocluster ($\delta>2.5$) the LAHs show elevated
flux out to $>100$ pkpc (Figure 13, 20). If LAHs trace the gas distribution of
the CGM, the former indicates that the large-scale environment (except for
protocluster environments) does not have a large impact on the matter
distribution out to $\sim 100$ pkpc and it is determined by individual halo
mass or other internal processes. Alternatively, LAEs could simply be poor
tracers of large-scale environments. Recently, Momose et al. (2021) found that
LAEs behind the foreground large-scale structure tend to be missed due to
absorption by the foreground structure (see also Shimakawa et al., 2017). It
may also be the case that even our large field of view of 1.2 deg diameter may
be insufficient to capture diverse environments including voids while
targeting a single protocluster. Other line emitters or continuum-selected
galaxies would be ideal, though would be expensive for current facilities. On
the other hand, the HS1549 protocluster is confirmed by overdensity of LAEs
and continuum-selected galaxies, with $\sim 200$ member galaxies
spectroscopically identified (Trainor & Steidel 2012; Mostardi et al. 2013; C.
Steidel et al. in prep.). Enlarging the size of cutout images, we confirmed
that flux higher than $>1\sigma$ level continues out to $\sim 500$ pkpc. In
the following, we investigate the possible cause of this emission in the
$\delta>2.5$ subsample.

First, we further divide the protocluster sample into a “core” group (LAEs
within a projected distance of $<500$ pkpc from the HLQSO but excluding the
HLQSO itself, $N=25$) and an “outskirt” group (the remainder, $N=29$), and
stacked them separately. This time, pixels with
SB$>10^{-17}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ around the HLQSO
and the associated bright nebula (Kikuta et al., 2019) are masked before
stacking in each cutout image to exclude their contribution. As shown in
Figure 15, the core sample (the yellow curve) clearly shows excess emission
which does not decrease toward the edge above the original orange curve, while
the outskirt stack (the green curve) shows just a mild bump around $\sim 15$
pkpc and no evidence of excess. This suggests the extended emission of
$\delta>2.5$ sample is solely produced by the LAEs at the core of the
protocluster.

Recalling that the protocluster LAEs tend to have higher
$M_{\mathrm{UV}},~{}L_{\mathrm{Ly\alpha}}$ and lower EW (Figure 2) and that
such galaxies generally have more extended Ly$\alpha$ SB profiles (Figure 13
and Figure 19), the overabundance of such LAEs leads to larger LAHs. Moreover,
the core region is sufficiently crowded that Ly$\alpha$ emission from
neighboring LAEs may overlap, thereby leading to an overestimation of the SB
profile. We evaluated this effect using the stacked Ly$\alpha$ images of
subsamples based on UV magnitude (those shown in the top row of Figure 10) as
follows. We embedded their images in a blank image, mimicking the observed
spatial distribution of LAEs (including the HLQSO) using IRAF task
“mkobjects”. When embedding LAEs with a certain $\mathrm{M_{UV}}$, we used the
stacked Ly$\alpha$ image of the appropriate $\mathrm{M_{UV}}$ subsample scaled
to match the observed UV magnitude. For example, we embedded the scaled
stacked image of the $-19.2<\mathrm{M_{UV}<-18.6}$ subsample at the location
of LAEs in the same UV magnitude range. Cutout images of the simulated
Ly$\alpha$ image were then created at the locations of the embedded LAEs and
stacked in the same manner as the real “core” LAEs. The result is also shown
in Figure 15 (the blue curve). The observed large LAH (the orange curve) is
not reproduced. Thus, we conclude that the combined effect of the
overabundance of bright LAEs and overlapping do contribute to the large LAHs,
but cannot fully explain the extent of the LAH of the core LAEs.

Kikuta et al. (2019) showed that there is a Mpc-scale diffuse Ly$\alpha$
emitting structure around the HLQSO. Such diffuse emission would explain the
remaining excess in the core region of the protocluster. The excess comes not
only from gas directly associated with LAEs but also from gas out of the LAEs;
the pixel value distribution of the core region ($<500$ pkpc from the HLQSO)
after masking $<10$ arcsec ($\sim 80$ pkpc) regions around the detected LAEs
and a $\sim 200$ pkpc $\times 260$ pkpc box covering the bright QSO nebula
clearly shows excess at $1\times 10^{-18}<$SB$<4\times
10^{-18}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ compared to that of
the outer region (with regions around LAEs masked in the same manner as the
core case). Previously, two studies reported similar large diffuse LAHs in
protoclusters at $z=2$–3; Steidel et al. (2011) observed HS1549 ($z=2.84$,
same as this study 111111Although their LBG sample is an order of magnitude
brighter in UV than our LAE sample, we did a consistency check with Steidel et
al. (2011) results; we confirmed that our stack of LBGs used also in their
sample gives a consistent Ly$\alpha$ SB profile as their work.), HS1700
($z=2.30$; see also Erb et al., 2011), and SSA22 ($z=3.09$; also observed by
Matsuda et al., 2012). For the HS1549 and SSA22, direct evidence of Mpc-scale
diffuse Ly$\alpha$ emission has been observed (Kikuta et al. 2019 and Umehata
et al. 2019, respectively), and the identification of filamentary structure
traced by six Ly$\alpha$ blobs in the HS1700 protocluster by Erb et al. (2011)
suggests that this protocluster also harbors such diffuse emission. In the
forming protocluster core at $z=2$–3, a large amount of cold gas can be
accreted through the filamentary structure (the cosmic web) penetrating the
cores (Kereš et al., 2005, 2009; Dekel et al., 2009). In addition to abundant
gas, Umehata et al. (2019) showed that an enhanced ionizing UV background due
to a local overdensity of star-forming galaxies and AGNs may play a crucial
role in boosting fluorescent Ly$\alpha$ emission to a detectable level. The
HS1549 protocluster also has tens of active sources (i.e., AGN, Ly$\alpha$
blobs (Kikuta et al., 2019, C. Steidel et al. in prep.), and SMGs (Lacaille et
al., 2019)) within a few arcmins from the HLQSO, producing sufficient UV
radiation to power the diffuse emission; in Supplementary Material S9 of
Umehata et al. (2019), they calculated the required UVB strength to boost SB
of optically thick gas. Fifty times stronger UVB than the cosmic average at
$z=3$, which is easily realized by the HLQSO alone in the area within several
pMpc from it, would boost SB to $\sim
10^{-18}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ level. To summarize,
the very large LAHs of stacked Ly$\alpha$ profiles reported previously and in
this work can be attributed to an overlap of crowded LAEs and diffuse
fluorescent Ly$\alpha$ emission within the forming protocluster core. The
$\delta^{2}$ dependence of LAH scale-length claimed in Matsuda et al. (2012)
should be revisited using new data targeting more protocluster fields at
similar redshift together with appropriate analysis methods as discussed in
Section 5.1.

Figure 15: Stacked Ly$\alpha$ SB profiles of the protocluster ($\delta>2.5$)
sample (orange), the core sample (yellow), the outskirt sample (green),
embedded and stacked image (blue, see text), and all LAEs (purple, same as
Figure 7). Downward triangles show 1$\sigma$ error levels after residual sky
subtraction.

###  5.3 Discovery of “UV Halos” and Its Implication to Low-Mass Galaxy
Evolution

As seen in Figure 12 and Figure 13 in Section 4.2, we have discovered “UV
halos” around UV/Ly$\alpha$-bright and/or low-EW LAEs. This has a significant
impact on our understanding of the origin of LAHs, and also of galaxy
evolution, because it provides direct evidence of star formation activity in
the outskirts of high-redshift low-mass galaxies.

To gain more insight on the latter point, we used the data products from the
TNG100 run of the IllustrisTNG simulation (e.g., Nelson et al., 2019;
Pillepich et al., 2018). We make median stacked SFR surface density profiles
of FOF (friends-of-friends) halos, which roughly represent collections of
gravitationally bound DM particles, at $z=3$ for 4 SFR bins;
$0.1<\mathrm{SFR}<1~{}\mathrm{M_{\odot}~{}yr^{-1}}$,
$1<\mathrm{SFR}<10~{}\mathrm{M_{\odot}~{}yr^{-1}}$,
$10<\mathrm{SFR}<100~{}\mathrm{M_{\odot}~{}yr^{-1}}$, and
$100~{}\mathrm{M_{\odot}~{}yr^{-1}}<\mathrm{SFR}$ and compare them with those
of our three UV-brightest subsamples after converting the simulation data to
UV flux density using the SFR-UV luminosity density conversion of Murphy et
al. (2011); Kennicutt & Evans (2012) relation (Figure 16). The SFR surface
density profiles were convolved with a Gaussian with 0$.\\!\\!\arcsec$77
FWHM121212As we see in Section 3, the PSF of our images are not Gaussian but
have a power-law tail with an index of $\sim-3$, but the PSF behavior at large
radii does not affect our result. to make a fair comparison. Here, SFR of TNG
galaxies denotes a total of all particles which belong to one FOF halo. We
also plot the prediction of Lake et al. (2015, convolved with a Gaussian with
1$.\\!\\!\arcsec$32 FWHM or 10.3 pkpc at $z=3.1$) for discussion in Section
5.4. The three UV-bright subsamples have median UV absolute magnitude of
$-19.62,-18.88,-18.31$ (Table 1), but these were derived using
1$.\\!\\!\arcsec$5 diameter apertures and thus may be underestimated. The
total magnitudes derived by integrating the UV SB profiles down to the radii
where emission is detected at more than $1\sigma$ significance are
$-21.34,-20.49,-20.10$, respectively, corresponding to SFR of
$17,7.7,5.4~{}\mathrm{M_{\odot}~{}yr^{-1}}$. The SFR of simulated galaxies
whose UV SB profiles match those of our LAEs is higher than that of our LAEs
(10–100 vs. 5.4–17), but there is considerable uncertainty in the conversion
between UV luminosity and SFR, as it depends on stellar age, dust attenuation,
and metal abundance, and not all galaxies in simulations would be selected as
LAEs. Reconciling this mismatch is beyond the scope of our paper. Rather, to
compare the profile shapes we normalized each curve in the bottom panel of
Figure 16. While the UV SB profiles of the two fainter subsamples seem to be
slightly more compact than the SFR density profiles of TNG galaxies, the UV-
brightest subsample has a remarkably similar shape as the SFR surface density
profiles of TNG galaxies with $1<\mathrm{SFR}<10$ and $10<\mathrm{SFR}<100$
subsamples.

Figure 16:  (Top) SFR surface density profiles of TNG100 galaxies converted to
UV SB profiles (thick colored solid curves, see text). Different colors show
different SFR bins. Dashed and dotted curves show contributions to the
profiles from the main halo and other halos, respectively. The thin purple
curve shows that predicted in Lake et al. (2015). Gray curves with error bars
are observed UV SB profiles of our three UV-brightest LAE subsamples. These
profiles are corrected for a slight redshift dimming effect by scaling
$(1+z)^{3}/(1+2.845)^{3}$. (Bottom) Curves of top panels are normalized to
better compare their shapes.

We further decompose the simulated FOF halos into the main halo and subhalos.
The decomposed SFR profiles demonstrate that the flattened outer part is
dominated by the contribution of satellites. Given the similarity of the
profiles, we suggest that the UV halo of the UV-brightest LAEs is also due to
such satellites. To characterize the subhalos around the central galaxy, we
extracted dark matter halo mass, gas mass, stellar mass, and SFR of those
within 50 pkpc ($\sim 6$ arcsec, 2D projected distance) from their main halos.
We only handle those halos with $M_{\mathrm{DM}}>7.5\times
10^{7}~{}\mathrm{M_{\odot}}$ ($=10\times$DM particle mass),
$M_{\mathrm{stellar}}>7.0\times 10^{6}~{}\mathrm{M_{\odot}}$
($=5\times$stellar particle mass), and non-zero SFR to avoid spurious objects.
The distribution of DM mass, stellar mass, gas mass, and SFR of satellites are
plotted in Figure 17. Satellites responsible for the UV halo of our LAEs would
be similar to those around central galaxies with $1<\mathrm{SFR}<10$ and
$10<\mathrm{SFR}<100$. On average, they have 1.9 and 2.3 satellites,
respectively, with median DM halo masses of $3.3\times
10^{9}~{}\mathrm{M_{\odot}}$ and $4.4\times 10^{9}~{}\mathrm{M_{\odot}}$, and
mean total halo SFR of 0.30 and 2.6 $\mathrm{M_{\odot}}~{}\mathrm{yr^{-1}}$
(i.e., $\sim 10\%$ of central galaxies). This suggests that the UV halo is
comprised of a few satellite galaxies around the main halo and not by an
intrinsically diffuse halo, unlike optical stellar halos of local galaxies
(e.g., D’Souza et al., 2014). Under this hypothesis, individual LAEs would not
have smooth UV SB profiles like those presented in Figure 12, but would be
more likely to exhibit stochastic shapes made by $\sim 2$ discrete satellites;
thus the term UV “halo” may not be appropriate, if the quoted simulations
represent the reality. The smooth profile of stacked UV SB is simply a
reflection of the radial and SFR distribution of satellite galaxies.

Figure 17:  The distribution of DM mass, stellar mass, gas mass, and SFR (from
left to right) of satellite halos for different main halo’s SFR bins, the
corresponding colors to which are shown in the third panel. Small colored
ticks below each panel indicate the medians of each bin.

Are such satellites observable? The SFR distributions of satellites around
galaxies with $0.1<\mathrm{SFR}<10$ in the simulations have medians in the
range $0.01<\mathrm{SFR/[M_{\odot}~{}\mathrm{yr^{-1}}]}<0.1$, although some
very low-mass objects may be affected by resolution effects (Figure 17). This
translates into UV absolute magnitude of $-13.3$ and $-15.8$, and apparent
g-band magnitude of 32.1 to 29.6 mag at $z=2.84$ (assuming no K-correction, or
equivalently, flat UV continua). Brighter satellites are thus well within the
reach of HST and JWST sensitivity in some deep fields or with the aid of
gravitational lensing (e.g., Alavi et al., 2016; Bouwens et al., 2022).

###  5.4 On the Origin of the LAHs

Lastly, we infer the origin of LAHs of LAEs through comparison with recent
numerical simulations. As introduced in Section 1, Ly$\alpha$ photons in the
halo regions are generated either by ex-situ (mostly from the host galaxy) and
transported by resonant scattering in neutral gas, or in-situ via
photoionization followed by recombination or collisional excitation. In-situ
photoionization is maintained by ionizing photons from star formation and/or
AGN activity in satellites, central galaxies, or other nearby sources of
ionizing UV such as QSOs. In-situ collisional excitation would be driven by
shocks due to fast outflows due to feedback or gravitational energy of
inflowing gas, but the former is considered to be effective in more massive
and energetic sources such as radio galaxies and Ly$\alpha$ blobs (e.g., Mori
et al., 2004). In-situ Ly$\alpha$ photons may also experience resonant
scattering, but its effect on the redistribution of photons is likely to be
relatively weak due to lower HI column densities than the central part by more
than a few dex (Hummels et al., 2019; van de Voort et al., 2019).

To predict Ly$\alpha$ SB profiles around galaxies, one needs to know physical
parameters such as hydrogen density, neutral fraction, gas kinematics, and
temperature which depend on SFR and AGN activity around a point of interest,
and ionizing/Ly$\alpha$ photon escape fraction, etc. Solving all these
quantities is practically impossible, but recent simulations are beginning to
reproduce observed Ly$\alpha$ SB profiles reasonably well. Among such studies
is Byrohl et al. (2021, hereafter B21); B21 presents full radiative transfer
calculations via post-processing of thousands of galaxies in the stellar mass
range $8.0<\log(M_{*}/M_{\odot})<10.5$ drawn from the TNG50 simulations of the
IllustrisTNG project. One of the advantages of B21 is the sample size, which
is far larger than those of previous studies. For example, Lake et al. (2015,
hereafter L15) performed radiative transfer modeling of 9 LAEs, obtaining
significantly diverse SB profiles. Those with massive neighbors have elevated
SB profiles both in UV and Ly$\alpha$, significantly boosting the average
profile. But the 9 galaxies modeled in the simulation would not be
representative of star-forming galaxies at $z\sim 3$ if not carefully
selected. Mitchell et al. (2021) calculated Ly$\alpha$ SB profiles of a single
galaxy at $z=3$–4 using the RAMSES-RT code, a radiation hydrodynamics
extension of the RAMSES code, and succeeded in reproducing SB profiles similar
to MUSE observations. The small sample size is somewhat mitigated by using all
available outputs between $z=4$ and $z=3$, but still there may remain biases
with respect to environment or evolutionary phase. For these reasons, we
compare our results primarily with those of B21.

####  5.4.1 Dominance of Star Formation in Central and Satellite Galaxies

In Figure 18, we plot B21 (convolved with a Gaussian with 0$.\\!\\!\arcsec$7
FWHM) and L15 (convolved with a Gaussian with 1$.\\!\\!\arcsec$32 FWHM)
results with our observations of subsamples based on UV magnitude. Although
B21 probed only up to $<50$ pkpc, the overall shapes of the
$9.5<\log(M_{*}/M_{\odot})<10.0$ stack and our UV-brightest subsample stack
match remarkably well, and the $9.0<\log(M_{*}/M_{\odot})<9.5$ stack and the
other UV subsamples stack also show fairly good agreement. We also highlight
the diversity of Ly$\alpha$ SB profiles here by drawing all the profiles of
galaxies with $9.0<\log(M_{*}/M_{\odot})<9.5$ in B21 in Figure 18 with thin
curves. Similarly to L15, bumps in some curves are caused by companion
galaxies. This demonstrates the difficulty of studying halo origins with a
small sample as discussed in Section 5.2. Future (observational and
theoretical) studies should keep this in mind before discussing the halo
dependence on physical properties. B21 concluded that scattering of Ly$\alpha$
photons originate from star formation in the central and satellite galaxies is
almost always dominant within 50 pkpc from the center ($\sim 50$% at $>20$
pkpc), with in-situ collisions and recombination explaining the remaining 30%
and 20%, respectively. The contribution from satellite star formation
dominates over that from the central galaxy beyond $\sim 40$ pkpc, and they
show that halos that have more massive neighbors within 500 pkpc can have very
extended LAHs compared to those residing in normal environments (Figure 12 of
B21). Similar conclusions are reached by L15 and Mitchell et al. (2021) as to
the importance of satellites. The dashed curve shows the L15 result which
extends to 100 pkpc. As we saw in Figure 16 (see also Figure 4 of L15), their
galaxies (whose mean stellar mass is $\sim 2.9\times
10^{10}\mathrm{~{}M_{\odot}}$) have more star formation outside of the host
halos and have an enhanced SB profile at outer regions. A comparably massive
galaxy sample is required to confirm their prediction. Our first detection of
satellites (Section 5.3) and reasonable agreement of both UV (Figure 16) and
Ly$\alpha$ (Figure 18) SB profiles with simulations suggest that star
formation in central and satellite galaxies are important Ly$\alpha$ sources
contributing to LAHs.

Rest-frame Ly$\alpha$ EW of each annulus of observed LAEs is plotted with the
right axis of Figure 12. Considering the outward diffusion of Ly$\alpha$ from
the central galaxies, this EW is always an upper limit for the EW of in-situ
Ly$\alpha$ emission. If low-mass satellites are responsible for the outer LAHs
as simulations suggest, then their expected dark matter halo masses are about
$10^{9-10}\mathrm{~{}M_{\odot}}$. EW0,Lyα of $\gtrsim 200$ Å observed at
$r\sim 30$ pkpc of the UV/Ly$\alpha$-brightest or EW-lowest subsamples can be
explained by halo star formation alone if these low-mass galaxies have average
EW0,Lyα of $\geq 200$ Å, even without scattered Ly$\alpha$ from central
galaxies. For other subsamples, the EW0,Lyα is $>240$ Å due to the faintness
of UV SB and the extended Ly$\alpha$ SB profiles at outer regions. Such high
EW0,Lyα is hard to explain by star formation alone (Schaerer, 2003);
scattering of Ly$\alpha$ photons produced elsewhere and in-situ
recombination/collisional excitation should dominate the Ly$\alpha$ photon
budget.

Figure 18:  Stacked Ly$\alpha$ SB profiles of galaxies at $z=3$ for five
different intervals of stellar mass from the Byrohl et al. (2021) simulation
(thick gray) along with the averaged SB profile predicted in Lake et al.
(2015) (dashed) and our results (colored with errorbars, UV subsamples from
the top left panel of Figure 19). The thick gray curves represent, from bottom
to top, profiles of galaxies with $8.0<\log(M_{*}/\mathrm{M_{\odot}})<8.5$,
$8.5<\log(M_{*}/\mathrm{M_{\odot}})<9.0$,
$9.0<\log(M_{*}/\mathrm{M_{\odot}})<9.5$,
$9.5<\log(M_{*}/\mathrm{M_{\odot}})<10$, and
$10<\log(M_{*}/\mathrm{M_{\odot}})<10.5$. To show the diversity of Ly$\alpha$
SB profiles, we show all profiles of galaxies with
$9.0<\log(M_{*}/\mathrm{M_{\odot}})<9.5$ with thin gray curves. The B21 and
L15 results are corrected for a slight redshift dimming effect by scaling
$(1+z)^{4}/(1+2.845)^{4}$ ($z=3.0$ for B12 and 3.1 for L15).

####  5.4.2 Non-Negligible Contribution from In-Situ Ly$\alpha$ Production

Processes other than star formation, e.g. QSO-boosted fluorescence or
collisional excitation via gravitational cooling, can still be important at
large radii, since they are predicted to make non-negligible contributions and
there are situations under which such processes become more important. In
Section 5.2.1 and Section 5.2.2, we specified conditions where fluorescence
could be dominant, namely regions near bright QSOs and/or near protocluster
cores. A number of simulations have suggested that stronger ionizing radiation
fields boost fluorescent Ly$\alpha$ emission from the CGM/IGM (Cantalupo et
al., 2005; Kollmeier et al., 2010, see also Appendix A5 of B21). A major
problem is that both processes are very hard to accurately predict even with
state-of-the-art simulations; recombination emissivity could be significantly
boosted without changing total hydrogen column density if there are many tiny
($\ll 1$ pkpc) clumps with locally increased density in the CGM/IGM regions
unresolved in current standard simulations but suggested from observations
(Rauch et al., 1999; Cantalupo et al., 2014, 2019; McCourt et al., 2018;
Hummels et al., 2019; van de Voort et al., 2019). The total Ly$\alpha$
luminosity from gravitational cooling should have a strong dependence on halo
mass (Goerdt et al., 2010). In addition, the emissivity of the collisional
process depends extremely sensitively on temperature exponentially in the
range $T=10^{4}$–$10^{5}$ K characteristic of cold accretion, and treatments
of the effect of self-shielding against the UVB may have a critical impact on
results (Faucher-Giguère et al., 2010; Kollmeier et al., 2010; Rosdahl &
Blaizot, 2012).

There remains a possibility that our main conclusion about the dominance of
central/satellite star formation may apply only to relatively massive LAEs,
since lower-mass halos would have less scattering media and less satellites.
For example, high-EW LAEs are efficient producers of Ly$\alpha$ photons. But
because they are on average less massive and should have lower HI gas (Rakic
et al., 2013; Turner et al., 2017), they have more compact LAHs despite
efficient Ly$\alpha$ production. Kakiichi & Dijkstra (2018) showed that
scattering of Ly$\alpha$ photons produced by central galaxies with realistic
HI distributions constrained by Ly$\alpha$ forest observations of LBGs results
in power-law-like Ly$\alpha$ SB profiles (see also Steidel et al., 2011). In
Figure 11, UV/Ly$\alpha$-faint and high-EW LAEs seem to deviate from power-law
fits at $r>25$ kpkc. This could be a hint of the dominance of other processes.
Similar conclusion of the dominance of scattered Ly$\alpha$ from central
galaxies and possible contribution from the other processes is reached by
recent observational studies (Lujan Niemeyer et al., 2022a, b). In this way,
much information is buried beyond the “flattening radius” around 20 pkpc,
outside of which contributions from central/satellites appear insufficient to
explain the observations. With a larger sample and deeper data, we can further
probe the behavior of LAHs e.g., by making EW-based subsamples with matched UV
magnitude, etc.

##  6 Summary

We have investigated the rest-frame UV continuum ($\lambda\sim 1225$ Å) and
Ly$\alpha$ radial surface brightness (SB) profiles of LAEs at $z=2.84$ through
stacking analyses of UV and Ly$\alpha$ images created from Subaru/HSC g-band
and the NB468 narrow-band images. The depth and wide field coverage, including
a known protocluster, enable us to study both SB profiles with unprecedented
depth because of the large sample ($N=3490$) at $z\sim 3$. Our major findings
are as follows:

  1. 1.

Stacking of 3490 LAEs yields a SB sensitivity of $\sim 1\times
10^{-20}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}arcsec^{-2}}$ in Ly$\alpha$ and
$\sim 1\times
10^{-33}\mathrm{~{}erg~{}s^{-1}~{}cm^{-2}~{}Hz^{-1}~{}arcsec^{-2}}$ (Figure
9). Our analyses reveal that systematic errors should be at the same levels at
most (Section 4.1 and that the choice of mesh size for local sky estimation
could have a large impact on the results (Figures 7 and 8).

  2. 2.

By dividing the LAEs into subsamples according to various photometric
properties, UV magnitude, Ly$\alpha$ luminosity, Ly$\alpha$ EW, projected
distance from a hyperluminous (HL) QSO residing at the center of the
protocluster (as a proxy for radiation field strength boosted by the HLQSO),
and LAE overdensity $\delta$ on a 1$.\\!\\!\arcmin$8 ($\sim 840$ pkpc) scale,
we study the dependence of the SB profiles on these quantities (Figures 11 and
12). To quantify the radial dependence of SB profiles, we fit 2-component
exponential functions (Equation 1) to observed profiles. For Ly$\alpha$ SB
profiles, we consistently obtain $r_{\mathrm{2,Ly\alpha}}$, the scale-length
of the more extended component, of $\sim 10$ pkpc for all subsamples. However,
we do not observe any clear trend of $r_{\mathrm{2,Ly\alpha}}$ with any
property probed here (Figure 13), whereas the scale-length of the compact
components (both of $r_{\mathrm{1,Ly\alpha}}$ and $r_{\mathrm{1,UV}}$) varies
monotonically with respect to UV magnitude, Ly$\alpha$ luminosity, and
Ly$\alpha$ EW.

  3. 3.

We find an exceptionally large exponential scale-length
$r_{\mathrm{2,Ly\alpha}}$ for LAEs in the inner core (those within 500 pkpc
from the HLQSO) of the protocluster, and a significant variation in
$r_{\mathrm{2,Ly\alpha}}$ with respect to the projected distance from the
HLQSO. These findings could be explained by enhanced ionizing background
radiation due to abundant active sources and cool gas at the forming
protocluster core and the past activity of HLQSO, respectively (Section 5.2.1
and Section 5.2.2).

  4. 4.

We for the first time identify extended UV components (i.e.,
$r_{\mathrm{2,UV}}$ inconsistent with zero), or “UV halos” around some bright
LAE subsamples, which provides direct evidence for the contribution of star
formation in halo regions and/or satellite galaxies to Ly$\alpha$ halos.
Comparison with cosmological hydrodynamical simulations suggests that UV halos
could be composed of 1–2 low-mass ($M_{\mathrm{DM}}\sim
10^{9.5}~{}\mathrm{M_{\odot}}$) galaxies (Section 5.3) with total SFR of $\sim
10\%$ of that of their central galaxies.

  5. 5.

Combining our results with predictions of recent numerical simulations, we
conclude that star formation in both the central galaxy and in satellites,
together with resonant scattering of Ly$\alpha$ photons are the dominant
factors determining the Ly$\alpha$ SB profiles at least within a few tens
pkpc. In outer regions (projected distances $\gtrsim 30$ pkpc) other
mechanisms such as fluorescence can also play a role especially in certain
situations like dense regions of the Universe and near zones of bright QSOs
(Section 5.4).

The low-mass satellite galaxies suggested by our deep UV stacked images will
be very important targets for revealing the role of minor mergers in galaxy
evolution and cosmic reionization, as they are believed to be promising
analogues of the main galaxy contributors to ionizing photon budget at $z>6$
(Robertson et al., 2013; Mason et al., 2015). In a year or so, deep surveys
with JWST will detect these galaxies within several arcsecs from LAE-class
galaxies, which are to be observed in coming programs. Finally, H$\alpha$
observations of LAHs open up a new pathway to study star formation in the CGM
and fluorescent clumps without the blurring effect of resonant scattering.
Observations of H$\alpha$ from $z=2.84$ has not been possible with ground-
based telescopes due to heavy atmospheric absorption and extremely bright
backgrounds, but now this is also becoming feasible thanks to JWST. In the
future, we will combine cross-analyses with a larger LAE sample and new
constraints from JWST to further constrain the origin of LAHs.

We are grateful to the anonymous referee for their careful reading and
constructive comments and suggestions, Chris Byrohl for providing simulation
data, and Haibin Zhang, Masafumi Yagi, Masayuki Umemura, Tadafumi Takata,
Kazuhiro Shimasaku, Yusei Koyama, and Kazuhiro Hada for fruitful discussions.
We thank Yukie Oishi and the HSC pipeline team for their helpful comments on
HSC data analyses. We would like to acknowledge all who supported our
observations at the Subaru Telescope, including the staff of the National
Astronomical Observatory of Japan, Maunakea Observatories, and the local
Hawaiian people who have been making efforts to preserve and share the
beautiful dark sky of Maunakea with us. We are honored and grateful for the
opportunity of observing the Universe from Maunakea, which has the cultural,
historical, and natural significance in Hawaii. Data analysis was carried out
on the Multi-wavelength Data Analysis System operated by the Astronomy Data
Center (ADC), National Astronomical Observatory of Japan and on analysis
servers at Center for Computational Astrophysics, National Astronomical
Observatory of Japan. SK acknowledges supports from the Japan Society for the
Promotion of Science (JSPS) KAKENHI grant Nos. 18J11477, 19H00697 and the
Course-by-Course Education Program of SOKENDAI. YM acknowledges supports from
the JSPS KAKENHI grant Nos. 25287043, 17H04831, 17KK0098. CCS acknowledges
support by US NSF grant AST-2009278. ZZ acknowledges support by US NSF grant
AST-2007499.

##  Appendix A Comparison of Lyman-alpha and UV SB profiles for each subsample

In Figures 19 and 20, we show Ly$\alpha$ (left) and UV (right) SB profiles of
each subsample in each row for easier comparison.

In the top-left panel of Figure 20, the curves appear to deviate from each
other beyond 2 arcsec, although the sensitivity at these angular separations
is not very high. To check whether this difference is significant, we stacked
Ly$\alpha$ images of 700 randomly selected LAEs (roughly corresponding to the
number of LAEs in each bin of projected distance subsample) from the whole
($N=3490$) LAE sample and derived its Ly$\alpha$ SB profile. We repeated this
500 times and derived the 5th and 95th percentile of the SB distribution in
each annulus. These are shown as a gray shaded region in Figure 20; in the
bottom panel, we also show the SB distribution with 1000 randomly selected
LAEs to see the difference between $\delta<1.0$ subsamples. The curves are
almost within the shaded regions, suggesting that the difference is apparently
marginal.

Figure 19: Ly$\alpha$ SB profiles (Left) and UV continuum SB profiles (Right)
derived by stacking analysis. From top to bottom, the LAE sample is grouped by
their UV magnitude, Ly$\alpha$ luminosity, and Ly$\alpha$ equivalent width in
the manner specified in Table 1. The error bars are slightly shifted
horizontally for display purpose. Downward triangles show 1$\sigma$ error
levels after residual sky subtraction.

Figure 20: Same as Figure 19, but here the LAEs are grouped by their projected
distances from the HLQSO (Top) and local environment (Bottom) in the manner
specified in Table 1. The gray shaded regions in the Top (Bottom) panel shows
5th and 95th percentile of the Ly$\alpha$ SB distribution of stacked images
created with randomly selected 700 (1000) LAEs.

## References

  * Abolfathi et al. (2018) Abolfathi, B., Aguado, D. S., Aguilar, G., et al. 2018, The Astrophysical Journal Supplement Series, 235, 42, doi: 10.3847/1538-4365/aa9e8a
  * Adelberger et al. (2005) Adelberger, K. L., Shapley, A. E., Steidel, C. C., et al. 2005, The Astrophysical Journal, 629, 636, doi: 10.1086/431753
  * Alavi et al. (2016) Alavi, A., Siana, B., Richard, J., et al. 2016, The Astrophysical Journal, 832, 56, doi: 10.3847/0004-637x/832/1/56
  * Ando et al. (2006) Ando, M., Ohta, K., Iwata, I., et al. 2006, The Astrophysical Journal, 645, L9, doi: 10.1086/505652
  * Bacon et al. (2021) Bacon, R., Mary, D., Garel, T., et al. 2021, Astronomy and Astrophysics, 647, 1, doi: 10.1051/0004-6361/202039887
  * Becker & Bolton (2013) Becker, G. D., & Bolton, J. S. 2013, Monthly Notices of the Royal Astronomical Society, 436, 1023, doi: 10.1093/mnras/stt1610
  * Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, Astronomy and Astrophysics Supplement Series, 117, 393, doi: 10.1051/aas:1996164
  * Bond et al. (2010) Bond, N. A., Feldmeier, J. J., Matković, A., et al. 2010, Astrophysical Journal Letters, 716, L200, doi: 10.1088/2041-8205/716/2/L200
  * Borisova et al. (2016) Borisova, E., Lilly, S. J., Cantalupo, S., et al. 2016, Astrophys. J., 830, 120, doi: 10.3847/0004-637X/830/2/120
  * Bouwens et al. (2022) Bouwens, R. J., Illingworth, G., Ellis, R. S., et al. 2022, The Astrophysical Journal, 931, 81, doi: 10.3847/1538-4357/ac618c
  * Byrohl et al. (2021) Byrohl, C., Nelson, D., Behrens, C., et al. 2021, Monthly Notices of the Royal Astronomical Society, 506, 5129, doi: 10.1093/mnras/stab1958
  * 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, The Astrophysical Journal, 628, 61, doi: 10.1086/430758
  * Cantalupo et al. (2019) Cantalupo, S., Pezzulli, G., Lilly, S. J., et al. 2019, Monthly Notices of the Royal Astronomical Society, 483, 5188, doi: 10.1093/mnras/sty3481
  * Cen & Zheng (2013) Cen, R., & Zheng, Z. 2013, The Astrophysical Journal, 775, 112, doi: 10.1088/0004-637X/775/2/112
  * Chen et al. (2020) Chen, Y., Steidel, C. C., Hummels, C. B., et al. 2020, Monthly Notices of the Royal Astronomical Society, 499, 1721, doi: 10.1093/mnras/staa2808
  * Chen et al. (2021) Chen, Y., Steidel, C. C., Erb, D. K., et al. 2021, Monthly Notices of the Royal Astronomical Society, 508, 19, doi: 10.1093/mnras/stab2383
  * Dekel et al. (2009) Dekel, A., Sari, R., & Ceverino, D. 2009, The Astrophysical Journal, 703, 785, doi: 10.1088/0004-637X/703/1/785
  * Dijkstra & Kramer (2012) Dijkstra, M., & Kramer, R. 2012, Monthly Notices of the Royal Astronomical Society, 424, 1672, doi: 10.1111/j.1365-2966.2012.21131.x
  * Dijkstra & Loeb (2009) Dijkstra, M., & Loeb, A. 2009, Monthly Notices of the Royal Astronomical Society, 400, 1109, doi: 10.1111/j.1365-2966.2009.15533.x
  * D’Souza et al. (2014) D’Souza, R., Kauffman, G., Wang, J., & Vegetti, S. 2014, Monthly Notices of the Royal Astronomical Society, 443, 1433, doi: 10.1093/mnras/stu1194
  * Erb et al. (2011) Erb, D. K., Bogosavljević, M., & Steidel, C. C. 2011, Astrophysical Journal Letters, 740, L31, doi: 10.1088/2041-8205/740/1/L31
  * Erb et al. (2018) Erb, D. K., Steidel, C. C., & Chen, Y. 2018, The Astrophysical Journal, 862, L10, doi: 10.3847/2041-8213/aacff6
  * Faucher-Giguère et al. (2010) Faucher-Giguère, C.-A., Kereš, D., Dijkstra, M., Hernquist, L., & Zaldarriaga, M. 2010, The Astrophysical Journal, 725, 633, doi: 10.1088/0004-637X/725/1/633
  * Feldmeier et al. (2013) Feldmeier, J. J., Hagen, A., Ciardullo, R., et al. 2013, The Astrophysical Journal, 776, 75, doi: 10.1088/0004-637X/776/2/75
  * Goerdt et al. (2010) Goerdt, T., Dekel, A., Sternberg, A., et al. 2010, Monthly Notices of the Royal Astronomical Society, 407, 613, doi: 10.1111/j.1365-2966.2010.16941.x
  * Gronke & Bird (2017) Gronke, M., & Bird, S. 2017, The Astrophysical Journal, 835, 207, doi: 10.3847/1538-4357/835/2/207
  * Haardt & Madau (2012) Haardt, F., & Madau, P. 2012, The Astrophysical Journal, 746, 125, doi: 10.1088/0004-637X/746/2/125
  * Haiman & Rees (2001) Haiman, Z., & Rees, M. J. 2001, The Astrophysical Journal, 556, 87, doi: 10.1086/321567
  * Hashimoto et al. (2017a) Hashimoto, T., Ouchi, M., Shimasaku, K., et al. 2017a, Monthly Notices of the Royal Astronomical Society, 465, 1543, doi: 10.1093/mnras/stw2834
  * Hashimoto et al. (2017b) Hashimoto, T., Garel, T., Guiderdoni, B., et al. 2017b, Astronomy & Astrophysics, 608, A10, doi: 10.1051/0004-6361/201731579
  * Hathi et al. (2016) Hathi, N. P., Le Fèvre, O., Ilbert, O., et al. 2016, Astronomy & Astrophysics, 588, A26, doi: 10.1051/0004-6361/201526012
  * Hayashino et al. (2004) Hayashino, T., Matsuda, Y., Tamura, H., et al. 2004, The Astronomical Journal, 128, 2073, doi: 10.1086/424935
  * Hummels et al. (2019) Hummels, C. B., Smith, B. D., Hopkins, P. F., et al. 2019, The Astrophysical Journal, 882, 156, doi: 10.3847/1538-4357/ab378f
  * Infante-Sainz et al. (2019) Infante-Sainz, R., Trujillo, I., & Román, J. 2019, Monthly Notices of the Royal Astronomical Society, 14, 1, doi: 10.1093/mnras/stz3111
  * Kakiichi & Dijkstra (2018) Kakiichi, K., & Dijkstra, M. 2018, Monthly Notices of the Royal Astronomical Society, 480, 5140, doi: 10.1093/mnras/sty2214
  * Kennicutt & Evans (2012) Kennicutt, R. C., & Evans, N. J. 2012, Annual Review of Astronomy and Astrophysics, 50, 531, doi: 10.1146/annurev-astro-081811-125610
  * Kereš et al. (2009) Kereš, D., Katz, N., Fardal, M., Davé, R., & Weinberg, D. H. 2009, Monthly Notices of the Royal Astronomical Society, 395, 160, doi: 10.1111/j.1365-2966.2009.14541.x
  * Kereš et al. (2005) Kereš, D., Katz, N., Weinberg, D. H., & Davé, R. 2005, Monthly Notices of the Royal Astronomical Society, 363, 2, doi: 10.1111/j.1365-2966.2005.09451.x
  * Kida et al. (2019) Kida, A., Kikuta, S., & Matsuda, Y. 2019, Proceedings of the International Astronomical Union, 15, 277, doi: 10.1017/S1743921319002631
  * Kikuta et al. (2019) Kikuta, S., Matsuda, Y., Cen, R., et al. 2019, Publications of the Astronomical Society of Japan, 71, 1, doi: 10.1093/pasj/psz055
  * Kollmeier et al. (2010) Kollmeier, J. A., Zheng, Z., Davé, R., et al. 2010, The Astrophysical Journal, 708, 1048, doi: 10.1088/0004-637X/708/2/1048
  * Kusakabe et al. (2019) Kusakabe, H., Shimasaku, K., Momose, R., et al. 2019, Publications of the Astronomical Society of Japan, 71, 1, doi: 10.1093/pasj/psz029
  * Kusakabe et al. (2022) Kusakabe, H., Verhamme, A., Blaizot, J., et al. 2022, Astronomy & Astrophysics, 660, A44, doi: 10.1051/0004-6361/202142302
  * Lacaille et al. (2019) Lacaille, K. M., Chapman, S. C., Smail, I., et al. 2019, Monthly Notices of the Royal Astronomical Society, 488, 1790, doi: 10.1093/mnras/stz1742
  * Lake et al. (2015) Lake, E., Zheng, Z., Cen, R., et al. 2015, Astrophysical Journal, 806, 46, doi: 10.1088/0004-637X/806/1/46
  * Laursen et al. (2009) Laursen, P., Sommer-Larsen, J., & Andersen, A. C. 2009, Astrophysical Journal, 704, 1640, doi: 10.1088/0004-637X/704/2/1640
  * Leclercq et al. (2017) Leclercq, F., Bacon, R., Wisotzki, L., et al. 2017, Astronomy & Astrophysics, 608, A8, doi: 10.1051/0004-6361/201731480
  * Lujan Niemeyer et al. (2022a) Lujan Niemeyer, M., Komatsu, E., Byrohl, C., et al. 2022a, The Astrophysical Journal, 929, 90, doi: 10.3847/1538-4357/ac5cb8
  * Lujan Niemeyer et al. (2022b) Lujan Niemeyer, M., Bowman, W. P., Ciardullo, R., et al. 2022b, The Astrophysical Journal Letters, 934, L26, doi: 10.3847/2041-8213/ac82e5
  * Mas-Ribas & Dijkstra (2016) Mas-Ribas, L., & Dijkstra, M. 2016, The Astrophysical Journal, 822, 84, doi: 10.3847/0004-637X/822/2/84
  * Mason et al. (2015) Mason, C. A., Trenti, M., & Treu, T. 2015, The Astrophysical Journal, 813, 21, doi: 10.1088/0004-637X/813/1/21
  * Matsuda et al. (2012) Matsuda, Y., Yamada, T., Hayashino, T., et al. 2012, Monthly Notices of the Royal Astronomical Society, 425, 878, doi: 10.1111/j.1365-2966.2012.21143.x
  * Mawatari et al. (2012) Mawatari, K., Yamada, T., Nakamura, Y., Hayashino, T., & Matsuda, Y. 2012, The Astrophysical Journal, 759, 133, doi: 10.1088/0004-637X/759/2/133
  * McCourt et al. (2018) McCourt, M., Oh, S. P., O’Leary, R., & Madigan, A.-M. 2018, Monthly Notices of the Royal Astronomical Society, 473, 5407, doi: 10.1093/mnras/stx2687
  * Mitchell et al. (2021) Mitchell, P. D., Blaizot, J., Cadiou, C., et al. 2021, Monthly Notices of the Royal Astronomical Society, 5775, 5757, doi: 10.1093/mnras/stab035
  * Miyazaki et al. (2012) Miyazaki, S., Komiyama, Y., Nakaya, H., et al. 2012, Proc.SPIE, 8446, 84460Z, doi: 10.1117/12.926844
  * Momose et al. (2021) Momose, R., Shimasaku, K., Nagamine, K., et al. 2021, The Astrophysical Journal Letters, 912, L24, doi: 10.3847/2041-8213/abf04c
  * Momose et al. (2014) Momose, R., Ouchi, M., Nakajima, K., et al. 2014, Monthly Notices of the Royal Astronomical Society, 442, 110, doi: 10.1093/mnras/stu825
  * Momose et al. (2016) —. 2016, Monthly Notices of the Royal Astronomical Society, 457, 2318, doi: 10.1093/mnras/stw021
  * Mori et al. (2004) Mori, M., Umemura, M., & Ferrara, A. 2004, The Astrophysical Journal, 613, L97, doi: 10.1086/425255
  * Mostardi et al. (2013) Mostardi, R. E., Shapley, A. E., Nestor, D. B., et al. 2013, The Astrophysical Journal, 779, 65, doi: 10.1088/0004-637X/779/1/65
  * Murphy et al. (2011) Murphy, E. J., Condon, J. J., Schinnerer, E., et al. 2011, The Astrophysical Journal, 737, 67, doi: 10.1088/0004-637X/737/2/67
  * Muzahid et al. (2021) Muzahid, S., Schaye, J., Cantalupo, S., et al. 2021, Monthly Notices of the Royal Astronomical Society, 508, 5612, doi: 10.1093/mnras/stab2933
  * Nelson et al. (2019) Nelson, D., Springel, V., Pillepich, A., et al. 2019, Computational Astrophysics and Cosmology, 6, 2, doi: 10.1186/s40668-019-0028-x
  * Nilsson et al. (2009a) Nilsson, K. K., Möller-Nilsson, O., Møller, P., Fynbo, J. P. U., & Shapley, A. E. 2009a, Monthly Notices of the Royal Astronomical Society, 400, 232, doi: 10.1111/j.1365-2966.2009.15439.x
  * Nilsson et al. (2009b) Nilsson, K. K., Tapken, C., Møller, P., et al. 2009b, Astronomy & Astrophysics, 498, 13, doi: 10.1051/0004-6361/200810881
  * Pillepich et al. (2018) Pillepich, A., Nelson, D., Hernquist, L., et al. 2018, Monthly Notices of the Royal Astronomical Society, 475, 648, doi: 10.1093/mnras/stx3112
  * Rakic et al. (2013) Rakic, O., Schaye, J., Steidel, C. C., et al. 2013, Monthly Notices of the Royal Astronomical Society, 433, 3103, doi: 10.1093/mnras/stt950
  * Rauch et al. (1999) Rauch, M., Sargent, W. L. W., & Barlow, T. A. 1999, The Astrophysical Journal, 515, 500, doi: 10.1086/307060
  * Rauch et al. (2008) Rauch, M., Haehnelt, M., Bunker, A., et al. 2008, The Astrophysical Journal, 681, 856, doi: 10.1086/525846
  * Robertson et al. (2013) Robertson, B. E., Furlanetto, S. R., Schneider, E., et al. 2013, Astrophysical Journal, 768, doi: 10.1088/0004-637X/768/1/71
  * Rosdahl & Blaizot (2012) Rosdahl, J., & Blaizot, J. 2012, Monthly Notices of the Royal Astronomical Society, 423, 344, doi: 10.1111/j.1365-2966.2012.20883.x
  * Rudie et al. (2012) Rudie, G. C., Steidel, C. C., Trainor, R. F., et al. 2012, The Astrophysical Journal, 750, 67, doi: 10.1088/0004-637X/750/1/67
  * Scarlata et al. (2009) Scarlata, C., Colbert, J., Teplitz, H. I., et al. 2009, Astrophysical Journal, 704, 98, doi: 10.1088/0004-637X/704/2/L98
  * Schaerer (2003) Schaerer, D. 2003, Astronomy and Astrophysics, 397, 527, doi: 10.1051/0004-6361:20021525
  * Shapiro et al. (2006) Shapiro, P. R., Iliev, I. T., Alvarez, M. A., & Scannapieco, E. 2006, The Astrophysical Journal, 648, 922, doi: 10.1086/506242
  * Shibuya et al. (2019) Shibuya, T., Ouchi, M., Harikane, Y., & Nakajima, K. 2019, The Astrophysical Journal, 871, 164, doi: 10.3847/1538-4357/aaf64b
  * Shibuya et al. (2014) Shibuya, T., Ouchi, M., Nakajima, K., et al. 2014, Astrophysical Journal, 785, doi: 10.1088/0004-637X/785/1/64
  * Shimakawa et al. (2017) Shimakawa, R., Kodama, T., Hayashi, M., et al. 2017, Monthly Notices of the Royal Astronomical Society: Letters, 468, L21, doi: 10.1093/mnrasl/slx019
  * Steidel et al. (2011) Steidel, C. C., Bogosavljević, M., Shapley, A. E., et al. 2011, The Astrophysical Journal, 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, Astrophysical Journal, 717, 289, doi: 10.1088/0004-637X/717/1/289
  * Takada et al. (2014) Takada, M., Ellis, R. S., Chiba, M., et al. 2014, Extragalactic science, cosmology, and Galactic archaeology with the Subaru Prime Focus Spectrograph, doi: 10.1093/pasj/pst019
  * Taniguchi & Shioya (2000) Taniguchi, Y., & Shioya, Y. 2000, The Astrophysical Journal, 532, L13, doi: 10.1086/312557
  * Trainor & Steidel (2013) Trainor, R., & Steidel, C. C. 2013, The Astrophysical Journal, 775, L3, doi: 10.1088/2041-8205/775/1/L3
  * Trainor & Steidel (2012) Trainor, R. F., & Steidel, C. C. 2012, The Astrophysical Journal, 752, 39, doi: 10.1088/0004-637X/752/1/39
  * Tumlinson et al. (2017) Tumlinson, J., Peeples, M. S., & Werk, J. K. 2017, Annual Review of Astronomy and Astrophysics, 55, 389, doi: 10.1146/annurev-astro-091916-055240
  * Turner et al. (2017) Turner, M. L., Schaye, J., Crain, R. A., et al. 2017, Monthly Notices of the Royal Astronomical Society, 471, 690, doi: 10.1093/MNRAS/STX1616
  * Umehata et al. (2019) Umehata, H., Fumagalli, M., Smail, I., et al. 2019, Science, 366, 97, doi: 10.1126/science.aaw5949
  * van de Voort et al. (2019) van de Voort, F., Springel, V., Mandelker, N., van den Bosch, F. C., & Pakmor, R. 2019, Monthly Notices of the Royal Astronomical Society: Letters, 482, L85, doi: 10.1093/mnrasl/sly190
  * Wisotzki et al. (2016) Wisotzki, L., Bacon, R., Blaizot, J., et al. 2016, Astronomy & Astrophysics, 587, A98, doi: 10.1051/0004-6361/201527384
  * Wisotzki et al. (2018) Wisotzki, L., Bacon, R., Brinchmann, J., et al. 2018, Nature, doi: 10.1038/s41586-018-0564-6
  * Xue et al. (2017) Xue, R., Lee, K.-S., Dey, A., et al. 2017, The Astrophysical Journal, 837, 172, doi: 10.3847/1538-4357/837/2/172
  * Yajima et al. (2013) Yajima, H., Li, Y., & Zhu, Q. 2013, Astrophysical Journal, 773, doi: 10.1088/0004-637X/773/2/151
  * Zheng et al. (2010) Zheng, Z., Cen, R., Trac, H., & Miralda-Escudé, J. 2010, Astrophysical Journal, 716, 574, doi: 10.1088/0004-637X/716/1/574
  * Zheng et al. (2011) Zheng, Z., Cen, R., Weinberg, D., Trac, H., & Miralda-Escudé, J. 2011, The Astrophysical Journal, 739, 62, doi: 10.1088/0004-637X/739/2/62