<figure><img src="image_1.png" /></figure>

# Crucial Factors for Ly$\alpha$ Transmission in the Reionizing Intergalactic Medium: Infall Motion, H II Bubble Size, and Self-shielded Systems

Hyunbae Park (박현 배)<sup>1,2</sup>①, Intae Jung<sup>3,4,5</sup>①, Hyunmi Song<sup>6</sup>, Pierre Ocvirk<sup>7</sup>, Paul R. Shapiro<sup>8</sup>①, Taha Dawoodbhoy<sup>8</sup>, Ilian T. Iliev<sup>9</sup>①, Kyungjin Ahn<sup>10</sup>①, Michele Bianco<sup>9</sup>, and Hyo Jeong Kim<sup>10</sup>

Received 2021 May 21; revised 2021 September 21; accepted 2021 October 9; published 2021 December 6

## Abstract

Using the CoDa II simulation, we study the Ly$\alpha$ transmissivity of the intergalactic medium (IGM) during reionization. At $z>6$, a typical galaxy without an active galactic nucleus fails to form a proximity zone around itself due to the overdensity of the surrounding IGM. The gravitational infall motion in the IGM makes the resonance absorption extend to the red side of Ly$\alpha$, suppressing the transmission up to roughly the circular velocity of the galaxy. In some sight lines, an optically thin blob generated by a supernova in a neighboring galaxy results in a peak feature, which can be mistaken for a blue peak. Redward of the resonance absorption, the damping-wing opacity correlates with the global IGM neutral fraction and the UV magnitude of the source galaxy. Brighter galaxies tend to suffer lower opacity because they tend to reside in larger H II regions, and the surrounding IGM transmits redder photons, which are less susceptible to attenuation, owing to stronger infall velocity. The H II regions are highly nonspherical, causing both sight-line-to-sight-line and galaxy-to-galaxy variation in opacity. Also, self-shielded systems within H II regions strongly attenuate the emission for certain sight lines. All these factors add to the transmissivity variation, requiring a large sample size to constrain the average transmission. The variation is largest for fainter galaxies at higher redshift. The 68% range of the transmissivity is similar to or greater than the median for galaxies with $M_{UV} \ge -21$ at $z \ge 7$, implying that more than a hundred galaxies would be needed to measure the transmission to 10% accuracy.

*Unified Astronomy Thesaurus concepts: Reionization (1383); Lyman-alpha galaxies (978); Hydrodynamical simulations (767); Radiative transfer simulations (1967); Cosmology (343)*

## 1. Introduction

How was the intergalactic medium (IGM) reionized at $z \gtrsim 6$? This question is connected to the unknown properties of early galaxies and is, therefore, one of the key questions of modern astronomy. Cosmic microwave background radiation observed by the Planck satellite constrains the midpoint of reionization to be at $z \sim 8$ (Planck Collaboration et al. 2020).¹¹Another powerful probe of cosmic history is the Ly$\alpha$ emission from high-$z$ star-forming galaxies. Those galaxies are expected to generate strong emission near Ly$\alpha$, which the neutral gas in the IGM can easily scatter. Thus, the statistics of the Ly$\alpha$ line strength is closely related to the global neutral fraction of the IGM (e.g., Malhotra & Rhoads 2004; Dijkstra 2014).

Today, numerous Ly$\alpha$ emitters (LAEs) have been observed up to $z \sim 9$ by narrowband imaging (e.g., Kashikawa et al. 2006;

Ouchi et al. 2008, 2010; Hu et al. 2010; Kashikawa et al. 2011; Konno et al. 2014; Matthee et al. 2015; Santos et al. 2016; Zheng et al. 2017; Ota et al. 2017; Inoue et al. 2018; Ouchi et al. 2018) and spectroscopic surveys (e.g., Deharveng et al. 2008; Adams et al. 2011; Finkelstein et al. 2013; Schenker et al. 2014; Cassata et al. 2015; Oesch et al. 2015; Zitrin et al. 2015; Song et al. 2016; Herenz et al. 2017; Stark et al. 2017; Hoag et al. 2017; Matthee et al. 2017; Jung et al. 2018; Shibuya et al. 2018; Jung et al. 2020; Khusanova et al. 2020; Morales et al. 2021; Taylor et al. 2021). The observed samples indeed show a steep decline in the number of LAEs found per galaxy above $z=6$, contrary to the much more gradual evolution at lower redshifts. This decline appears more prominent for relatively fainter galaxies ($M_{UV} \ge -20$) (e.g., Fontana et al. 2010; Stark et al. 2011; Pentericci et al. 2011; Ono et al. 2012; Maller et al. 2012; Curtis-Lake et al. 2012; Treu et al. 2013; Tilvi et al. 2014; Zheng et al. 2017; Endsley et al. 2021, Jung et al. 2021 in preparation). On the other hand, the number of bright LAEs tend to change less with redshift and are often found clustered in overdense regions (Castellano et al. 2018; Jung et al. 2019; Harikane et al. 2019; Higuchi et al. 2019; Jung et al. 2020; Tilvi et al. 2020; Hu et al. 2021). In most recent surveys, several LAEs have been detected at $z \sim 7$, constraining the neutral fraction at the highest redshifts ever (Hoag et al. 2019; Mason et al. 2019; Jung et al.

<figure><img src="image_2.png" /></figure>

<footnote n="1"><sup>1</sup>Lawrence Berkeley National Laboratory, CA 94720-8139, USA; hyunbae.park@ipmu.jp</footnote>

<footnote n="2"><sup>2</sup>Kavli IPMU (WPI), UTIAS, The University of Tokyo, Kashiwa, Chiba 277-8583, Japan</footnote>

<footnote n="3"><sup>3</sup>Astrophysics Science Division, Goddard Space Flight Center, Greenbelt, MD 20771, USA</footnote>

<footnote n="4"><sup>4</sup>Department of Physics, The Catholic University of America, Washington, DC 20064, USA</footnote>

<footnote n="5"><sup>5</sup>Center for Research and Exploration in Space Science and Technology, NASA/GSFC, Greenbelt, MD 20771, USA</footnote>

<footnote n="6"><sup>6</sup>Department of Astronomy and Space Science, Chungnam National University, Daejeon 34134, Republic of Korea</footnote>

<footnote n="7"><sup>7</sup>Universite de Strasbourg, CNRS, Observatoire astronomique de Strasbourg, UMR 7550, F-67000 Strasbourg, France</footnote>

<footnote n="8"><sup>8</sup>Department of Astronomy, University Texas, Austin, TX 78712-1083, USA</footnote>

<footnote n="9"><sup>9</sup>Astronomy Centre, Department of Physics & Astronomy, Pevensey III Building, University of Sussex, Falmer, Brighton, BN1 9QH, UK</footnote>

<footnote n="10"><sup>10</sup>Department of Earth Sciences, Chosun University, Gwangju 61452, Republic of Korea; kjahn@chosun.ac.kr</footnote>

<footnote n="11">¹¹We note that this statement is somewhat model dependent as the electron-scattering optical depth measured by Planck only constrains the integrated column density of free electrons. For a more detailed analysis on this matter, we refer the readers to, e.g., Heinrich & Hu (2021).</footnote>

<footnote>Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.</footnote>
2020; Morales et al. 2021) although the results still have a large uncertainty and have not converged yet.

In order to constrain the reionization history from observed Lyα emission statistics precisely, one crucially needs to understand the IGM opacity and how it evolves during reionization. Transmission of Lyα photons in the IGM depends on the neutral hydrogen density, velocity, and temperature of the intergalactic gas. Ultraviolet (UV) radiation from star-forming galaxies forms expanding H II regions around them, where some Lyα emission can be transmitted owing to low H I density. Thus, capturing the statistics of these H II regions, which can be up to tens of comoving megaparsecs (cMpc) large, is key to modeling the IGM transmissivity during the reionization era (e.g., Furlanetto & Oh 2005). Additionally, small-scale velocity and density structures below the megaparsec scale around the source galaxy are also important factors for the IGM transmissivity (Bolton & Haehnelt 2013; Mesinger et al. 2015; Kakiichi et al. 2016).

The IGM opacity of Lyα photons can be classified into two components: resonance and damping-wing opacity. At the line resonance, the cross section is so high that even a tiny trace of neutral fraction in highly ionized gas (well below 1%) can make the IGM opaque to Lyα photons. Thus, photons on the blue side of Lyα typically face high opacity at some point in the rest frame of the IGM and are scattered off the line of sight. On the other hand, photons on the red side of the resonance are much more likely to be transmitted. The IGM opacity for these photons mostly comes from the damping-wing cross section, which extends far away from the resonance. This damping-wing cross section is much smaller than the resonance cross section, but the H I density in neutral regions can be high enough to create a substantial optical depth. Therefore, this damping-wing opacity is expected to depend on the global neutral fraction of the IGM during reionization.

Our goal in this study is to clarify what aspects of the IGM density/velocity field modulate the resonance and damping-wing absorption and quantify how they depend on the progress of reionization. To this end, we shall analyze the Cosmic Dawn II (CoDa II) simulation of Ocvirk et al. (2020), which covers a suitable dynamic range to compute the IGM opacity during reionization and the relevant physics. This work parallels Gronke et al. (2021) in that both works are based on similar IGM transmission calculations from the same data set. Whereas Gronke et al. (2021) mainly focused on the transmission of the resonance to understand the recently observed blue peaks from some LAEs at $z \gtrsim 6.5$ (Songaila et al. 2018; Matthee et al. 2018; Meyer et al. 2021), we focus on the red side that is relevant to the evolution of the LAE number density during reionization.

The rest of the paper is as follows. In Section 2, we describe the CoDa II simulation and how we calculate the IGM transmission. In Section 3, we present our main results. In Section 4, we discuss our results in light of the current understanding of LAEs. For the numerical calculations, we adopt the cosmology parameters used by the CoDa II simulation for self-consistency. The parameters are $\Omega_{m,0} = 0.307$, $\Omega_{b,0} = 0.048$, and $h = 0.678$ (Planck Collaboration et al. 2014).

## 2. Method

### 2.1. Cosmic Dawn II Simulation

We compute the IGM transmission of Lyα photons during the Epoch of Reionization based on the simulated galaxies

and IGM of the CoDa II simulation. CoDa II is a fully coupled radiation-hydrodynamic simulation of the reionization era in a cubic volume of $[64 h^{-1} \text{cMpc}]^3$ on a 4096³ regular grid, which resolves atomically cooling halos down to $10^8 M_{\odot}$ with 200 dark matter particles. The simulation was initialized at $z=150$ and run down to $z=5.8$ using the RAMSES-CUDATON code (Ocvirk et al. 2016). The global ionization fraction of the IGM reaches 99% at $z \approx 6.2$ in the simulation. The simulated galaxies reproduced the observed UV luminosity function at $z \gtrsim 6$. The details of this simulation are elaborated on in Ocvirk et al. (2020). Here, we briefly list the features that are relevant to our work.

In CoDa II, photoionization is computed self-consistently by casting rays from stellar particles, which form from gas cells that satisfy certain star formation conditions. The CoDa II simulation captures a statistically meaningful number of H II regions with its large volume and resolves small-scale structures around galaxies needed to compute the transmission accurately. Figure 1 illustrates this with the H I density in the IGM around the UV-brightest galaxies at $z=6, 7$, and 8 from the full-box scale in the left column down to the virial radius scale in the right column. We also note that CoDa II uses the exact speed of light in its radiative transfer calculation and, thus, is free from possible problems arising from using the reduced speed-of-light approximation (Deparis et al. 2019; Ocvirk et al. 2019).

The UV luminosity of each galaxy is computed from the age of the stellar particles based on stellar synthesis models. We name galaxies according to their ranking in UV luminosity at 1600 Å (i.e., $M_{UV} \equiv M_{AB1600}$) within the snapshot to which it belongs. For example, we refer to the 100th UV-brightest galaxy as galaxy #0100. The brightest sample galaxy (galaxy #0001) of the $z=7$ snapshot is $M_{UV} = -23.1$ in UV magnitude, for example, and there are hundreds of galaxies up to $M_{UV} \sim -19$, which makes CoDa II suitable for understanding LAEs observed at $z \gtrsim 7$. Note that galaxy #n at $z=7$ and galaxy #n at $z=8$ are generally not the same because the ranking in UV magnitude changes over time.

### 2.2. IGM Transmission Calculation

We calculate the IGM transmissivity for photons emitted around the Lyα frequency from galaxies in the $z=6, 7$, and 8 snapshots. In the simulation, the IGM is almost fully ionized at $z=6$, while it is 50% and 13.2% ionized at $z=7$ and 8, respectively. Throughout this work, we assume neutral hydrogen is the only opacity source given that the dust and metals from star formation are unlikely to exist in meaningful amounts in the IGM during reionization.

The transmissivity, $T \equiv e^{-\tau}$, for a frequency $\nu_e$ at the rest frame of a source is given by the optical depth

$$ \tau(\nu_e) = \int_{s_{\min}}^{s_{\max}} n_{\text{HI}}(\mathbf{r}) \sigma_{\alpha}(\nu_a(\mathbf{r}), T(\mathbf{r})) a(z) ds, \quad (1) $$

where $s$ is the comoving distance between the source and the IGM, $a(z)$ is the scale factor at redshift $z$, $\mathbf{r} = \mathbf{r}_s + s\hat{\gamma}$ is the comoving coordinate of the IGM along the line-of-sight direction $\hat{\gamma}$ from the source location $\mathbf{r}_s$, $n_{\text{HI}}$ is the neutral hydrogen number density, $\sigma_{\alpha}(\nu, T)$ is the ensemble-averaged Lyα cross section of gas with temperature $T$, and

$$ \nu_a(\mathbf{r}) = \nu_e \left[ 1 - \frac{\nu_{\text{pe}}(\mathbf{r}) + sH(z)a(z)}{c} \right] $$
<figure><img src="image_3.png" /><figcaption>**Figure 1.** xy plane of the neutral hydrogen density map centered on the most massive galaxy in the CoDa II simulation at z = 8 (top), 7 (middle), and 6 (lower panels). We zoom in from the entire simulation box of (64 h⁻¹ cMpc)² to the (8 h⁻¹ cMpc)² and (1 h⁻¹ cMpc)² subboxes from left to right. The map is color-coded according to the neutral hydrogen density, but the blue/red color coincides with neutral/ionized regions in most cases. In the middle and right column panels, the partially transparent arrows describe the peculiar motion of gas with respect to the galaxy. The white circles in the left panels describe the maximum integration distance in Equation (1), smax = 24 h⁻¹ Mpc, for the local optical depth τL described in Section 2.2. The yellow circles in the right panels show the virial radius of the galaxies at the center, which is the minimum integration distance, smin = Rvir.</figcaption></figure>

describes the frequency shift in the rest frame of the IGM due to the IGM peculiar motion vpe and the cosmic expansion rate H(z). Here, vpe is calculated with respect to the peculiar motion of the galaxy, which we obtain by taking the density-weighted average of the velocity field with the virial radius of the host halo. We define the virial radius as the distance within which the mean density becomes 200 times the cosmic mean (Rvir ≡ r200).

The Lyα cross section is given by

$$ \sigma_{\alpha}(\nu, T) = 5.889 \times 10^{-14} T_{4}^{-0.5} \phi(x) \text{ cm}^{2}, \quad (2) $$

where T4 ≡ T/[10⁴ K] and φ(x) is the Voigt function as a function of the dimensionless frequency x ≡ [ν − να]/ΔνD.

Here, να ≈ 2.46 × 10¹⁵ Hz is the resonance frequency of Lyα, and

$$ \Delta\nu_D = \nu_\alpha \sqrt{\frac{2k_B T}{m_p c^2}} = 4.28 \times 10^{-5} T_4^{0.5} \nu_\alpha $$

is the thermal broadening frequency (Doppler width), kB is the Boltzmann constant, and mP is the proton mass. The Voigt function is given by

$$ \phi(x) = \frac{a_V}{\pi} \int_{-\infty}^{\infty} dy \frac{e^{-y^2}}{[y - x]^2 + a_V^2}, \quad (3) $$
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Lyα scattering cross section as a function of wavelength offset in velocity units ($v_α$). The blue, black, and red lines describe the cross section for three different gas temperatures, $T = 10^2$, $10^4$, and $10^6$ K, respectively.</figcaption></figure>

where $a_V \approx 4.7 \times 10^{-4} T_4^{-0.5}$ is the Voigt parameter. For fast computation, it is convenient to use the analytic fitting formula provided by Tasitsiomi (2006).

The shape of $\sigma_\alpha(\nu, T)$ is shown in Figure 2 for $T = 10^2$, $10^4$, and $10^6$ K. The profile of $\sigma_\alpha(\nu)$ has a core at the center with the damping wing extending outward. The FWHM is roughly 2, 20, and 200 km s$^{-1}$ for $T = 10^2$, $10^4$, and $10^6$ K, respectively. The gas temperature is approximately $T = 10^4$ K for the ionized IGM and $T \lesssim 10^2$ K for the neutral IGM. Some gas is shock-heated to $\sim 10^6$ K by supernova feedback, but such gas rarely exists outside the virial radius, where we shall calculate the IGM transmissivity.

For each galaxy, we calculate $\tau(\nu_e)$ for 2000 random lines of sight to account for the sight-line variation. For each sight line, we evaluate Equation (1) in two segments to obtain $\tau = \tau_L + \tau_{DW}$. Here, $\tau_L$ is the local contribution within a distance of $24 h^{-1}$ cMpc, which is far enough from the source galaxy to include the surrounding H II bubble. $\tau_{DW}$ is the damping-wing absorption by distant IGM beyond that distance.

For the local contribution $\tau_L$, we capture the local IGM dynamics and the inhomogeneity of reionization by sampling $n_{HI}(r)$, $T(r)$, and $v_{pe}(r)$ along each sight line by interpolating the mesh data of the simulation. In this work, we regard the IGM as the gas outside the virial radius ($R_{vir}$) and integrate to $24 h^{-1}$ cMpc from the galaxy by setting $s_{min} = R_{vir}$ and $s_{max} = 24 h^{-1}$ cMpc in Equation (1). We describe this integration range by white circles in Figure 1. Within this range, we do not consider the time evolution of IGM and calculate the opacity from the still snapshot of the target redshift. The circles in the left panels of Figure 1 describe the integration range for $\tau_L$, which is large enough to capture the surroundings of the ionized bubble at $z=7$ and 8 created by the central galaxy.

For the large-scale opacity contribution from the neutral gas in the IGM, $\tau_{DW}$, we integrate in Equation (1) from $s_{min} = 24 h^{-1}$ cMpc to $z=6$, where the reionization has ended in the simulation. For this large-scale contribution, we only consider the damping-wing opacity from neutral IGM and ignore the peculiar motion. For the neutral hydrogen density,

we use the globally averaged value given by

$$ \begin{aligned} \bar{n}_{\text{HI}} &= \bar{x}_{\text{HI}}(z)[1+z]^3 \Omega_b \frac{3H_0^2}{8\pi G} \frac{X}{m_p} \\ &= 9.62 \times 10^{-5} \bar{x}_{\text{HI}}(z) \left[ \frac{\Omega_b}{0.048} \right] \left[ \frac{h}{0.678} \right]^2 \left[ \frac{z+1}{8} \right]^3 \text{cm}^{-3}, \end{aligned} \quad (4) $$

where $\bar{x}_{\text{HI}}$ is the global average neutral fraction and $X=0.76$ is the mass fraction of hydrogen in baryons. In the damping wing, the Voigt profile can be approximated as $\phi(x) = a_V \pi^{-0.5} x^{-2}$, giving

$$ \tau_{\text{DW}}(v_\alpha) = 5.7 \left[ \frac{h}{0.678} \right] \left[ \frac{\Omega_b}{0.048} \right] \left[ \frac{\Omega_m}{0.307} \right]^{-1} \int_{s_{\text{min}}}^{s_{\text{max}}} \times \left[ \frac{ds}{h^{-1} \text{cMpc}} \right] \left[ \frac{z+1}{8} \right] \bar{x}_{\text{HI}}(z) \left[ \frac{s + \Delta s(v_\alpha)}{h^{-1} \text{cMpc}} \right]^{-2}, \quad (5) $$

where

$$ \begin{aligned} \Delta s(v_\alpha) = 0.64 [h^{-1} \text{cMpc}] \\ \times \left[ \frac{v_\alpha}{100 \text{ km s}^{-1}} \right] \left[ \frac{\Omega_m}{0.307} \right]^{-0.5} \left[ \frac{1+z}{8} \right]^{-0.5} \end{aligned} \quad (6) $$

accounts for the initial frequency offset of emitted photons from Lyα in the velocity unit: $v_α \equiv -c\Delta\nu/v_α$. Here, we have assumed $H \approx H_0 \sqrt{\Omega_m(1+z)^3}$, given that $z \gg 1$ during reionization. We calculate $\bar{x}_{\text{HI}}$ by interpolating in redshift the globally averaged values within each snapshot. In the case of $\tau_{\text{DW}}$, the integration starts from $s_{\text{min}} = 24 h^{-1}$ cMpc and ends at $s_{\text{max}} = 263 (482) h^{-1}$ cMpc for galaxies at $z=7$ (8). We find that roughly 10% (20%) of photons are scattered due to this large-scale damping-wing opacity at $z=7$ (8).

## 3. Results

Figure 3 shows the IGM transmissivity $T = e^{-\tau}$ as a function of frequency in the velocity unit ($v_α$) for galaxies #0001, #0014, and #0250 in the $z=6, 7$, and 8 snapshots. Note that the positive side of the x-axis ($v_α > 0$) is the red side of Lyα. The transmissions for different sight lines are shown as multiple partially transparent curves to depict the sight-line variation. The blue solid line shows the median value at each frequency, and the cyan lines bracket the 68% range around the median.

### 3.1. Resonance Absorption

The cross section of Lyα is large enough to make even the mostly ionized IGM opaque at the resonance. Photons that escape the circumgalactic medium (CGM) at a frequency bluer than the resonance are redshifted to the line center in the rest frame of the IGM at a certain point and encounter this high opacity of the resonance absorption. Notably, all the transmission curves in Figure 3 show a rapid drop at $v_α = 100-300$ km s$^{-1}$ on the red side of the Lyα, implying that the photons are blueshifted to the resonance in the rest frame of IGM.

This red-side resonance absorption occurs for two reasons. First, the peculiar velocity field of the IGM is dominated by the gravitational infall motion toward the source galaxy. Second,
<figure><img src="image_5.png" /><figcaption>**Figure 3.** IGM transmission as a function of frequency. The upper, middle, and lower panels show the results for galaxy #0001, #0014, and #0250, respectively, in the $z=8$ (left), 7 (center), and 6 (right) snapshots. The dashed lines denote the large-scale absorption ($\tau_{DW}$) coming from the IGM farther than $s > 24 h^{-1}$ cMpc from the source galaxy. Fifty out of 2000 individual sight lines are shown as partially transparent black solid lines. The blue solid line marks the median of the results from 2000 sight lines binned for each $v_\alpha$, and the cyan lines bracket the 68% range.</figcaption></figure>

the rising IGM density toward the galaxy keeps the IGM optically thick at the resonance frequency despite the rising ionizing intensity from the source galaxy. Without the infall motion and the overdensity of the IGM, the transmission cutoff should appear at or on the blue side of $Ly\alpha$ (e.g., Mason & Gronke 2020).

The velocity map around galaxies in Figure 1 clearly shows the infall motion from all directions. The infall motion extends to at least several comoving Mpc from the galaxy (see middle column panels) and can extend up to tens of Mpc (see Figure 6 of Iliev et al. 2008). There also are outflow motions generated by star formation feedback in hot low-density blobs shown as dark red spots. However, the outflowing gas is mostly confined within the virial radius, which we take as the CGM in this work.

In order to describe how the infall motion affects the IGM transmission, we plot gas density, H I density, peculiar velocity, and optical depth toward the $-y$ direction from galaxy #0001 at $z=7$ in Figure 4. The plotted quantities are baryon density, ionizing photon density, H I number density, peculiar velocity, the optical depth of IGM to photons emitted at three different frequencies at the galaxy center, and the local density map around the line of sight.

The peculiar velocity shown as the black solid line in panel d is negative (i.e., infalling) up to $\sim 5 h^{-1}$ cMpc from the source. The infall velocity becomes as large as $\approx 300$ km s$^{-1}$ when it peaks near $s = R_{vir}$ and gradually decreases toward larger distances. The infall velocity at a distance $s$ is similar to the circular velocity at that distance:

$$v_{\text{pe},m}(s) \equiv \sqrt{\frac{GM_h}{s}}, \qquad (7)$$

where $M_h$ is the total galaxy mass (see the green dotted line). This indicates that the gravitational force from the galaxy is the dominant cause of the infall motion. The peculiar velocity profile is also similar to the analytical models of Santos (2004) and Dijkstra et al. (2007).

The blueshift in the rest frame of the IGM due to the infall motion results in a high resonance opacity for some photons on the red side of $Ly\alpha$. The dashed lines in panel d show that the frequency change of the photons emitted at $v_{\alpha,i} = -200, 0$, and $200$ km s$^{-1}$ due to cosmic expansion along the line of sight. When the dashed line intersects with the solid line (i.e., $v_\alpha + v_{\text{pe}} = 0$), the photons face a sudden increase in the optical depth due to resonance opacity (see panel e). The
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Physical quantities of the IGM relevant to Lyα transmission along the −y direction from the center of galaxy #0001 at z = 7. Panel a: gas density normalized by the cosmic mean. The range within the virial radius of the source galaxy is plotted as a dotted line and is not included in the optical depth calculation. This also applies to panels b, c, and d. Panel b: ionization photon density. Panel c: neutral hydrogen number density. Panel d: peculiar velocity of IGM w.r.t. to the galaxy is shown as the black solid line, and the green dotted lines are the model fit, $v_{\text{pe,m}} = \sqrt{GM_h/s}$. The blue, black, and red dashed lines show the frequency of the photons that were emitted at $v_{\alpha,i} = -200, 0, \text{and } 200 \text{ km s}^{-1}$, respectively. We plot minus the photon frequency, $-v_{\alpha}$, so that the intersection between the dashed and solid lines (i.e., $v_{\alpha} + v_{\text{pe}} = 0$) is where the photons encounter resonance in the rest frame of the IGM. Panel e: optical depth of the IGM for four different frequencies. The line colors correspond to the same frequencies as for the dashed lines of panel d plus the magenta line corresponds to $v_{\alpha,i} = 600 \text{ km s}^{-1}$. Panel f: the neutral hydrogen density map around the line of sight. The color scheme is the same as the one used in Figure 1. The half-circle in the black dashed line indicates the virial radius. The contours lines describe the location where $v_{\alpha} + v_{\text{pe}} = 0$ for the frequencies plotted in panels d and e. The thickness of the contours depicts the thermally broadened line-absorption profile of photoionized gas with T = $10^4$ K. The vertical dotted lines extending across panels d, e, and f connect the locations where $v_{\alpha} + v_{\text{pe}} = 0$.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Probability distribution of the maximum infall velocity $v_{in}$ for 2000 sight lines for galaxies #0001, #0014, and #0250 shown in black, blue, and red, respectively. $v_{in}$ is obtained by identifying the location of the cutoff in the Ly$\alpha$ transmission. The galaxy masses are shown in the legend.</figcaption></figure>

photons that started at $v_{\alpha,i} = 0$ and $200 \text{ km s}^{-1}$ would never have encountered the resonance cross section without the infall motion.

Another reason for the red-side absorption is that the infalling gas near the virial radius remains opaque to Ly$\alpha$ photons despite the intense ionizing radiation from the source. From $s = 10R_{vir} \approx 1.8 \text{ h}^{-1} \text{ cMpc}$ to $s = R_{vir} \approx 0.18 \text{ h}^{-1} \text{ cMpc}$, the gas density in panel a rises by 10 times the cosmic mean while the ionizing radiation intensity ($\Gamma$) rises roughly by a factor of a hundred (panel b). The gas temperature remains close to 20,000 K throughout the interval. The 10 times increase in density cancels the hundred times increase in the intensity because the neutral hydrogen density of highly ionized gas ($n_{HI} = x_{HI}n_{H}$) goes as $n_{H}(n_{HI}/\Gamma) \approx n_{H}^{2}/\Gamma$ in the ionization equilibrium for constant gas temperature. Therefore, the neutral hydrogen density does not fall toward the galaxy as shown in panel c, keeping the IGM optically thick to Ly$\alpha$.

The transmission cutoff frequency is set by the maximum infall velocity outside the virial radius $R_{vir}$. In most cases, the infall velocity peaks around $s = R_{vir}$ (panel d), where the infall velocity is close to the circular velocity at the virial radius, $V_c \equiv \sqrt{GM_h/R_{vir}}$. Here, the galaxy mass, $M_h$, is the total mass that includes both baryonic and dark matter within $R_{vir}$. For each transmission curve, we measure the maximum infall velocity, $v_{in}$, by identifying the location of the transmission cutoff, which coincides with the peak of the transmission curve slope. We plot in Figure 5 the probability distribution of $v_{in}$ based on the results from the 2000 sightlines. The $v_{in}$ distribution has a fairly large scatter but shows a clear galaxy mass dependence. $v_{in}$ is distributed around 300, 200, and $100 \text{ km s}^{-1}$ for galaxies with their masses 1.09, 0.36, and $0.09 \times 10^{12} M_{\odot}$, respectively.

We show the complete statistics of $v_{in}$ from the entire galaxy sample as a function of halo mass in Figure 6 and as a function of UV magnitude in Figure 7. The error bar for each galaxy denotes the $1\sigma$ sight-line variation calculated from 2000 sight lines for each galaxy. $v_{in}$ is similar to $V_c$ throughout the mass range of our sample, although it is slightly lower at high mass. For example, $v_{in}$ prefers $V_c$ at $M_h = 10^{11} M_{\odot}$ while it prefers $0.85V_c$ at $M_h = 10^{12} M_{\odot}$.

We show in Figure 7 that $v_{in}$ also correlates with the UV magnitude. This is due to the strong correlation between $M_{UV}$ and $M_h$, which hardly evolves in time (see Figure 8). We find that the correlation coefficient between $M_{UV}$ and $M_h$ is above 95% for the galaxy samples used in this work. The relation between $v_{in}$ and $M_{UV}$ is well described by a fitting function

$$v_{in,fit} \approx [160 \text{ km s}^{-1}]e^{-0.19[M_{UV}+20]}, \quad (8)$$

with a $1\sigma$ scatter of $\sim 40 \text{ km s}^{-1}$ (see the cyan dashed line in Figure 7).

Nonvanishing transmission on the blue side of Ly$\alpha$ often appears at $z = 6$, while it is mostly zero at $z = 7$ and 8. At $z = 7$, only bright galaxies with $M_{UV} \lesssim -23$ can transmit a small fraction of blue-side photons. According to Gronke et al. (2021), this is because the typical H I density in ionized IGM is too high ($\gtrsim 10^{-9} \text{ cm}^{-3}$) to transmit the blue-side photons at $z > 6.5$. We refer readers to their work for a detailed analysis on this subject.

We note that the residual IGM neutral fraction at the postoverlap phase ($z < 6.2$) is a factor of $\sim 5$ lower than the estimate of Fan et al. (2006). We checked for a subsample of halos that increasing the calculated optical depths by a similar factor did not affect the transmission redward of the cutoff.

## 3.2. Damping-wing Absorption

On the red side of the transmission cutoff, photons do not encounter resonance absorption, leading to a much higher transmissivity. Here, the opacity comes from the damping-wing absorption from highly neutral segments of sight lines. The opacity contribution from each segment depends sensitively on the distance to the segment from the source: $d\tau_{DW} \propto s^{-2}$, as shown in Equation (5).

The red side of the transmission curves in Figure 3 show the damping-wing opacity for different sight lines, halos, and redshifts. This red-side transmission decreases toward high redshift due to increasing damping-wing opacity. For example, the median transmission at $v_{\alpha} = 400 \text{ km s}^{-1}$ for galaxies #0001, #0016, and #0150 is nearly 100% at $z = 6$, while it drops to 80, 53, and 62% at $z = 7$ and to 44, 42, and 37% at $z = 8$, respectively. We find a substantial sight-line-to-sight-line variation in $\tau_{DW}$ for individual galaxies at $z > 6$. The $1\sigma$ range of the transmission bracketed by the cyan lines is up to 30%, implying that the ionization field surrounding the source galaxies is highly anisotropic.

### 3.2.1. H II Bubble Size Estimation

If IGM is fully ionized inside the H II bubble and fully neutral outside, the damping-wing opacity, $\tau_{DW}$, is given by integrating Equation (5) from the bubble size $s_{min} = R_b$ to the end of reionization $s_{max} = s(z = 6)$. Then, $\tau_{DW}$ would depend only on the bubble size $R_b$ and the reionization history $\bar{x}_{HI}$ and one can constrain $\bar{x}_{HI}$ and $R_b$ by measuring $\tau_{DW}$ from observations.

Often, an approximate expression for Equation (5) is used for convenience. In that case, one assumes that the opacity comes entirely from the neighborhood of the source galaxy before the redshift term in the integral changes significantly (Dijkstra 2014). The redshift-dependent terms can then be
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Maximum infall velocity $v_{in}$ as a function of galaxy mass $M_h$. Error bars for each data point show the 1σ range computed for individual galaxies from 2000 sight lines.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Maximum infall velocity $v_{in}$ as a function of UV magnitude $M_{UV}$. Error bars for each data point show the 1σ range computed for individual galaxies from 2000 sight lines. The cyan dashed line is a fitting function given by Equation (8).</figcaption></figure>

taken out from the integral, giving

$$ \tau_{\text{DW, app}}(v_{\alpha}) = 5.7 \left[ \frac{z_g + 1}{8} \right] \bar{x}_{\text{HI}}(z_g) \left[ \frac{R_b + \Delta s(v_{\alpha})}{h^{-1} \text{cMpc}} \right]^{-1} \times \left[ \frac{h}{0.678} \right] \left[ \frac{\Omega_b}{0.048} \right] \left[ \frac{\Omega_m}{0.307} \right]^{-1}, \quad (9) $$

where $z_g$ is the redshift of the source galaxy. In this case, the opacity does not depend on the reionization history along the line of sight, and one can trivially convert the measured $\tau_{\text{DW}}$ to $\bar{x}_{\text{HI}}$, given a relation between $\bar{x}_{\text{HI}}$ and $R_b$ computed from reionization models (e.g., Furlanetto & Oh 2005).

In this work, we define $R_b$ to be the value that matches the damping-wing optical depth $\tau_{\text{DW}}$ with the computed IGM transmissivity at $v_{\alpha} = 400 \text{ km s}^{-1}$ (hereafter $T_{400}$). $v_{in}$ is smaller than $400 \text{ km s}^{-1}$ for all sight lines considered in this work. Thus, the optical depth at $400 \text{ km s}^{-1}$ is purely from the damping-wing opacity.

In Figure 9, we show $R_b$ for galaxies #0001 and #0100 at $z=7$ as the line contours overplotted on the H I density maps. We computed $R_b$ for 360 equally spaced directions ($\hat{n}$) on the

xy plane and connected the locations of $R_b \hat{n}$ to create a closed contour. The green solid line contour is from the original expression of Equation (5), while the dashed line contour is from the approximate expression of Equation (9).

We note that the original and approximate expressions generally give similar results, although the approximate value occasionally overshoots in some directions. The overshoot tends to be larger when $R_b$ is larger (i.e., $\tau_{\text{DW}}$ is small). This is because assuming the $[1+z]\bar{x}_{\text{HI}}(z)$ term in Equation (5) to be constant overestimates $\tau_{\text{DW}}$ by ~0.1 at $z=7$, due to the decline of $\bar{x}_{\text{HI}}$ toward low redshift. For large H II bubbles like the one surrounding galaxy #0001, the error in $\tau_{\text{DW}}$ can bias $R_b$ more because $\tau_{\text{DW}}$ is generally smaller. However, galaxy #0001 is a rare extreme case, and most of the galaxies are in smaller H II bubbles like galaxy #0100, where the error is insignificant.

The $R_b$ contour generally follows the H II bubble outline with some deviations that make the contour spikier than the actual bubble shape. The deviation between the apparent and effective bubble sizes largely owes to the bubble geometry being more complex than what is assumed in Equation (5): “fully ionized inside the H II bubble and fully neutral outside.” When another
<figure><img src="image_10.png" /><figcaption>Figure 8. Trend of $M_{UV}$ with $M_h$. Black dots are for individual galaxies in the $z=7$ snapshot. The red solid and blue dashed lines are the median trend for $z=7$ and 8, respectively.</figcaption></figure>

large H II bubble is located in the vicinity of the H II bubble surrounding the source galaxy, for example, $R_b$ would be larger than what Equation (5) gives toward the direction of the neighboring bubble.

Another important factor ignored in Equation (5) is self-shielded systems, which are small and dense pockets of neutral gas in mini halos within H II regions. When a sight line encounters such a clump, it can pick up a substantial amount of opacity even inside a H II region. The $-y$ direction ($\phi = 270^\circ$) from galaxy #0001 (upper panels of Figure 9) is an example that is affected by a clump. This clump is described in detail in Figure 4, where it is shown as a blue spot in panel f at $s \approx 2.3 h^{-1}$ cMpc. The baryon density in this neutral clump reaches 30 times the cosmic mean at the peak (panel a). The optical depth for the photon emitted at $v_{\alpha,i} = 600$ km s$^{-1}$ rises to 0.3 at this clump (magenta line in panel e). $R_b$ for this sight line is $0.33 h^{-1}$ cMpc, although the neutral region actually starts at $s \approx 4 h^{-1}$ cMpc.

We present the $R_b$ statistics of the entire sample with $M_{UV} < -18.7$ at $z=7$ and 8 in Figure 10. We also give the corresponding transmissivity at $v_\alpha = 400$ km s$^{-1}$ on the right-hand side y-axis. We also tabulate the median values of $T_{400}$ for $M_{UV} = -21, -20, \text{and}, -19$ in Table 1. The scattered data points with vertical bars in Figure 10 highlight the large galaxy-to-galaxy and sight-line-to-sight-line variation in the damping-wing opacity; the $1\sigma$ range depicted by the gray shape shows is larger than the median value.

$R_b$ is expected to evolve with redshift, which is considered to be the main cause of reduced LAE visibility at high redshifts. The data points of Figure 10 shows this expected redshift dependence of $R_b$ and $T_{400}$ clearly: the median values of $R_b$ and $T_{400}$ are around $5 h^{-1}$ Mpc and 73% at $z=7$, while $2.5 h^{-1}$ Mpc and 37% at $z=8$, respectively.

The $M_{UV}$ dependence is not easily visible from the data points, but the median curve reveals a weak dependence. At $z=7$, the median $R_b$ is $5.5 h^{-1}$ Mpc for galaxies $M_{UV} = -21$, whereas it is about 1 Mpc smaller for the galaxies with $M_{UV} = -19$. Similarly, at $z=8$, $R_b$ is about $3 h^{-1}$ Mpc for $M_{UV} = -21$ and $2 h^{-1}$ Mpc for $M_{UV} = -19$. The corresponding $T_{400}$ drops from 73% to 70% at $z=7$ in the same $M_{UV}$ interval, while from 37% to 30% at $z=8$.

### 3.3. Self-shielded Systems

Self-shielded systems can add a large opacity to certain sight lines to create a large discrepancy between the effective and actual bubble size as mentioned in Section 3.2. In this section, we estimate the quantitative impact of these systems on transmissivity by removing these systems in the transmissivity calculation. For this, we first identify highly neutral segments of sight lines with $n_{HI} > 10^{-5}$ cm$^{-3}$. We then assume all those segments shorter than $0.3 h^{-1}$ cMpc are the compact neutral clumps within the H II region and lower the H I density of those segments to $10^{-5}$ cm$^{-3}$ when evaluating Equation (1). This way, we exclude the damping-wing opacity of the self-shielded system like the one in Figure 4.

In Figure 11, we compare $T_{400}$ before and after removing the self-shielded systems for 2000 sight lines of galaxies #0001, #0021, and #0100. The data points that are above the $x=y$ line are the sight lines that are affected by the neutral clumps. The figure shows that some sight lines have substantially higher transmission when clumps are removed.

We find that removing the self-shielded clumps has a small impact on the global statistics of the transmission. Removing the clumps increases the median value of $T_{400}$ for galaxies #0001, #0021, and #0100 increase from 80.1% to 80.9%, from 82.9% to 83.5%, and from 62.6% to 63.8%, respectively. Also, only 0.78%, 0.69%, and 2.3% of the sight lines experience more than a 10% increase in $T_{400}$, respectively.

However, if we limit the statistics to low-transmission sight lines with $T_{400} < 0.3$, the impact appears much stronger, finding that $T_{400}$ increases by more than 10% for 33%, 100%, and 12% of the sight lines, respectively. This result suggests that unusual nondetections of Ly$\alpha$ in UV bright galaxies can be partially attributed to the self-shielded systems blocking the sight line.

We note that our method of removing self-shielded systems is not perfectly accurate. For example, if the edges of two adjacent H II bubbles are closer than $0.3 h^{-1}$ cMpc, our method can erroneously regard the thin H I region between the edges as an isolated neutral clump. However, such complications on the bubble edge would have small impacts on the opacity compared to the impact of dense neutral clumps within the edge of a H II region, which is generally closer to the source. We leave more rigorous and comprehensive analyses of the self-shielded systems for future studies.

### 3.4. Transmissivity of the Ly$\alpha$ Emission Line

We adopt a model spectral energy distribution for the Ly$\alpha$ emission line profile from the CGM to estimate the quantitative impact of IGM absorption. The adopted model is

$$F_{\text{in}}(v_{\alpha}) \propto \exp\left(-\left[\frac{v_{\alpha} - V_c}{V_c/2.355}\right]^2\right). \quad (10)$$

This model assumes that the profile is given by a Gaussian which the peak location and the FWHM are given by the circular velocity of the galaxy $V_c = \sqrt{GM_h/R_{vir}}$ as motivated by recent observations of LAEs at $z<6$ (Yang et al. 2016; Verhamme et al. 2018).$^{12}$ Here, we only model the ratio between the intrinsic and transmitted flux for an arbitrary

<footnote>$^{12}$ We discuss how varying the intrinsic profile affects our results in Section 4.4 based on the experiment results in the Appendix.</footnote>
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Effective H II bubble size $R_b$ for directions on the xy plane is shown on top of the H I density field around galaxies #0001 (upper panel) and #0100 (lower panel) in the $z=7$ snapshot. The green solid line contour describes the $R_b$ calculated from the original expression (Equation (5)), while the blue dashed line is from the approximate expression (Equation (9)). The left panels show the maps in a 32 $h^{-1}$ cMpc box, while the right panels show in a 2 $h^{-1}$ cMpc box. The black dotted lines with angles in degree mark the sight lines described in Figure 12.</figcaption></figure>

normalization for the intrinsic flux. Also, we do not consider any emission from the blue side of Ly$\alpha$ as the transmission is nearly zero there at $z=7$ and 8. We note that, in reality, the spectral shape of the Ly$\alpha$ emission entering the IGM is expected to vary from sight line to sight line due to the complex dynamics of the CGM (e.g., Verhamme et al. 2012; Song et al., in preparation). We take the product of the transmission curve and the intrinsic spectral shape to obtain the IGM transmitted flux: $F_{\text{out}}(v_\alpha) = T(v_\alpha)F_{\text{in}}(v_\alpha)$. Then, the transmitted fraction of the line emission is given by

$$X_{\text{Ly}\alpha} \equiv \frac{\int dv_{\alpha} T(v_{\alpha}) F_{\text{in}}(v_{\alpha})}{\int dv_{\alpha} F_{\text{in}}(v_{\alpha})}. \qquad (11)$$

### 3.4.1. Notable Cases

We first visually inspect the transmission curves and pick several individual that give useful insights into IGM transmissivity. In Figure 12, we plot $F_{\text{in}}$, $F_{\text{out}}$, and $T$ for several notable sight lines on the xy planes of galaxies #0001 and #0100 at $z=7$. We show these sight lines as black dotted lines in Figure 9.

In Figures 12(a) and (b), We compare two sight lines in almost opposite directions ($\phi = 5^\circ$ versus $190^\circ$) to galaxy #0100 at $z=7$. Although the former case has a larger distance to the neutral region ($R_b = 11$ versus $3.7 h^{-1}$ cMpc), it has a lower transmitted fraction $X_{\text{Ly}\alpha}$ (35% versus 46%) because it has a higher infall velocity (183 versus 114 km s$^{-1}$). This comparison highlights that the infall velocity is another crucial factor for the IGM transmission and its sight-line variation.
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Effective H II bubble size $R_b$ for sample galaxies in the $z = 7$ (left) and 8 (right) snapshots. Vertical bars on the data points describe the sight-line variation of individual galaxies. On the right y-axis, we show the corresponding IGM transmission at $v_\alpha = 400$ km s$^{-1}$. The black solid line marks the median of $R_b$ at given $M_{UV}$. The gray shade brackets the $1\sigma$ range of the combined sight-line-to-sight-line and galaxy-to-galaxy variation. Individual data points are shown only for relatively rare bright galaxies to avoid crowding the figure with too many data points.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 11.** Transmission at $v_\alpha = 400$ km s$^{-1}$ before and after removing the self-shielded systems in the opacity integral of Equation (1) as explained in Section 3.3. The results for 2000 sight lines of galaxies #0001, #0021, and #0100 are shown in the left, middle, and right panels, respectively.</figcaption></figure>

In Figures 12(b) and (c), we compare two sight lines that are only 3° apart ($\phi = 190°$ versus $193°$). Despite the small difference in direction, the latter case has a much stronger damping-wing absorption ($R_b = 3.7$ versus $0.26$ h$^{-1}$ cMpc), which results in a more than four times lower transmission ($X_{Ly\alpha} = 46\%$ versus $11\%$). This dramatic difference is due to a self-shielded system that is located in the direction of $\phi = 193°$ from the galaxy (lower-right panel of Figure 9).

Figure 12(d) shows an interesting case in which the transmitted spectrum is double peaked even though the intrinsic spectrum is single peaked. This “pseudo” blue peak occurs due to a leakage feature in the transmission curve, which appears as a sharp peak in the inset panel. This sight line is in the $\phi = 24°$ direction from galaxy #0001 (upper-right panel of Figure 9), which goes through an extremely low H I density blob with $n_{HI} \lesssim 10^{-10}$ cm$^{-3}$ near the virial radius shown as a dark red spot on the H I density map. We find that the temperature of this blob is $\sim 10^6$ K, implying that the blob was created by a supernova explosion in a neighboring galaxy. This double-peak feature appears for other adjacent sight lines that go through the optically thin spot. We find similar features for sight lines between $\phi = 15°$ and $30°$ for this galaxy on the xy plane. From

<table><caption>**Table 1**
Median IGM Transmissivity at $v_\alpha = 400$ km s$^{-1}$</caption><thead><tr><th></th><th></th><th></th><th></th></tr></thead><tbody><tr><td>$M_{UV}$</td><td>-21</td><td>-20</td><td>-19</td></tr><tr><td>$T_{400}(z = 7)$</td><td>73.1%</td><td>70.6%</td><td>70.3%</td></tr><tr><td>$T_{400}(z = 8)$</td><td>37.2%</td><td>33.3%</td><td>30.0%</td></tr></tbody></table>

visual inspections of several other galaxies, we find a few percent of sight lines exhibit similar features. We leave a more comprehensive analysis of this pseudo-double-peak feature for future studies.

### 3.4.2. Global Statistics

We calculate $X_{Ly\alpha}$ for 2000 sight lines for galaxies with $M_{UV} \lesssim -19$. We show the full statistics of the results in Figure 13, where we show $X_{Ly\alpha}$ as a function of $M_{UV}$ with the galaxy-to-galaxy and sight-line-to-sight-line variations. Similarly to Figure 10, the error bars on the data points show the 68% range for individual galaxies. The solid line marks the median value for the corresponding $M_{UV}$, and the gray shade covers the 68% range with the galaxy-to-galaxy and sight-line-to-sight-line variations.
<figure><img src="image_14.png" /><figcaption>**Figure 12.** Spectral shape of the model’s intrinsic Lyα emission from the circumgalactic medium (dotted line) and transmitted emission through the intergalactic medium (solid line) along with the transmission curve in the inset panel for several sight lines on the xy planes of galaxies #0001 and #0100 at z = 7. Panels a, b, and c are for the azimuthal angles of φ = 5°, 190°, and 193° from galaxy #0100, respectively, and Panel d is for 24° from galaxy #0001, respectively (see Figure 9). The maximum infall velocity $v_{\text{in}}$, H II bubble size $R_b$, and Lyα line transmissivity $X_{\text{Ly}\alpha}$ for each sight line are given in the corresponding panel.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 13.** Lyα transmissivity $X_{\text{Ly}\alpha}$ for sample galaxies in the z = 6 (left), 7 (middle), and 8 (right panel) snapshots. Each data point with an error describes the sight-line variation of an individual halo. The black solid line marks the median of $X_{\text{Ly}\alpha}$ at given $M_{\text{UV}}$. The gray shades bracket the 1σ range of the combined sight-line-to-sight-line and galaxy-to-galaxy variation. Individual data points are shown only for relatively rare bright galaxies to avoid crowding the figure with too many data points.</figcaption></figure>
<table><caption>**Table 2**
Median and 68% Range of IGM Transmissivity $X_{Ly\alpha}$</caption><thead><tr><th>M<sub>UV</sub></th><th>-22</th><th>-21</th><th>-20</th><th>-19</th></tr></thead><tbody><tr><td>z = 6</td><td>58.8<sup>+17.9</sup><sub>-9.9</sub>%</td><td>59.3<sup>+19.8</sup><sub>-22.8</sub>%</td><td>55.7<sup>+22.8</sup><sub>-23.1</sub>%</td><td>49.8<sup>+26.1</sup><sub>-23.5</sub>%</td></tr><tr><td>z = 7</td><td>44.5<sup>+12.5</sup><sub>-12.1</sub>%</td><td>38.4<sup>+19.3</sup><sub>-16.7</sub>%</td><td>32.2<sup>+18.8</sup><sub>-14.8</sub>%</td><td>27.5<sup>+19.9</sup><sub>-15.1</sub>%</td></tr><tr><td>z = 8</td><td>N/A</td><td>19.2<sup>+10.6</sup><sub>-8.8</sub>%</td><td>12.3<sup>+10.3</sup><sub>-7.9</sub>%</td><td>7.3<sup>+9.2</sup><sub>-5.2</sub>%</td></tr></tbody></table>

variations combined. We list the 1σ range for $M_{UV} = -22, -21, -20$, and $-19$ at $z = 6, 7$, and 8 in Table 2.

First, $X_{Ly\alpha}$ decreases toward high $z$. For example, the median at $M_{UV} = -21$ is 59.3%, 38.4%, and 19.2% at $z = 6, 7$, and 8, where the IGM is 99.9%, 50%, and 13.2% ionized, respectively. This is in accordance with the common expectation that $X_{Ly\alpha}$ is a useful probe of the IGM neutral fraction.

It is also notable that $X_{Ly\alpha}$ is well below 100% at $z = 6$, where the IGM is nearly 100% ionized. The IGM absorption at $z = 6$ should be entirely due to the resonance scattering by infalling gas given that the IGM is nearly fully ionized. Thus, we take the $z = 6$ results as the fully ionized limit and calculate the relative change in $X_{Ly\alpha}$ at higher redshifts. That is, we take the ratio of the median at $z = 7$ or 8, $\langle X_{Ly\alpha} \rangle_{z=7}$ or 8, to that at $z = 6$, $\langle X_{Ly\alpha} \rangle_{z=6}$, in order to separate out the damping-wing opacity from the neutral IGM. For example, the median transmission fraction of galaxies with $M_{UV} = -21$ is suppressed to 65% and 32% at $z = 7$ and 8, respectively, compared to that at $z = 6$. We list this ratio for $M_{UV} = -22, -21, -20$, and $-19$ in Table 3. Both Figure 13 and Table 3 show that UV-fainter galaxies show a stronger suppression in $X_{Ly\alpha}$ than the brighter ones do. For example, the median transmission fraction of $M_{UV} = -21$ galaxies drops to 35% from $z = 6$ to 8, while that of $M_{UV} = -19$ galaxies drops to 15%.

The $M_{UV}$ dependence of $X_{Ly\alpha}$ is often attributed to the difference in H II bubble size between bright and faint galaxies. Interestingly, we find that the bubble size only partially explains the dependence. According to Table 1, the IGM transmissivity at a fixed frequency of $v_\alpha = 400 \text{ km s}^{-1}$ changed from 73.1% to 70.3% between $M_{UV} = -21$ and $-19$ at $z = 7$. This difference is smaller than the difference in $X_{Ly\alpha}$ between the same $M_{UV}$ (65% versus 55%). This extra difference in $X_{Ly\alpha}$ is attributed to the difference in the velocity offset of transmitted flux (Figure 7), caused by the $M_{UV}$ dependence of the IGM infall velocity $v_{in}$ and the peak location of our intrinsic line profile. That is, brighter galaxies tend to transmit redder photons due to the stronger infall motion, and these redder photons are subject to smaller damping-wing opacity for a given H I density.

We illustrate our understanding of the $M_{UV}$ dependence of $X_{Ly\alpha}$ in Figure 14, where we show the median (transmitted) emission line profile for bright ($M_{UV} < -22$) and faint ($M_{UV} > -19$) galaxy groups at $z = 6, 7$, and 8. We first construct the median transmission curve for each group using the median of $v_{in}$ and $R_b$ for each group. Then, multiply the median transmission by the intrinsic emission profile obtained from Equation (10) using the median circular velocity. The line flux relative to $z = 6$ drops to 70% (37%) for the bright group and to 58% (12%) for the faint groups at $z = 7$ (8). The figure shows both bubble size and wavelength of the transmitted photons create an $M_{UV}$ dependence in the damping-wing opacity as we explained above.

We also highlight the large scatter in the $X_{Ly\alpha} - M_{UV}$ relation. Due to the complicated nature of $Ly\alpha$ radiative processes,

<table><caption>**Table 3**
IGM Transmissivity Relative to $z = 6$</caption><thead><tr><th>M<sub>UV</sub></th><th>-22</th><th>-21</th><th>-20</th><th>-19</th></tr></thead><tbody><tr><td>⟨X<sub>Lyα</sub>⟩<sub>z=7</sub>/⟨X<sub>Lyα</sub>⟩<sub>z=6</sub></td><td>77%</td><td>65%</td><td>58%</td><td>55%</td></tr><tr><td>⟨X<sub>Lyα</sub>⟩<sub>z=8</sub>/⟨X<sub>Lyα</sub>⟩<sub>z=6</sub></td><td>N/A</td><td>32%</td><td>22%</td><td>15%</td></tr></tbody></table>

galaxies at the same redshift with the same brightness can transmit a significantly different number of $Ly\alpha$ photons from sight line to sight line and from galaxy to galaxy. As $X_{Ly\alpha}$ decreases toward low brightness and high redshift, the variation in $X_{Ly\alpha}$ becomes larger relative to the median. The 1σ range of $X_{Ly\alpha}$ for $M_{UV} = -21$ at $z = 7$ is $59.3^{+19.8}_{-22.8}\%$, while that for $M_{UV} = -20$ at $z = 8$ is $\bar{X}_{Ly\alpha} = 7.3^{+9.2}_{-5.2}\%$ (see Table 2). The ratio of the 1σ range to the median is $43/59 \approx 0.7$ in the former case, while $14.4/7.3 \approx 2$ in the latter case. This means that one needs to measure $Ly\alpha$ line strengths of at least 50 (400) galaxies to constrain $\langle X_{Ly\alpha} \rangle$ to 10% in the former (latter) case assuming perfectly accurate measurements.

## 4. Summary and Discussion

We have computed the IGM transmission for the $Ly\alpha$ emission line from star-forming galaxies during reionization from the Cosmic Dawn II simulation to understand what aspects of the IGM determine the transmission. Thanks to the star formation physics and the subsequent radiative transfer implemented in a large box with high resolution, we were able to compute the IGM transmissivity for galaxies with $M_{AB1600} \gtrsim -22$ while accounting for the small-scale physics around the virial radii of the galaxies. We find that the IGM infall motion, H II bubble size, and self-shielded neutral gas systems are the critical factors. We discuss the implications of our results for the current understandings of LAEs below.

### 4.1. Significance of Infall Motion for the IGM Transmissivity

Previous works have repeatedly reported the impact of the infalling IGM on the $Ly\alpha$ transmission (Santos 2004; Dijkstra et al. 2007; Iliev et al. 2008; Dayal et al. 2011; Laursen et al. 2011; Bolton & Haehnelt 2013; Sadoun et al. 2017; Weinberger et al. 2019; Garel et al. 2021). Although galactic outflow from supernova explosions may extend beyond the virial radius for some occasions, the outflow is confined within the virial radius in most sight lines (e.g., Shen et al. 2013). Thus, the peculiar motion of IGM outside the virial radius is dominated by the gravitational field of the source galaxy. The gravitationally infalling IGM suppresses a significant portion of the red-peak emission even after reionization has ended, requiring higher intrinsic $Ly\alpha$ luminosity to explain the observed $Ly\alpha$ luminosity function (Dayal et al. 2011). The absorption by infalling IGM extends to $v_\alpha \sim 100 \text{ km s}^{-1}$ in Laursen et al. (2011) and Garel et al. (2021), for example, while it goes to $v_\alpha \gtrsim 200 \text{ km s}^{-1}$ in our work and Weinberger et al. (2019). This difference likely resulted from larger simulation box sizes of the latter two studies, which involve more massive galaxies with stronger infall motion.

The truncation frequency of the IGM transmission is given by the maximum infall velocity along the line of sight. In most sight lines, the infall velocity peaks around the virial radius of the galaxy. If one uses a distance larger than the virial radius to define the IGM, the maximum infall velocity would decrease. In
<figure><img src="image_16.png" /><figcaption>**Figure 14.** The median emission-line profile (solid line) calculated by multiplying the median transmission (dotted line) to the intrinsic profile (dashed line). The red and blue lines show the results for the bright ($M_{UV} < -22$) and faint ($M_{UV} > -19$) galaxy groups, respectively. The left, middle, and right panels show the results for $z = 6, 7$, and 8 galaxies, respectively. The median transmission curve is constructed using the median infall velocity $v_{\text{in}}$ and effective H II bubble size $R_{\text{b}}$ for each group, the values of which are shown on the plot. The intrinsic profile is given by Equation (10), where the circular velocity $V_c$ is given by the median $M_h$ of each group. The integrated line flux at $z = 7$ or 8 with respect to $z = 6$ is also shown on the plot for each group in percent.</figcaption></figure>

this case, some of the absorption by the infalling IGM in this work would be regarded as the absorption by the CGM. This is illustrated in Figure 6 of Weinberger et al. (2019), where they define gas at $R_{\text{vir}} < s < 5R_{\text{vir}}$, $5R_{\text{vir}} < s < 10R_{\text{vir}}$, and $s > 10R_{\text{vir}}$ as the “inner CGM,” “outer CGM,” and “IGM,” respectively. All three components are regarded as the IGM in this work.

We note that AGNs, which the CoDa II simulation does not include, can result in a much higher transmission owing to much stronger outflow and ionizing radiation. Jets from AGNs can disrupt the infall motion well beyond the virial radius, and the H I density is low enough to be transparent at resonance ($<10^{-10} \text{ cm}^{-3}$) within the proximity zones. Indeed, low-$z$ observations show that luminous LAEs are often associated with AGNs (Konno et al. 2016). Neighboring galaxies within the proximity zone can also receive a boost in the Ly$\alpha$ transmission as was reported by Bosman et al. (2020). Also, the X-ray from AGNs can change the reionization morphology by partially ionizing the surrounding H I region (Kakiichi et al. 2017).

The relation between $v_{\text{in}}$ and the galaxy mass (or UV magnitude) does not evolve much with redshift as the infall motion is purely gravitational in nature. However, $v_{\text{in}}$ can still affect the IGM transmission by setting the minimum wavelength of transmitted photons. Stronger infall motion results in the transmission of redder photons, which experience less damping opacity. As a result, Ly$\alpha$ emission from brighter galaxies, which are likely more massive, are less subject to increasing neutral IGM fraction toward higher redshifts.

The infall motion varies depending on the viewing angle due to the complexity of density structures around the source galaxy. This variation adds to the statistical uncertainty of the transmitted flux from a single observation. We find the infall velocity typically has a scatter of $\sim 40 \text{ km s}^{-1}$ from sight line to sight line (Figure 7). The sight-line difference in the infall velocity can overwhelm that in the damping-wing opacity in determining the total IGM transmissivity (e.g., Figure 12(a) versus (b)).

## 4.2. Observational Signature of the Infall Motion

Truncation of the red peak by the infalling IGM can make the spectral shape of the transmitted flux blue-skewed in frequency even if the intrinsic emission is symmetric

(Figure 12). Such blue-skewed spectra are often observed from recent high-$z$ spectroscopic LAE surveys (see, e.g., Figure 4 of Jung et al. 2020).

The infall motion is a potential explanation for the offset between the observed line and the systemic redshift of the source galaxy. The absorption by infalling opaque IGM will force the peak to be on the red side of $v_{\alpha} = v_{\text{in}} \approx V_c$ regardless of the intrinsic emission. Given the correlation between the UV magnitude and the infall velocity (Figure 7), the offset is expected to correlate with the UV magnitude, and this is supported by existing observations (e.g., see Figure 10 of Hashimoto et al. 2019).

The pseudo-blue-peak feature described in Figure 12(d) is an interesting possible outcome of infalling IGM. Given that the IGM is unlikely to transmit photons from the blue side of Ly$\alpha$ at $z \gtrsim 6.5$ (Figure 3; also see Gronke et al. 2021), the currently existing double-peak cases at $z \sim 6.5$ without exact systemic redshifts (Matthee et al. 2018; Songaila et al. 2018; Meyer et al. 2021) may be this case. The line shape of our case is the most similar to NEPLA4 of Songaila et al. (2018, see their Figure 12), although the overall shape is about two times narrower in our case. The peak separation of double-peaked LAEs ranges from $100$ to $300 \text{ km s}^{-1}$, which is narrower than the reported value of $\sim 750 \text{ km s}^{-1}$ from $z = 2-3$ LAEs by Kulas et al. (2012). The pseudo-blue-peak scenario can explain the narrow separation because the two peaks come from incomplete resonance absorption of a single red peak. We also note that more recent observations with improved spectral resolution by Matthee et al. (2021) do find lots of narrowly separated double peaks at $z = 2-3$. If any of the reported $z > 6$ double-peaked LAE falls into the pseudo-double-peak case, its systemic redshift will be found blueward of both peaks in future observation.<fnref n="13" />

## 4.3. Observational Impacts of Self-shielded Systems

Before reionization, mini halos below $10^8 M_{\odot}$ could accrete gas from the cold IGM and develop dense gas clumps inside their potential wells (Shapiro et al. 2003, 2004). Upon

<footnote n="13">¹³ Matthee et al. (2018) also suggested a model that produces a double-peak feature redward of the systemic redshift in their Figure 9. However, the suggested H I density configuration is different from our case.</footnote>
reionization, these clumps were exposed to the ionizing background radiation and started evaporating. However, the evaporation took up to several hundred megayears to finish for relatively massive ones (Iliev et al. 2005; Park et al. 2016; D’Aloisio et al. 2020) and, thus, some of the clumps are expected to have existed within H II regions and attenuated Lyα emission for incoming sight lines (Bolton & Haehnelt 2013; Mesinger et al. 2015; Kakiichi et al. 2016). Previous works relied on analytical or subgrid models of these clumps to demonstrate their impact on large-scale Lyα transmission.

We have revisited this subject with the CoDa II simulation, which directly resolves these dense clumps at (proper) kiloparsec scales while simultaneously reproducing the H II bubbles as large as several (proper) Mpc. We report that the self-shielded systems do not seem to have a large impact on the average transmission, but they seem to be responsible for a large fraction of sight lines with unusually low transmission. Thus, Lyα nondetections from UV-bright galaxies can be at least partially attributed to these dense clumps blocking the sight line. We note that the spatial resolution of CoDa II may have underestimated the opacity from these systems in this work given that high-resolution simulation studies (e.g., Nasir et al. 2021) suggest there are a lot of self-shielded absorbers smaller than the cell size of CoDa II ($16 h^{-1}$ ckpc).

The statistics of these clumps depend on multiple parameters including the ionizing background radiation intensity, local overdensity, timing of reionization, X-ray preheating, and baryon dark matter streaming motion (Mesinger et al. 2015; Park et al. 2016; D’Aloisio et al. 2020; Cain et al. 2020; Park et al. 2021). The first three parameters are especially relevant to relatively massive ($\gtrsim 10^6 M_\odot$) systems that can survive the ionizing background radiation for more than $10^7$ yr. The IGM opacity would also depend on these parameters if the unresolved self-shielded systems add a significant amount of opacity to the IGM. Thus, it is worth exploring the Lyman-limit opacity of these absorbers in future studies dedicated to these small-scale structures.

### 4.4. $M_{UV}$ Dependence of IGM Transmissivity

A number of past observations have reported that Lyα emission from relatively brighter galaxies tends to be less affected by the neutral IGM during reionization than fainter ones do (e.g., Fontana et al. 2010; Stark et al. 2011; Pentericci et al. 2011; Ono et al. 2012; Mallerly et al. 2012; Curtis-Lake et al. 2012; Treu et al. 2013; Tilvi et al. 2014; Endsley et al. 2021, I. Jung et al. 2021, in preparation). We confirm this trend from our results in Figure 13 and Table 3). We illustrate our results by calculating the “median emission profile” for bright ($M_{UV} < -22$) and faint ($M_{UV} > -19$) galaxies in Figure 14, which is obtained by applying the median $v_{in}$ and $R_b$ to the intrinsic emission profile created for the median $M_h$.

Based on Figure 14, we conclude that the $M_{UV}$ dependence of the transmission is driven by two factors. First, bright galaxies tend to reside in H II regions larger than where fainter galaxies reside ($R_b = 4.7$ versus $3.5 h^{-1}$ cMpc at $z=7$) because they can ionize larger volumes and are more likely to form in overdense regions, where numerous neighboring galaxies power the H II bubble together (Finkelstein et al. 2019; Endsley et al. 2021). Second, bright galaxies tend to transmit in longer wavelengths because they have a larger velocity offset in their intrinsic emission profile (Mason et al. 2018; Garel et al. 2021;

<table><caption>**Table 4**
68% Range to Median Ratio of $X_{L\alpha}$</caption><thead><tr><th>$M_{UV}$</th><th>-22</th><th>-21</th><th>-20</th><th>-19</th></tr></thead><tbody><tr><td>$z=6$</td><td>0.47</td><td>0.72</td><td>0.82</td><td>1.00</td></tr><tr><td>$z=7$</td><td>0.55</td><td>0.94</td><td>1.04</td><td>1.27</td></tr><tr><td>$z=8$</td><td>N/A</td><td>1.01</td><td>1.48</td><td>1.97</td></tr></tbody></table>

also compare $V_c$ in Figure 14) and higher IGM infall velocity ($v_{in} = 248$ versus $138 \text{ km s}^{-1}$ at $z=7$) suppressing the transmission up to a longer wavelength. Often the former is considered to be the major reason for the $M_{UV}$ dependence, but we find the latter is comparably responsible for the dependence. We also find that the $M_{UV}$ dependence grows stronger toward high $z$. We find the median transmission relative to $z=6$ is 69% versus 57% at $z=7$ while it is 37% versus 11% at $z=8$. This highlights the importance of dealing with the dependence carefully in high-$z$ LAE surveys to avoid biasing the IGM neutral fraction estimation.

The infall velocity $v_{in}$ and the circular velocity $V_c$ stay nearly constant between $z=6$ and 8, and thus do not directly contribute to the IGM transmission evolution. However, the strength of the $M_{UV}$ dependence depends on the wavelength of the transmitted flux set by those two parameters. In the Appendix, we demonstrate the flux-weighted transmission and its $M_{UV}$ dependence can be strengthened or weakened by modest amounts if we shrink or expand the intrinsic profiles, respectively. However, our main conclusion that the IGM transmission for brighter galaxies is less suppressed by the neutral IGM remains solid in this analysis.

### 4.5. Variation in IGM Transmissivity and Its Implication for Observations

While an apparent decline in Lyα visibility at $z>6$ has been interpreted as a sign of an increasing IGM neutral fraction, the complex nature of Lyα radiative processes in the IGM limits the measurement accuracy of $\bar{x}_{HI}$ based on the Lyα emission strength. The IGM infall velocity and size of H II regions vary from sight line to sight line and galaxy to galaxy even for galaxies with the same mass or brightness, and self-shielded systems occasionally block unfortunate sight lines. These uncertainties make it difficult to pin down the average IGM transmissivity at a given redshift with a small number of measurements. Moreover, the intrinsic Lyα emission from the CGM appears to also depend on not only the viewing angle but also various characteristics of galaxies (e.g., Tang et al. 2021; Matthee et al. 2021) as theoretically expected (e.g., Verhamme et al. 2012; Smith et al. 2019, H. Song et al. 2021, in preparation).

We show in Table 4 the ratio between the 68% range and the median value of $X_{L\alpha}$ based on Table 2. The ratio is larger for fainter galaxies at higher redshifts. For galaxies with $M_{UV} \ge -21$ at $z \ge 7$, the ratio is $\ge 1$, implying that we need at least a hundred measurements of Lyα emission strength to constrain the mean IGM transmissivity to 10% accuracy with those galaxies. In practice, the required sample size is likely much larger because most of the measurements would yield only nondetections, and there will be extra uncertainty in the intrinsic emission from the CGM.

Recent state-of-art LAE surveys are starting to constrain the neutral fraction of the IGM at $z \sim 7.5$, but their constraints have
not converged yet: Hoag et al. (2019), Jung et al. (2020), and Whitler et al. (2020) give $\bar{x}_{HI} = 88_{-10}^{+5}\%$, $49_{-19}^{+19}\%$, and $55_{-13}^{+11}\%$, respectively. According to our findings, these studies are based on less than a hundred galaxies, which is not enough to suppress the uncertainty from the IGM transmission variation.

For LAEs at $z \lesssim 7$, ground-based experiments are expected to provide a large-enough sample soon (e.g., SILVERRUSH; Ono et al. 2021). For higher redshifts, future space experiments like the Nancy Grace Roman Space Telescope (NGRST), the James Webb Space Telescope (JWST), and Euclid will substantially improve the constraints with high-quality observations.

Because brighter galaxies have a relatively smaller transmission variation, surveying bright LAEs in a wide field may be more effective than surveying faint galaxies in a narrow field. In this regard, NGRST and Euclid may be more efficient than JWST in constraining the IGM transmissivity at high redshifts, although a definitive conclusion requires accounting for the variation in the intrinsic emission from the CGM.

We thank A. Inoue, M. Ouchi, N. Yoshida, K. Moriwaki, J. Matthee, C. Mason, and A. Smith for helpful comments on this paper. H.P. was supported by the World Premier International Research Center Initiative (WPI), MEXT, Japan and JSPS KAKENHI grant No. 19K23455. H.S. was supported by the Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education (2020R111A1A01069228). Numerical computations were carried out on the gfarm computing cluster of the Kavli Institute for Physics and Mathematics of the Universe and a high-performance computing cluster at the Korea Astronomy and Space Science Institute. T.D. was supported by the National Science Foundation Graduate Research Fellowship Program under Grant No. DGE-1610403. P.R.S. was supported in part by US NSF grant AST-1009799, NASA grant NNX11AE09G, NASA/JPL grant RSA Nos. 1492788 and 1515294, and supercomputer resources from NSF XSEDE grant TG-AST090005 and the Texas Advanced Computing Center

(TACC) at the University of Texas at Austin. I.T.I. was supported by the Science and Technology Facilities Council (grant Nos. ST/I000976/1 and ST/T000473/1) and the South-east Physics Network (SEPNet). K.A. was supported by NRF-2016R1D1A1B04935414 and NRF-2016R1A5A1013277. I.J. acknowledges support from NASA under award number 80GSFC21M0002.

## Appendix

### Dependence of IGM Transmissivity on Intrinsic Emission Profile

We test how sensitive our results are to our choice of the intrinsic emission profile (Equation (10)) by varying the profile in our calculation. We apply the median transmission to bright and faint galaxies in the same way as in Figure 14 in Section 4.4, except that we shrink the profile by reducing the circular velocity term $V_c$ by 33% in one case and inflating the profile by increasing $V_c$ by 50% in the other case. The median infall velocity ($v_{in}$) and effective bubble size ($R_b$) are kept the same as in Figure 14. The results for shrunk and inflated profiles are shown in Figures 15 and 16, respectively.

The flux-weighted transmission, $\langle X_{Lyo} \rangle$, with respect to $z=6$ is 69% (37%) and 57%(11%) for bright and faint galaxies at $z=7$ (8), respectively, in the fiducial case. These fractions decrease to 58% (27%) and 38% (5%) at $z=7$ (8), respectively, for the shrunken intrinsic profile and increase to 71% (41%) and 58% (15%) at $z=7$ (8), respectively, for the inflated profile. Wider emission profiles result in the emission transmitted at longer wavelengths on average, which leads to higher flux-weighted transmission because the transmission increases toward longer wavelengths for a given $R_b$. We also note that the transmission from the faint galaxies is more sensitive to the profile shape, increasing three times from the shrunken to the inflated profile at $z=8$ while that from the bright galaxies increases 1.5 times. Despite these quantitative changes in the relative transmissivity values, the general trend of fainter galaxies being more suppressed by the damping-wing opacity remains throughout our experiment.

<figure><img src="image_17.png" /><figcaption>**Figure 15.** Same as Figure 14 except that we reduce $V_c$ by 33% in Equation (10) to shrink the intrinsic emission profile.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>**Figure 16.** Same as Figure 14 except that we increase $V_c$ by 50% in Equation (10) to inflate the intrinsic emission profile.</figcaption></figure>

## ORCID iDs

Hyunbae Park (박현배) https://orcid.org/0000-0002-7464-7857

Intae Jung https://orcid.org/0000-0003-1187-4240

Paul R. Shapiro https://orcid.org/0000-0002-0410-3045

Ilian T. Iliev https://orcid.org/0000-0002-5174-1365

Kyungjin Ahn https://orcid.org/0000-0003-3974-1239

## References

Adams, J. J., Blanc, G. A., Hill, G. J., et al. 2011, ApJS, 192, 5

Bolton, J. S., & Haehnelt, M. G. 2013, MNRAS, 429, 1695

Bosman, S. E. I., Kakiichi, K., Meyer, R. A., et al. 2020, ApJ, 896, 49

Cain, C., D'Aloisio, A., Iršič, V., McQuinn, M., & Trac, H. 2020, ApJ, 898, 168

Cassata, P., Tasca, L. A. M., Le Fèvre, O., et al. 2015, A&A, 573, A24

Castellano, M., Pentericci, L., Vanzella, E., et al. 2018, ApJL, 863, L3

Curtis-Lake, E., McLure, R. J., Pearce, H. J., et al. 2012, MNRAS, 422, 1425

D'Aloisio, A., McQuinn, M., Trac, H., Cain, C., & Mesinger, A. 2020, ApJ, 898, 149

Dayal, P., Maselli, A., & Ferrara, A. 2011, MNRAS, 410, 830

Deharveng, J.-M., Small, T., Barlow, T. A., et al. 2008, ApJ, 680, 1072

Deparis, N., Aubert, D., Ocvirk, P., Chardin, J., & Lewis, J. 2019, A&A, 622, A142

Dijkstra, M. 2014, PASA, 31, e040

Dijkstra, M., Lidz, A., & Wyithe, J. S. B. 2007, MNRAS, 377, 1175

Endsley, R., Stark, D. P., Charlot, S., et al. 2021, MNRAS, 502, 6044

Fan, X., Carilli, C. L., & Keating, B. 2006, ARA&A, 44, 415

Finkelstein, S. L., D'Aloisio, A., Paardekooper, J.-P., et al. 2019, ApJ, 879, 36

Finkelstein, S. L., Papovich, C., Dickinson, M., et al. 2013, Natur, 502, 524

Fontana, A., Vanzella, E., Pentericci, L., et al. 2010, ApJL, 725, L205

Furlanetto, S. R., & Oh, S. P. 2005, MNRAS, 363, 1031

Garel, T., Blaizot, J., Rosdahl, J., et al. 2021, MNRAS, 504, 1902

Gronke, M., Ocvirk, P., Mason, C., et al. 2021, MNRAS, 508, 3697

Harikane, Y., Ouchi, M., Ono, Y., et al. 2019, ApJ, 883, 142

Hashimoto, T., Inoue, A. K., Mawatari, K., et al. 2019, PASJ, 71, 71

Heinrich, C., & Hu, W. 2021, PhRVD, 104, 063505

Herenz, E. C., Urrutia, T., Wisotzki, L., et al. 2017, A&A, 606, A12

Higuchi, R., Ouchi, M., Ono, Y., et al. 2019, ApJ, 879, 28

Hoag, A., Bradač, M., Huang, K., et al. 2019, ApJ, 878, 12

Hoag, A., Bradač, M., Trenti, M., et al. 2017, NatAs, 1, 0091

Hu, E. M., Cowie, L. L., Barger, A. J., et al. 2010, ApJ, 725, 394

Hu, W., Wang, J., Infante, L., et al. 2021, NatAs, 5, 485

Iliev, I. T., Shapiro, P. R., McDonald, P., Mellema, G., & Pen, U.-L. 2008, MNRAS, 391, 63

Iliev, I. T., Shapiro, P. R., & Raga, A. C. 2005, MNRAS, 361, 405

Inoue, A. K., Hasegawa, K., Ishiyama, T., et al. 2018, PASJ, 70, 55

Jung, I., Finkelstein, S. L., Dickinson, M., et al. 2019, ApJ, 877, 146

Jung, I., Finkelstein, S. L., Dickinson, M., et al. 2020, ApJ, 904, 144

Jung, I., Finkelstein, S. L., Livermore, R. C., et al. 2018, ApJ, 864, 103

Kakiichi, K., Dijkstra, M., Ciardi, B., & Graziani, L. 2016, MNRAS, 463, 4019

Kakiichi, K., Graziani, L., Ciardi, B., et al. 2017, MNRAS, 468, 3718

Kashikawa, N., Shimasaku, K., Malkan, M. A., et al. 2006, ApJ, 648, 7

Kashikawa, N., Shimasaku, K., Matsuda, Y., et al. 2011, ApJ, 734, 119

Khusanova, Y., Le Fèvre, O., Cassata, P., et al. 2020, A&A, 634, A97

Konno, A., Ouchi, M., Nakajima, K., et al. 2016, ApJ, 823, 20

Konno, A., Ouchi, M., Ono, Y., et al. 2014, ApJ, 797, 16

Kulas, K. R., Shapley, A. E., Kollmeier, J. A., et al. 2012, ApJ, 745, 33

Laursen, P., Sommer-Larsen, J., & Razoumov, A. O. 2011, ApJ, 728, 52

Malhotra, S., & Rhoads, J. E. 2004, ApJL, 617, L5

Mallery, R. P., Mobasher, B., Capak, P., et al. 2012, ApJ, 760, 128

Mason, C. A., Fontana, A., Treu, T., et al. 2019, MNRAS, 485, 3947

Mason, C. A., & Gronke, M. 2020, MNRAS, 499, 1395

Mason, C. A., Treu, T., de Barros, S., et al. 2018, ApJL, 857, L11

Matthee, J., Sobral, D., Darvish, B., et al. 2017, MNRAS, 472, 772

Matthee, J., Sobral, D., Gronke, M., et al. 2018, A&A, 619, A136

Matthee, J., Sobral, D., Hayes, M., et al. 2021, MNRAS, 505, 1382

Matthee, J., Sobral, D., Santos, S., et al. 2015, MNRAS, 451, 400

Mesinger, A., Aykutalp, A., Vanzella, E., et al. 2015, MNRAS, 446, 566

Meyer, R. A., Laporte, N., Ellis, R. S., Verhamme, A., & Garel, T. 2021, MNRAS, 500, 558

Morales, A., Mason, C., Bruton, S., et al. 2021, ApJ, 919, 120

Nasir, F., Cain, C., D'Aloisio, A., Gangolli, N., & McQuinn, M. 2021, arXiv:2108.04837

Ocvirk, P., Aubert, D., Chardin, J., Deparis, N., & Lewis, J. 2019, A&A, 626, A77

Ocvirk, P., Aubert, D., Sorce, J. G., et al. 2020, MNRAS, 496, 4087

Ocvirk, P., Gillet, N., Shapiro, P. R., et al. 2016, MNRAS, 463, 1462

Oesch, P. A., van Dokkum, P. G., Illingworth, G. D., et al. 2015, ApJL, 804, L30

Ono, Y., Ouchi, M., Mobasher, B., et al. 2012, ApJ, 744, 83

Ono, Y., Itoh, R., Shibuya, T., et al. 2021, ApJ, 911, 78

Ota, K., Iye, M., Kashikawa, N., et al. 2017, ApJ, 844, 85

Ouchi, M., Harikane, Y., Shibuya, T., et al. 2018, PASJ, 70, S13

Ouchi, M., Shimasaku, K., Akiyama, M., et al. 2008, ApJS, 176, 301

Ouchi, M., Shimasaku, K., Furusawa, H., et al. 2010, ApJ, 723, 869

Park, H., Shapiro, P. R., Ahn, K., Yoshida, N., & Hirano, S. 2021, ApJ, 908, 96

Park, H., Shapiro, P. R., Choi, J.-h., et al. 2016, ApJ, 831, 86

Pentericci, L., Fontana, A., Vanzella, E., et al. 2011, ApJ, 743, 132

Planck Collaboration, Ade, P. A. R., Aghanim, N., et al. 2014, A&A, 571, A16

Planck Collaboration, Aghanim, N., Akrami, Y., et al. 2020, A&A, 641, A6

Sadoun, R., Zheng, Z., & Miralda-Escudé, J. 2017, ApJ, 839, 44

Santos, M. R. 2004, MNRAS, 349, 1137

Santos, S., Sobral, D., & Matthee, J. 2016, MNRAS, 463, 1678

Schenker, M. A., Ellis, R. S., Konidaris, N. P., & Stark, D. P. 2014, ApJ, 795, 20

Shapiro, P. R., Iliev, I. T., & Raga, A. C. 2004, MNRAS, 348, 753

Shapiro, P. R., Iliev, I. T., Raga, A. C., & Martel, H. 2003, in AIP Conf. Ser. 666, The Emergence of Cosmic Structure, ed. S. H. Holt & C. S. Reynolds (Melville, NY: AIP), 89

Shen, S., Madau, P., Guedes, J., et al. 2013, ApJ, 765, 89

Shibuya, T., Ouchi, M., Harikane, Y., et al. 2018, PASJ, 70, S15

Smith, A., Ma, X., Bromm, V., et al. 2019, MNRAS, 484, 39

Song, M., Finkelstein, S. L., Livermore, R. C., et al. 2016, ApJ, 826, 113
Songaila, A., Hu, E. M., Barger, A. J., et al. 2018, ApJ, 859, 91

Stark, D. P., Ellis, R. S., Charlot, S., et al. 2017, MNRAS, 464, 469

Stark, D. P., Ellis, R. S., & Ouchi, M. 2011, ApJL, 728, L2

Tang, M., Stark, D., Chevallard, J., et al. 2021, MNRAS, 503, 4105

Tasitsiomi, A. 2006, ApJ, 645, 792

Taylor, A. J., Cowie, L. L., Barger, A. J., Hu, E. M., & Songaila, A. 2021, ApJ, 914, 79

Tilvi, V., Malhotra, S., Rhoads, J. E., et al. 2020, ApJL, 891, L10

Tilvi, V., Papovich, C., Finkelstein, S. L., et al. 2014, ApJ, 794, 5

Treu, T., Schmidt, K. B., Trenti, M., Bradley, L. D., & Stiavelli, M. 2013, ApJL, 775, L29

Verhamme, A., Dubois, Y., Blaizot, J., et al. 2012, A&A, 546, A111

Verhamme, A., Garel, T., Ventou, E., et al. 2018, MNRAS, 478, L60

Weinberger, L. H., Haehnelt, M. G., & Kulkarni, G. 2019, MNRAS, 485, 1350

Whitler, L. R., Mason, C. A., Ren, K., et al. 2020, MNRAS, 495, 3602

Yang, H., Malhotra, S., Gronke, M., et al. 2016, ApJ, 820, 130

Zheng, Z.-Y., Wang, J., Rhoads, J., et al. 2017, ApJL, 842, L22

Zitrin, A., Labbé, I., Belli, S., et al. 2015, ApJL, 810, L12