UV & Ly$\alpha$ halos of Ly$\alpha$ emitters across environments at z=2.84

TL;DR

This study analyzes UV and Lyα halos around Lyα emitters at z=2.84, revealing environmental influences, the presence of UV halos likely from satellite galaxies, and the role of star formation and fluorescence in shaping Lyα halos.

Contribution

First to identify UV halos around LAEs and analyze their dependence on environment and satellite contributions at high redshift.

Findings

LAHs have scale-lengths of ~1.56 and 10.4 pkpc.

LAEs in protocluster centers have larger LAHs.

UV halos likely originate from satellite galaxies.

Abstract

We present UV and 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 $\sim1\times10^{-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\sim3$. Fitting of the two-component exponential function gives the scale-lengths of $1.56\pm0.01$ and $10.4\pm0.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.

Tables

Table 1. Table 1: Subsample used in our stacking analysis and the result of fitting performed in Section 4.2 .

quantity	criteria	median	$N$	$C_{1,\mathrm{Ly}\alpha }$	$r_{1,\mathrm{Ly}\alpha }$	$C_{2,\mathrm{Ly}\alpha }$	$r_{2,\mathrm{Ly}\alpha }$	$C_{3,\mathrm{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$
		-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$
		-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.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$
		-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$
		-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.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.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.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.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$
		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$
	Å	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$
	Å	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	Å	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$
	Å	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 Å	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$
	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$
	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	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	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$
	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$
		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.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.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.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$
		-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$

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

criteria	$C_{1,\mathrm{UV}}$	$r_{1,\mathrm{UV}}$
(1)	(2)	(3)
All	$8.146\pm 0.072$	$1.486\pm 0.019$
	$26.830\pm 0.144$	$1.632\pm 0.013$
	$13.896\pm 0.178$	$1.285\pm 0.027$
	$8.049\pm 0.179$	$1.203\pm 0.046$
	$5.552\pm 0.309$	$0.495\pm 0.214$
	–	–
	$25.078\pm 0.157$	$1.538\pm 0.014$
	$13.285\pm 0.164$	$1.278\pm 0.026$
	$8.606\pm 0.210$	$1.147\pm 0.050$
	$6.118\pm 0.220$	$1.042\pm 0.073$
	–	–
Å	$21.208\pm 0.142$	$1.710\pm 0.016$
Å	$11.772\pm 0.141$	$1.614\pm 0.028$
Å	$8.182\pm 0.170$	$1.353\pm 0.045$
Å	$6.173\pm 0.264$	$0.690\pm 0.102$
160 Å	–	–
pMpc	$8.451\pm 0.168$	$1.395\pm 0.043$
pMpc	$8.716\pm 0.154$	$1.458\pm 0.039$
pMpc	$7.798\pm 0.163$	$1.495\pm 0.047$
pMpc	$8.229\pm 0.149$	$1.475\pm 0.040$
pMpc	$7.512\pm 0.148$	$1.625\pm 0.046$
	$9.423\pm 0.471$	$1.815\pm 0.123$
	$8.319\pm 0.212$	$1.391\pm 0.054$
	$8.050\pm 0.143$	$1.445\pm 0.038$
	$7.828\pm 0.129$	$1.507\pm 0.036$
	$8.431\pm 0.133$	$1.530\pm 0.034$

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

criteria	$C_{1,\mathrm{UV}}$, $C_{2,\mathrm{UV}}$	$r_{1,\mathrm{UV}}$, $r_{2,\mathrm{UV}}$
(1)	(2)	(3)
	$26.730\pm 0.389$	$1.280\pm 0.048$
	$1.969\pm 0.495$	$4.845\pm 0.532$
	$22.466\pm 1.150$	$0.874\pm 0.110$
	$5.483\pm 0.495$	$3.072\pm 0.289$
Å	$20.804\pm 0.531$	$1.343\pm 0.071$
	$1.818\pm 0.660$	$4.570\pm 0.679$
Å	$12.102\pm 0.169$	$1.407\pm 0.050$
	$0.250\pm 0.080$	$10.586\pm 2.201$

Equations

$$\mathrm{PSF}\ast \left[C_1\exp \left(-\frac{r}{r_1}\right)+C_2\exp \left(-\frac{r}{r_2}\right)\right]$$

$$\mathrm{PSF}\ast C_3{r}^{-\alpha },$$

Full text

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.

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.

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$.

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).

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).

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.

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:

[TABLE]

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.

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.

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).

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.

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.

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.

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.

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:

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.