Lyman-alpha Emitting Galaxies and Ionized Bubbles

Modelling of Lyman-alpha Emitting Galaxies and Ionized Bubbles at the Era of Reionization

Hidenobu Yajima, Kazuyuki Sugimura, Kenji Hasegawa
Frontier Research Institute for Interdisciplinary Sciences, Tohoku University, Sendai 980-8578, Japan
Astronomical Institute, Tohoku University, Sendai 980-8578, Japan
Graduate School of Science, Nagoya University, Furo-cho, Chikusa-ku, Nagoya, Aichi 464-8602, Japan
E-mail: (HY)
Accepted ?; Received ??; in original form ???

Understanding emitting galaxies (LAEs) can be a key to reveal cosmic reionization and galaxy formation in the early Universe. Based on halo merger trees and radiation transfer calculations, we model redshift evolution of LAEs and their observational properties at . We consider ionized bubbles associated with individual LAEs and IGM transmission of photons. We find that luminosity tightly correlates with halo mass and stellar mass, while the relation with star formation rate has a large dispersion. Comparing our models with the observed luminosity function by Konno et al., we suggest that LAEs at have galactic wind of and Hi column density of . Number density of bright LAEs rapidly decreases as redshift increases, due to both lower star formation rate and smaller Hii bubbles. Our model predicts future wide deep surveys with next generation telescopes, such as JWST, E-ELT and TMT, can detect LAEs at with a number density of for the flux sensitivity of . By combining these surveys with future 21-cm observations, it could be possible to detect both LAEs with and their associated giant Hii bubbles with the size at .

radiative transfer – line: profiles – galaxies: evolution – galaxies: formation – galaxies: high-redshift
pagerange: Modelling of Lyman-alpha Emitting Galaxies and Ionized Bubbles at the Era of ReionizationReferencespubyear: 2008

1 Introduction

One of the major challenges in today’s astronomy is revealing cosmic reionization history with galaxy evolution. Recent CMB observations suggested cosmic reionization occurred at (Komatsu et al. 2011; Planck Collaboration et al. 2014, 2016). Gunn-Peterson tests by the observations of high-redshift QSOs indicated cosmic reionization completed at (e.g., Fan et al. 2006). However, the ionization history of integer-galactic medium (IGM) has not been understood yet. Recent observations of high-redshift galaxies at are gradually unveiling the cosmic star formation history (Ouchi et al. 2010; Bouwens et al. 2012, 2015; Oesch et al. 2015, 2016), and have allowed us to speculate the cosmic reionization history (Robertson et al. 2015). Yet, even considering the above observational constraints, various reionization histories remain viable (Cen 2003; Yajima & Khochfar 2015). One of the main uncertainties is low-mass galaxy formation with the halo mass less than . Due to their higher number density and the higher escape fraction of ionizing photons from them (Razoumov & Sommer-Larsen 2010; Paardekooper et al. 2013; Yajima et al. 2011, 2014), low-mass galaxies can be responsible for main ionizing sources. However, the detection sensitivities of current observations are not sufficient to constrain the formation of low-mass galaxies. Therefore it is important to theoretically investigate the formation of low-mass galaxies and their contribution to cosmic reionization.

Most of high-redshift low-mass galaxies are observed as emitters (LAEs; Hu & McMahon 1996; Steidel et al. 2000; Iye et al. 2006; Gronwall et al. 2007; Ouchi et al. 2008, 2010; Bond et al. 2011; Ciardullo et al. 2012; Yamada et al. 2012). Ouchi et al. (2010) indicated that LAEs were hosted in halos with the halo mass by the clustering analysis (see also, Gawiser et al. 2007). Verhamme et al. (2008) suggested that LAEs were not dust-enriched well yet so that photons escaped from galaxies against dust attenuation (see also, Yajima et al. 2014). These suggest LAEs are likely to be in the early phase of galaxy evolution (e.g., Mori & Umemura 2006). In addition, LAEs are one of major populations of galaxies contributing to the cosmic star-formation rate density in the early Universe (Ciardullo et al. 2012).

LAEs also have been used as a tool to investigate the early Universe (Iye et al. 2006; Vanzella et al. 2011; Ono et al. 2012; Shibuya et al. 2012; Finkelstein et al. 2013; Zitrin et al. 2015; Song et al. 2016b). Number density of LAEs rapidly decreases at (Ono et al. 2012; Konno et al. 2014), indicating that fluxes from some galaxies were significantly attenuated due to neutral hydrogen in the IGM (e.g., Kashikawa et al. 2006). Meanwhile, LAEs themselves could provide sufficient ionizing photons into the IGM (Yajima et al. 2009, 2014), and some of them could make giant Hii bubbles that allowed photons to reach us. In this case, galaxies can be observed as LAEs. In practice, the most distant LAE has been observed even at (Zitrin et al. 2015). Thus, LAEs can be the key objects in understanding the galaxy formation and cosmic reionization.

In this work, we investigate the evolution of LAEs with their associated Hii bubbles. By modeling both LAEs and Hii bubbles simultaneously, we estimate luminosity functions (LFs), number density of LAEs, size distribution of Hii bubbles, and the relation between flux and the size of Hii bubble. These estimations are useful for future observational missions. James Webb Space Telescope (JWST) is going to be launched at 2018, and aims to observe galaxies at . Its high sensitivity of spectroscopy will make it possible to detect flux from galaxies if they distribute in giant Hii bubbles. Later on, 30-m class ground telescopes, the European Extremely Large Telescope (E-ELT), the Thirty Meter Telescope (TMT), and the Giant Magellan Telescope (GMT), will investigate galaxies at statistically. In addition, several 21-cm observational missions are on going, e.g., the LOw Frequency ARray (LOFAR; Harker et al. 2010), and Murchison Widefield Array (MWA; Lonsdale et al. 2009). In future, Square Kilometer Array phase-2 (SKA-2; Dewdney et al. 2009) will perform large-scale 21-cm tomography. Since 21-cm emission comes from Hi gas, detections of holes of 21-cm signal indicate giant Hii bubbles around galaxies. Therefore, a combination of 21-cm and galaxy observations will provide fruitful information about the cosmic reionization and galaxy formation (e.g., Lidz et al. 2009).

Theoretically large-scale simulations showed that galaxies made patchy ionization structures in an inside-out fashion (e.g., Iliev et al. 2006, 2012; Mellema et al. 2006; Trac & Cen 2007; Ocvirk et al. 2016), whereas low-density void regions would be ionized first, in the so-called outside-in fashion, if AGNs were the main ionizing sources. However, AGNs are unlikely to be main ionizing sources because of the rapid decrease of observed number density at redshift (e.g., Richards et al. 2006), although the contribution of faint AGNs is still under the debate (Madau & Haardt 2015; Yoshiura et al. 2016; Khaire et al. 2016). Thus, in this paper, we model LAEs assuming the reionization model with the inside-out fashion. In order to understand cosmic reionization theoretically, we need high spatial resolution to follow the detailed physical processes of interstellar medium, as well as large volume to consider large-scale inhomogeneity of ionization structure. Even current numerical simulations cannot resolve such a wide dynamic range. For that reason, in this work, we semi-analytically investigate the cosmic reionization and galaxy formation based on simple structure formation models and radiative transfer calculations.

The distributions of LAEs at can provide us valuable information about cosmic reionization and galaxy formation (e.g., see the review by Dijkstra 2014). Furlanetto et al. (2004) investigated the IGM damping wing absorption of flux from star-forming galaxies, and suggested that the line could be a tool to investigate neutral fraction of IGM. McQuinn et al. (2007a) calculated the large-scale inhomogeneous ionization structure of IGM by combining -body simulations with post-processing radiation transfer calculations. They calculated the spatial distribution of LAEs after IGM transmission using simple line profiles of gaussian or single frequency with velocity offset, and showed that the neutral IGM significantly changed the distribution of LAEs (see also, Mesinger et al. 2015). Numerical modeling of large-scale IGM ionization structure always suffers from the expensive calculation cost of radiation transfer. Therefore, Mesinger & Furlanetto (2008) developed a semi-numerical method to calculate the ionization structure of IGM, and investigated reionization structures with the volume of in comoving unit. This semi-numerical method significantly reduces the calculation amount, but simulates the ionization structure accurately (Zahn et al. 2011).

These previous works investigated the impacts of neutral IGM on the distributions of LAEs with detailed calculations of the IGM ionization structure. In this work, we focus on the following two points: (1) What are physical properties of observable LAEs at ?, (2) Can LAEs be observed even at ? We investigate these points statistically based on a large galaxy sample and line profile models with a wide range of galactic inflow/outflow velocities, which have not been understood well in previous studies. The IGM transmission sensitively depends on the line profiles that are related with physical properties of gas in galaxies, e.g., velocity field, Hi column density (Dijkstra et al. 2006; Verhamme et al. 2006; Laursen et al. 2009; Yajima et al. 2012b). Some previous works showed the impacts of neutral IGM on the line profiles. Using the semi-numerical method, Dijkstra et al. (2011) calculated the IGM transmission to flux from a galaxy with at , and suggested such a galaxy could be observed even when IGM was not fully ionized. They assumed that the galaxy had an expanding gas shell with the outflow velocity of . Mesinger et al. (2015) investigated the redshift evolution of LAEs fraction in Lyman break galaxies (LBGs) by combining the large-scale ionization structures calculated by the semi-analytical models and simple Gaussian line profiles with the velocity offset and . Choudhury et al. (2015) also used a Gaussian profile with a velocity offset depending on the redshift, , and estimated the equivalent width distributions to reproduce the observations. Thus most previous works have used simplified models for the line profiles, while the ionization structures have calculated in detail by numerical simulations. Meanwhile, Jensen et al. (2013) estimated outflow velocities from hydrodynamics simulations, and claimed the cosmic reionization rapidly completed at redshifts between and . However, it is generally difficult to accurately model galactic outflow even with state-of-the-art hydrodynamics simulations because the limited numerical resolution does not allow us to determine the strength of supernova feedback. Therefore, as a complementary work to previous works, here we investigate the physical properties of LAEs, in particular, the typical outflow/inflow velocity by using wide-range velocity parameter sets from to . Moreover, in simulating intrinsic line profiles, we utilize two spherical symmetric gas models: (1) expanding uniform clouds, (2) expanding shells.

In addition, we make relations among luminosities, sizes of Hii bubbles, star formation rates, and stellar masses. Recent observations have detected some distant galaxy candidates at by the LBG selection (e.g., Oesch et al. 2015, 2016). They showed the stellar mass and star formation rates of the candidates by SED fittings. Therefore, the relations between luminosity and the physical properties can be a powerful tool to find out the objects which are likely to emit the flux for future deep spectroscopy observations.

Then, using modeled LAEs, we also study observability of LAEs at by next generation telescope, e.g., JWST. The detection of emission lines from distant galaxies is quite important to determine the redshift and confirm the authenticity of candidates. Recent observations have allowed us to detect a distant LAE up to (Zitrin et al. 2015). However, as shown in recent LAEs observations (e.g., Ono et al. 2012), number density of LAEs rapidly decreases toward higher redshifts at . This implies the difficulty of observations of LAEs at due to the IGM attenuation when the neutral fraction of the Universe is high. Therefore, it is important to study whether we can observe LAEs even at by JWST and other upcoming telescopes. We will show number densities of observable LAEs at based on simple expanding gas shell/cloud models.

In this work, we model star formation histories of galaxies based on halo merger histories, and ionized bubbles as expanding isolated ones associated with individual galaxies. Therefore, we do not consider expansion of ionized bubbles due to the overlap of them, i.e., we focus on field galaxies, not clustering galaxies that are likely to make the overlapped giant Hii bubbles. McQuinn et al. (2007b) presented cosmological simulations of different reionization scenarios with post-processing radiation transfer simulations. Depending on physical parameters, i.e., star formation efficiency and escape fraction of ionizing photons, various ionization structure and history can be considered under the same value of volume-weighted mean ionization fraction. They showed the overlaps of Hii regions were delayed in the model that assumed low-mass galaxies contributed to cosmic reionization significantly. In addition, recent CMB observation obtained the Thomson scattering optical depth of IGM , corresponding to the reionization redshift of under the assumption that all IGM were instantaneously ionized at . Since the Gunn-Peterson tests using QSOs indicated that the reionization was completed at (e.g., Fan et al. 2006), Hii bubbles might not be overlapped well at . Thus, in this work, we simply follow evolutions of isolated expanding Hii bubbles without the effect of overlap. With a simple prescription, we also discuss the effect of the overlap of Hii bubbles on the IGM transmission.

Tilvi et al. (2011) modeled LAEs at based on halo merger histories from cosmological -body simulations. They assumed that luminosity of each galaxy was simply proportional to SFR and all photons escaped from galaxies, and all ionizing photons were absorbed by interstellar gas. On the other hand, in this work, we focus on LAEs at , and consider non-zero escape fraction of ionizing photons to reproduce the observed Thomson scattering optical depth of IGM, and IGM transmission depending on line profiles.

We use the cosmological parameters, , , and (Komatsu et al. 2011; Planck Collaboration et al. 2016).

The paper is organized as follows. We describe our models of star formation and Hii bubbles in §2. In §3, we present the results, which include the properties, LFs of LAEs and redshift evolution of number density of LAEs. We discuss the relation between the flux of galaxies and the sizes of Hii bubbles in §4, and summarize in §5.

2 Model

In the current standard picture of structure formation, halos grow via minor/major mergers. In our model, the star formation rate (SFR) in a halo is assumed to be proportional to the growth rate of the halo as follows:


The relation between stellar and halo mass can be written as . The abundance matching analysis by Behroozi et al. (2013) indicated that for . In addition, Behroozi & Silk (2015) showed the weak dependence of on halo mass and redshift at the mass range , although increases with at (see also, Moster et al. 2013; Kravtsov et al. 2014). Therefore we here assume . Thus, we estimate SFR from the growth rate of halo mass with a constant tuning parameter , i.e., .

In order to estimate the growth of halos, we use halo merger trees based on an extended Press-Schechter formalism (Somerville & Kolatt 1999; Khochfar & Burkert 2001, 2006). The halo merger trees include 5000 realizations with the halo mass range of at . In this work, we allow star formation only for halos with the mass , because star formation in less massive halos can be significantly suppressed due to UV background or internal stellar feedback (e.g., Okamoto et al. 2008; Hasegawa & Semelin 2013). In deriving statistical properties, e.g., cosmic SFR density, stellar mass density and LF, we sum the contribution of each merger tree with normalization factors that reproduce the halo mass function of Sheth & Tormen (2002). We determine the tuning parameter by using observed cosmic SFR density and stellar mass density at . Figure 1 shows our modeled SFR and stellar mass density with observations. We choose the parameter by the least square fitting to the four points of the observations. The observed SFR and stellar mass densities consider only galaxies with and , respectively. Therefore, we consider only the galaxies satisfying the above criteria for the fitting. The best fit value is . Note that, so far there is no available data about at and can change with halo mass and redshift. However, for simplicity, we assume is constant. Following the above way, we derive star formation history from each halo merger tree. Although we consider only bright galaxies in deriving the parameter via the comparisons with the observations, we follow star formation of all halos with in studying cosmic reionization and LAEs in this work. Black dash lines in the panels of Figure 1 represent the total SFR and stellar mass densities as a function of redshift. The top panel of Figure 1 shows that low-mass faint galaxies have non-negligible contribution to the cosmic SFR density and are expected to contribute to cosmic reionization as recent studies indicated (e.g., Yajima et al. 2009, 2011, 2014; Robertson et al. 2013; Cai et al. 2014; Wise et al. 2014).

We also derive UV LFs by converting SFR to UV flux with (Madau et al. 1999). Figure 2 shows our modeled LFs at and with those from the recent observation by Bouwens et al. (2015). The observation indicated that a part of UV flux are absorbed by dust with the escape fraction . Our modeled LFs match the observation well with the same UV escape fraction.

Galaxies ionize the IGM as star formation proceeds. In a one-zone approximation, we estimate time evolution of cosmic ionization degree (Barkana & Loeb 2001),


where is the volume fraction of Hii, is the present-day hydrogen number density (), is the intrinsic ionizing photon emissivity per unit volume, is the case-B recombination rate, is a clumpiness factor of IGM, and is escape fraction of ionizing photons, which is a free parameter here. We estimate from star formation history of each halo by using a population synthesis code starburst99 (Leitherer et al. 1999) with the Salpeter initial mass function with the metallicity of . The reionization history is shown in the upper panel of Figure 3. As increases, the IGM is ionized earlier. We assume the recombination rate for () and , as suggested by numerical simulations (Pawlik et al. 2009; Jeon et al. 2014).

The cosmic reionization history is regulated by . Free electrons produced by the cosmic reionization contribute to the Thomson scattering optical depth () of CMB photons, defined as


where is the redshift at the time of recombination. In this work, we assume the single ionization fraction of helium is the same as the one for hydrogen at , and the double ionization takes place at (e.g., Wyithe et al. 2010; Inoue et al. 2013). Recent simulations show that indeed the fraction of Heii is close to the Hii fraction at high redshifts, although the ionization fraction of helium is slightly lower than the one for hydrogen (Ciardi et al. 2012). The top panel of Figure 3 represents the ionization history of hydrogen gas with the different . The bottom panel of Figure 3 shows with different . The Thomson scattering optical depth increases with redshift in a way depending on due to the different ionization histories. We find nicely reproduces the CMB observation (Planck 2016). Therefore, in this work, we adopt with no redshift evolution. In fact, Yajima et al. (2014) showed is constant with redshift and by cosmological simulations with radiative transfer calculations.

Figure 1: Upper panel: Star formation rate density. Red line represents our modeled star formation history based on the halo merger trees. Triangle symbols show the observation by Bouwens et al. (2015). The observed star formation rate densities are estimated by integrating the luminosity functions in the range of . Our model (red line) also consider only the galaxies brighter than the same limiting magnitude. Lower panel: Stellar mass density. Red line is the cumulated stellar mass of our model considering only the galaxies with . Triangle symbols show the observation by Song et al. (2016a), who integrated the derived stellar mass functions for galaxies with . Black dash lines represent the star formation rate and stellar mass densities considering all galaxies without the thresholds equivalent to the observations.
Figure 2: UV luminosity functions. Open squares and circles represent observed LBGs at and by Bouwens et al. (2015).
Figure 3: Upper panel: Ionization history. Different lines show ionization histories with different escape fractions of ionizing photons. Lower panel: Thomson scattering optical depth. Different lines are estimated with different ionization histories presented in the upper panel. Yellow shade represents the estimation by the CMB observation (Planck Collaboration et al. 2016).

3 Results

3.1 Evolution of ionized bubbles around galaxies

We follow the growth history of Hii bubbles around individual galaxies with the star formation efficiency and estimated above. Sizes of Hii bubbles evolve with ionizing photon emissions, recombination, and cosmic expansion as follows (Cen & Haiman 2000):


where is the Hubble constant at specific redshifts, and is intrinsic ionizng photon emissivity of each galaxy. The ionizing fronts can propagate up to the Strmgren radius during the recombination time scale (Spitzer 1978). The recombination time scale is , which is longer than the typical time scale over which SFR changes more than factor . Therefore, Hii bubbles do not reach the equilibrium state, and we need to consider SFR history to estimate the sizes of Hii bubbles at given redshifts. We estimate probability distribution function (PDF) of the sizes of ionized bubbles () as shown in Figure 4. In our model, higher mass halos tend to possess larger ionized bubbles. Due to the decrease in the number density of halos on the high-mass end of a halo mass function (Sheth & Tormen 2002), the PDF of rapidly decreases at larger . Although the ionizing front does not reach the size of Strmgren sphere , it can be used as a rough indicator of sizes of Hii bubbles. As redshift decreases, the IGM density decreases while the number density of massive halos with higher ionizing photon emissivity increases. Therefore, the tail of PDF in the large- end shifts to larger at lower redshift. Future 21 cm observations, e.g., SKA-2, is supposed to probe the IGM ionization structure with the angular resolution of . Therefore, at , the tail of PDF at large can be observationally investigated in future. The halo number density monotonically increases as the halo mass decreases at a fixed redshift. Since the size of Hii bubble is positively related with halo mass as will be shown in Section 3.3, it seems that the PDFs in the small- end monotonically increase as decreases. However, there are peaks in the PDFs, below which they decrease as decreases. This is caused by the threshold of halo mass for star formation imposed in this work.

Note that, in this work, we do not take the effect of overlap of Hii bubbles into account. The galaxy clustering causes the overlap of Hii bubbles that can extend the tails of PDFs to larger size. Since the overlapping effect can be larger as redshift decreases, the slopes of PDFs at lower redshifts can be changed to be more shallower. We will discuss the impacts of the overlap of Hii bubbles on luminosity functions in Section 3.2.

Figure 4: Probability distribution functions of sizes of Hii bubbles associated with individual galaxies. Vertical dash lines are corresponding to the viewing angle of 1 arcmin.

3.2 luminosity functions

Absorption of ionizing photons by interstellar medium within galaxies results in emissions via recombination processes, while escaped photons cause the cosmic reionization as shown in the previous section. The luminosity () is estimated by


where is the escape fraction of photons from galaxies, is the energy of a photon. The escape fraction of photons can be lower than because the path length of photons until escape can be longer due to multiple scattering process. However, if dust mainly distributes in Hi gas clumps, does not become lower than because photons are scattered by hydrogen on surface of the clumps before interacting with dust. Ciardullo et al. (2012) indicated that for observed LAEs at . In addition, cosmological simulations of Yajima et al. (2014) showed that was and similar to at . Therefore, in this work, we assume that .

Next, we estimate IGM transmission as a function of wavelength. As in Cen & Haiman (2000), we divide the paths along which photons travel from galaxies to us into the two parts, i.e., outside and inside ionized bubbles, and separately estimate each contribution. The transmission outside ionized bubbles is estimated as follows:


where . Here, is redshift when the cosmic reionization completes, which we set , is redshift when photons pass through ionizing front, is redshift of galaxy, is the scattering cross section for Hi gas, and is the size of the cosmological horizon at . We assume the outside of the bubble is completely neutral, i.e., . The is estimated by (Verhamme et al. 2006)


We set in this work. Here, is the Voigt function,


where , is the line-center frequency, is the Doppler width, , is the natural line width. We here use the fitting formula of given by Tasitsiomi (2006).

Even inside ionized bubbles, a tiny fraction of neutral hydrogens exist. We estimate the IGM transmission from ionizing front to virial radius with the neutral fraction under the ionization equilibrium state,


Note that, the optical depth outside the ionized bubble is dominant in our work. The IGM transmission is mostly almost zero at , where is the wavelength of line center.

Considering the IGM transmission, we derive LFs, and compare them with the observation of Konno et al. (2014). Depending on the shape of line profile, the IGM transmission significantly changes. Even with recent deep spectroscopies, however, it is difficult to determine intrinsic line profiles, i.e., before the IGM extinction (e.g., Ouchi et al. 2010). The intrinsic line profile depends on the physical nature of galaxies, e.g., Hi column density and velocity field (e.g., Verhamme et al. 2006). In this work, we calculate intrinsic line profiles by radiation transfer simulations using the code developed in Yajima et al. (2012b), and study the physical nature of LAEs through the comparison with the observation (Konno et al. 2014). For the purpose, we employ the following two models for the internal velocity structure of a galaxy: (1) expanding cloud, (2) expanding shell. Both cloud and shell models have only two parameters: (1) Hi column density , (2) outflow velocity . The intrinsic line profiles of the cloud model are calculated based on spherically outflowing gas with the following velocity structure,


where is the edge of the spherical cloud and is the outflow velocity at . Here uniform gas density in the cloud is assumed. Although the simulated profiles do not depend on , the source size is likely to increase with . Since the surface brightness decreases with the source size, the detectability of flux depends on the choice of . Nevertheless, in this work, we assume that LAEs are compact and no flux is lost based on theoretical motivation: the typical sizes of galaxies become (e.g., Mo et al. 1998) where is a halo spin parameter and expected to be less than (e.g., Bullock et al. 2001). Therefore we consider in which the source size may not be subject to the flux loss. The uncertainty of the surface brightness will be discussed in Sec. 4.4.

In the shell model, the velocity field is monochromatic, i.e.,


Our radiation transfer calculations are carried out in 100 spherical shells. All shells have same Hi density in the cloud model, while only most outside shell has Hi gas in the shell model.

The cloud and shell models approximate galactic outflows due to stellar feedback. When most of gas are being evacuated due to starburst, Hi gas structure may be close to the shell model. On the other hand, in the case of smooth star formation history, stars keep being formed in static Hi gas at near galactic centers, while a part of gas are evacuated from galaxies. In this case, the gas structure may be closer to the cloud model. The star formation and gas structure in high-redshift galaxies are still under the debate (e.g., Kimm & Cen 2014; Hopkins et al. 2014; Yajima et al. 2017). Therefore we here study LAEs by using these two models.

Figure 5 shows the modeled line profiles with various outflowing velocities. The line profiles of the cloud model get the asymmetric shape with a peak at redder wavelength for the outflow velocity field, because photons at bluer wavelength are scattered by Hi gas due to the Doppler shift. As increases, the line profiles are extended, and the peak frequency shift farther from the line center frequency. In the case of the shell model, the line profiles become more complicated because of the back scattering effect (Verhamme et al. 2006). Observable photons scattered by the shell at the far side from an observer make a bump at redder wavelength. In addition, when the Hi column density is low and the expanding velocity is large, most of photons can directly escape from the shell, resulting in the peaks at the line center frequency as shown in middle and right panels of the shell model. For the inflow velocity field, the line profile becomes the mirror symmetric shape to the one for outflow with the same absolute value of the velocity.

Figure 5: line profiles from expanding spherical gas clouds or shells. Upper panels represent the line profiles of the cloud model that assume the velocity field as , where is the radius of clouds and is the velocity at . Lower panels present the line profiles of the shell model in which the velocity is constant, i.e., . The different panels from left to right show the profiles considering different Hi column densities: (panel A), (panel B), and (panel C). The different lines represent the different expanding velocities at . The is normalized to be unity when it is integrated over the wavelength. The of panel (B) and (C) is artificially reduced by a factor 2 and 10.

With the simulated intrinsic line profiles, we can estimate LFs for given and . In this work, we infer typical and of LAEs by comparing modeled LFs with the observed one by Konno et al. (2014). Figure 6 show the LFs. Here we calculate the LFs for three Hi column density models, (model A), (model B), and (model C). The column density of model A is corresponding to Damped Lyman- Systems (DLAs: Wolfe et al. 2005). Yajima et al. (2012a) showed that DLAs distributed at lines of sight passing star-forming regions in high-redshift star-forming galaxies by combining cosmological SPH simulations with radiation transfer calculations. These column densities are optically thick to ionizing photons, hence not consistent with . However, recent simulations showed ionizing photons mostly escape along ionized holes created by radiative and SNe feedback (e.g., Yajima et al. 2009, 2011; Kimm & Cen 2014). Thus, of 20 can be considered as the fraction of viewing angle along which star forming regions are not covered by Hi gas. For simplicity, we do not take account of the effect of such holes on line profiles. Note that, however, line profiles somewhat change due to the holes, clumpiness or other detailed structure in Hi gas (Dijkstra & Kramer 2012).

The shaded regions in Figure 6 represent the LFs using different line profiles with the velocity range from to . The best fit velocities for the cloud model are (model A), (model B), and (model C), while in the case of the shell model they are (model A), (model B), and (model C), as summarized in Table 1. Only model A with the velocity range can reproduce the observed LF well. Therefore we suggest LAEs are likely to have high Hi column density with and outflowing Hi gas with velocity . These outflow velocities are consistent with recent observations of [Cii] emissions from LAEs at in which they showed the velocity offsets between and [Cii] lines (Pentericci et al. 2016). Note that, the cloud model has low velocity component at a galactic center, unlike the shell model. Therefore, the velocity can not be compared directly between the models. The profile of the cloud model monotonically shifts to redder one as the outflow velocity increases in model A as shown in Figure 5. On the other hand, as the Hi column density decreases, the profile moves back to the line center frequency at specific velocity of . This is because photons can escape from the cloud before shifting to longer wavelength due to lower optical depth. In this work, the best fit velocities for model B and C roughly correspond to those producing the profiles shifted farthest away.

In the case of the shell model, even for model B and C, the wavelength of the bump made by the back scattering effect shifts to longer one as the outflow velocity increases. However, the flux at the line center frequency becomes high as the outflow velocity increases. Therefore, for lower column density models, flux is efficiently attenuated when the outflow velocity is higher than specific values.

In addition, the width of line profile becomes smaller as the Hi column density decrease. Even with the best fit velocities, the number densities of LAEs in the model B and C are smaller than the observation because of the narrower line profiles resulting in lower IGM transmission. Hence, in order to get higher IGM transmission to reproduce the LF, LAEs are likely to have the column density higher than . Thus, in this work, we consider model A as our fiducial model. Verhamme et al. (2008) also suggested similar column densities and outflow velocities for LAEs at via comparisons of their modeled line profiles based on the shell model with the observations (but see, Konno et al. 2016). On the other hand, local LAE analogues are likely to have lower Hi column densities (Henry et al. 2015; Verhamme et al. 2015).

The IGM transmissions are presented in Figure 7. Even at more than half of intrinsic flux can be attenuated due to residual neutral IGM outside of Hii bubbles. This is roughly consistent with observations that suggested LAEs fraction in LBGs rapidly decreases from to (e.g., Ono et al. 2012; Konno et al. 2014). Bolton & Haehnelt (2013) suggested number densities of Lyman-limit systems or DLAs increased with redshift, resulting in the decreased number density of LAEs due to the IGM attenuation (see also, Choudhury et al. 2015; Mesinger et al. 2015). In addition, recently Sadoun et al. (2017) suggested that the decreased number density of LAEs can be explained by considering only infalling IGM near virial radius, without the modeling of redshift evolution of the neutral fraction of whole IGM. They assumed that the IGM ionization structure was determined by external UVB, while the IGM transmission in our models basically results from the neutral IGM outside Hii bubbles created by galaxies themselves.

As discussed in Sec. 3.3, the size of Hii bubble increases with stellar (or halo) mass. Therefore, the transmission increases with stellar mass (see Eq.  6). As redshift increases, the typical size of Hii bubbles decreases (Figure 4), and the mean IGM density increases. Thus, the transmission decreases as redshift increases. At , the transmission of even massive galaxies with is less than , while even low-mass galaxies have higher transmission than at . In addition, at higher redshifts, the transmission of the shell model becomes higher than the cloud model. As shown in Figure 5, the profiles of the shell model have a bump at longer wavelength due to the back scattering effect. photons at this bump can penetrate the IGM even when a Hii bubble is small at high redshifts.

Emergent line profiles are shown in Figure 8. As Hi column density decreases, intrinsic line profiles become narrower and peak positions shift to shorter wavelength. IGM transmission increases with wavelength because photons originally with long wavelength are redshifted and cease to be scattered by the IGM, while flux near the line center is reduced efficiently by the IGM scattering. FWHMs of the emergent line profiles of the cloud model are (model A), (model B) and (model C). The emergent line profiles of the shell model have somewhat larger FWHMs. Therefore it is difficult to distinguish the different column density models in the current spectroscopic observation with the resolution of (e.g., Shibuya et al. 2012). Future high-dispersion spectroscopies with , e.g., Prime Focus Spectrograph on Subaru or JWST, will be able to reveal the detailed shape of profile.

Line profiles in inflowing gas models result in the underproduction of observable LAEs, because most of photons are scattered by the IGM. This is consistent with the observation by Ouchi et al. (2010), which indicated LAEs at are likely to have outflowing gas by the composite spectrum. In addition, Shibuya et al. (2014) measured outflow velocities of individual LAEs at , and indicated that LAEs were likely to have outflow with .

Next we estimate the redshift evolution of LF based on the model A with (cloud model) and (shell model). Figure 9 shows the modeled LFs at and . Note that, here we use the same line profile for all halos and redshift. At higher redshifts, typical SFR is smaller due to lower halo mass and halo growth rate. In addition, as redshift increases, typical size of Hii bubbles decreases, resulting in the lower IGM transmission. As a result, the LF rapidly shifts to the fainter side at higher redshifts.

The LFs of the shell models at faint-end are somewhat larger than the cloud model. As explained above, the bump at longer wavelength due to the back scattering effect causes the higher IGM transmission, resulting in the higher LFs. We also compare our modeled LF at with the observed LF at by Ouchi et al. (2008). At this redshift, the cosmic reionization is thought to be completed, therefore the IGM transmission is likely to be . Black solid line represents the LF at without the IGM attenuation. The modeled LF nicely matches the observation at , but somewhat larger at . Bright LAEs are likely to reside in massive halos. As halos grow, interstellar gas could be dust enriched via type-II supernovae, resulting in the decrease of (e.g., Yajima et al. 2014). The lower might explain the lower LF at the bright end. In addition, Mesinger (2010) suggested that the IGM was not completely ionized even at . The residual neutral IGM could decrease number of observed LAEs.

Figure 6: luminosity functions (LFs) at . Square symbols represent the observed LF of LAEs at by Konno et al. (2014). Upper and lower panels present modeled LFs using line profiles of the cloud and shell models, respectively. Different panels from left to right show modeled LFs based on line profiles to different Hi column densities. Magenta and cyan shades show the range of LFs considering line profile with different outflow and inflow velocity with the range . Black solid lines are best fitted ones to the observation. Black dash line in the left panel shows the LF before considering IGM transmission.
Figure 7: IGM transmission as a function of stellar mass of galaxies. Filled and open circles represent median values of the cloud and shell models, respectively. Different colors indicate different redshifts.
Figure 8: line profiles. Red solid lines show emergent line profiles of a halo of at z=7.3, which evolves to a halo of at . Intrinsic luminosity of the halo is . Left and right panels represent the line profiles of the cloud and shell models. Emergent luminosities after IGM transmission are shown in the panels. Black dash lines are the line profiles before considering IGM transmission.
Figure 9: luminosity functions at and . Open triangles and squares show observed LAEs at (Ouchi et al. 2008) and (Konno et al. 2014). Solid and dash lines represent the LFs based on the cloud and shell models. Solid line at shows the LF without the IGM attenuation.

3.3 Relation between Ly luminosity and stellar mass

Figure 10 shows the sizes of Hii bubbles and luminosities considering the IGM transmission. properties are calculated by using the profile of model A of the cloud model. The shades represent the range of in the sample. We see that tightly correlate with stellar mass, as , while the relation with SFR shows a large dispersion.

SFR rapidly increase by major merger. However, is not so sensitive to the short-time fluctuation of SFR because of the longer time-scale for reaching the ionization equilibrium state. As a result, the relation between and SFR shows the large dispersion. High-redshift galaxies have been observed as so-called Lyman-break galaxies (LBGs: Bouwens et al. 2012), via the Lyman-break technique so far. Our results indicate LBGs at with similar UV brightness can have different fluxes due to the scatter of IGM transmission. At , our sample is limited by the stellar mass of galaxies . In our model, since the stellar mass is simply proportional to halo mass as , this means that there is no progenitors with at in our halo sample, which is constructed to have the mass range from to at .

In contrast to the SFR- relation, both luminosity and tightly correlate with . The - relation does not change with redshift significantly, while the - becomes smaller as redshift increases. This is because Hubble constant (i.e., expanding velocity of IGM) becomes large at higher redshift. Therefore, although decreases as redshift increases due to higher IGM density, the IGM transmission does not decreases significantly. Thus, luminosity does not depend sensitively on redshift. Note that, however, we have not taken the overlaps of Hii bubbles into account. When galaxies are clustered, giant Hii bubbles can form by the overlaps of individual Hii bubbles. This can cause scatters in the relations between and or .

The detection sensitivity of recent observations of LAEs at was corresponding to luminosity of (Ono et al. 2012; Shibuya et al. 2012; Finkelstein et al. 2013; Vanzella et al. 2011; Konno et al. 2014; Zitrin et al. 2015). In our model, median and minimum stellar masses producing at are and , respectively. Therefore, by considering the relation , we suggest that the observed LAEs at should be hosted in halos with , and the median halo mass is .

Figure 10: Upper panel: Sizes of Hii bubbles as a function of stellar mass and SFR. Green, red and blue lines represent median values to galaxies at and , respectively. The shades show quartiles at each bin. Lower panel: luminosity using the cloud model as a function of stellar mass and SFR.

3.4 Redshift evolution of number density of observable LAEs

So far the line has been used as the most strong tool to confirm the redshift of distant galaxy candidates (e.g., Iye et al. 2006; Finkelstein et al. 2013; Zitrin et al. 2015). However, it is widely thought that LAEs at are difficult to be detected because of the IGM opacity. Here, we estimate the number density of LAEs () with higher flux than specific thresholds. Figure 11 shows the number density of LAEs with and . The detection limits of current observations with a reasonable integration time are (e.g., Shibuya et al. 2012). As explained in Sec. 3.2, the number density of bright LAEs monotonically decreases with increasing redshift. Given the detection limits of , wide field surveys of are able to detect LAEs up to . This is consistent with recent observed LAEs at (Finkelstein et al. 2013; Oesch et al. 2015; Zitrin et al. 2015). The LAEs with at are quite rare, with .

Spectroscopies of next generation telescopes, e.g., JWST, are supposed to achieve the sensitivity of with a reasonable integration time. If galaxies at have outflowing gas with , the number density of LAEs with at is a few . As shown in Figure 10, bright LAEs are hosted in massive halos. The median halo and stellar mass of LAEs with at are and , respectively. It was suggested that the observed LBG at , GNz11, had the stellar mass of (Oesch et al. 2016), which is corresponding to in our model. Hence, it will be challenging to detect flux from GNz11 even by future spectroscopies with the line sensitivity of .

Different line profile models predict different number densities of observable LAEs. The IGM transmission becomes more sensitive to the intrinsic line profile models, since the typical size of Hii bubbles gets smaller with increasing redshift. As a result, the difference of among model A, B and C becomes larger at higher redshift, and more than order unity at . Future observation would also allow us to discriminate intrinsic line profiles, which in turn provide information about Hi column density and outflow velocity, by comparing the observed number density of LAEs with the theoretical models. On the other hand, it is difficult to distinguish the cloud and shell models from the number density alone as shown in the figure.

Even next generation telescopes, e.g., GMT, E-ELT, TMT, will be difficult to have the sensitivity of . However, if future telescopes somehow achieve such a high sensitivity, wide field surveys of would be able to reach LAEs at .

In this work, we consider isolated Hii bubbles associated with individual LAEs in the estimation of IGM transmission. However, at lower redshift, Hii bubbles can be overlapped each other (e.g., Iliev et al. 2012; Hasegawa & Semelin 2013; Ocvirk et al. 2016). The overlapped Hii bubbles can enhance the IGM transmission. This effect will be investigated in Hasegawa et al. (in preparation) by combining large-scale -body with small scale radiative-hydrodynamics simulations.

Figure 11: Number density of LAEs as a function of redshift. Red, blue and green lines show the number densities of LAEs with obtained assuming line profiles of model A, B and C, respectively. Colored solid and dash lines represent the number density using the line profiles of the cloud and shell models. Black dash and dot lines are the number density of LAEs with and obtained assuming the line profile of model A of the cloud model.
Model (cloud) (shell)
A 180 130
B 190 60
C 110 40

NOTES.—For each Hi column density, the outflow velocity is chosen to reproduce the observed luminosity function of LAE at (Konno et al. 2014). (cloud) and (shell) represent best-fitted outflow velocities in the expanding cloud and shell models. Escape fractions of ionizing, UV and photons are same for all models.

Table 1: Model parameters

4 Discussion

4.1 Impact of overlap of Hii bubbles on radiation properties

So far, we consider only isolated Hii bubbles around LAEs. However, due to the clustering of galaxies, Hii bubbles overlap each other (e.g., McQuinn et al. 2007b; Zahn et al. 2011; Iliev et al. 2012). This leads to the expansion of the size of Hii bubbles, resulting in the increase of IGM transmission. Recently, Castellano et al. (2016) observed two LAEs in a galaxy overdensity region at that might imply the clustering galaxies made a giant Hii bubble. Here, considering the overlap effect, we artificially expand sizes of all Hii bubbles associated with galaxies as follows:


where is the boost factor to expand the original size of Hii bubbles. For example, the condition of roughly considers that galaxies with similar luminosities distribute in an overlapped Hii bubble, since the size of the Hii bubbles is proportional to the power of one-third of luminosities. In practice, can depend on original sizes of Hii bubbles. Small Hii bubbles associated with low-mass galaxies near massive ones can be merged into giant Hii bubbles, resulting in high . On the other hand, the overlap between giant Hii bubbles is not frequent, leading to low .

As the size of Hii bubbles increases by the overlap effect, even photons with the line center frequency can penetrate IGM. This is likely to affect the outflow velocity to reproduce the observation. Figure 12 represents the best-fit outflow-velocities reproducing the observed LF at as a function of . Due to the expansion of Hii bubbles, the modeled LFs shift to brighter side, i.e., right side along x-axis in Figure 6. Therefore, as increases, the best-fit outflow-velocities monotonically decrease. All models show that the velocities become smaller than at . However, in the range of , the velocities do not go below , i.e., the LFs with do not exceed the observed LF. Therefore, if typical value of is smaller than 4, LAEs are likely to have outflowing gas.

The two-point correlation functions of LAEs derived from recent LAE surveys at showed the typical correlation length was in comoving scale and corresponding halo mass was (Ouchi et al. 2010, 2017). The correlation length of LAEs does not change with redshift significantly (Ouchi et al. 2010). In our model, the Hii bubble size around halos of is . The fact that the Hii bubble size is typically smaller than might indicate that does not largely exceed unity. The detailed overlap effect should be investigated by future cosmological simulations.

The velocities of model B and C of the cloud model does not change largely even if increases. As seen in Figure 5, the shell model with a large outflow velocity produces the strong peaks at the line center frequency. Therefore, even if is small, the velocities of the shell model do not increase largely.

Figure 12: Outflow velocities matching the observed luminosity function as a function of the boost factor of Hii bubble size. The IGM transmission is calculated to Hii bubble extended artificially, . Solid and dash lines represent the expanding cloud and shell models, respectively. Different colors indicate different Hi column densities in the calculations of line profiles.

4.2 Condition for galactic outflow from LAEs

As shown in Section 3.2, high-redshift LAEs are likely to have galactic wind with . Here we roughly derive the condition for making gas outflow with in a spherical gas cloud model. Supernova (SN) feedback can be responsible for causing strong outflow in high-redshift low-mass galaxies (e.g., Kimm & Cen 2014; Kimm et al. 2015). In the assumption that supernovae give feedback to all gas uniformly, we estimate the total kinetic energy of outflowing gas as


Here is gas mass after star formation, is stellar mass, is initial gas mass, is a star formation efficiency with respect to the initial gas mass, is the conversion efficiency from total supernova energy to kinetic one of gas, is the released energy for each supernova (e.g., Hamuy 2003), and is the number of supernova. For Salpeter-like IMF, . By using the star formation efficiency , we simply estimate as


where is the escape velocity of a halo.

Next we obtain as a function of . Recently Kim & Ostriker (2015) showed the final momentum produced by a single SN with the energy of as follows: (see also, Cioffi et al. 1988; Thornton et al. 1998). In this work, we ignore the weak density dependence and assume . Since the final momentum is almost linearly proportional to the injected SN energy (Cioffi et al. 1988; Kimm & Cen 2014), we approximate the total momentum produced by multiple SNe as follows: , where is the total energy of supernovae. Therefore, using , we derive . Note that can not exceed unity according to the energy conservation. For that reason, we set at because the above expression gives . Thus, Equation 14 is written by using as follows:


We find that the condition of is required to cause the galactic outflow with from a halo with at which is minimum halo mass to produce the observable luminosity . We can also convert the tuning parameter in our star formation model to , as . This is roughly similar to the value estimated above. As halo mass increases, higher is required to cause galactic outflow. Here we have estimated assuming all gas has same outflow velocity in the simple spherical cloud model. For example, in the case of disk galaxies, only a part of gas can be evacuated along the normal direction to galactic disk as shown in numerical simulations (e.g., Agertz et al. 2011). In this case, strong outflow can be caused even with smaller because piled gas mass can be lower.

4.3 Relation between Ly luminosity and size of ionized bubble

LAEs can be responsible for ionizing sources of cosmic reionization (e.g., Yajima et al. 2009, 2011, 2014). Figure 13 shows as a function of the size of associated Hii bubble . The luminosity steeply increases with due to higher IGM transmission for large Hii bubbles. The relation between and does not significantly change with redshift. IGM gas density increases with redshift as , resulting in lower IGM transmission at higher redshift. However, galaxies tend to have higher ionizing photon emissivity and intrinsic luminosity for a fixed size of Hii bubbles at higher redshift. The combination of these effects leads to the weak redshift dependence in the - relation.

The yellow shaded region represents the viewing angle of arcmin and the flux of for . This region corresponds to the LAEs with associated Hii bubbles that are detectable both as LAEs and holes in 21-cm signal by future galaxy observations by JWST and 21-cm tomography by SKA-2, respectively. Future 21-cm observations will be able to probe giant Hii bubbles around bright LAEs with at . Note that, the overlap of Hii bubbles is not considered in Figure 13. The effect of the overlap can shift our results to larger side of the Hii bubble.

The differential brightness temperature caused by galaxies shows inner positive and outer negative ring-like structure (e.g., Chen & Miralda-Escudé 2004; Yajima & Li 2014). The detailed structure depends on SED. If galaxies host X-ray sources like AGNs, the positive region is extended due to partial photo-ionization heating. Therefore, future 21-cm observations may also give us information about nature of X-ray and UV sources in bright LAEs.

Figure 13: luminosity as a function of size of Hii bubble. Triangle, circle and square symbols represent median values at and , respectively. Error bars show quartiles. The line profile of model A is used. Horizontal dash line shows flux of for . Vertical dash line represents the viewing angle of 1 arcmin for .

4.4 Limitations in current models

Here we have investigated properties of LAEs and ionized bubbles at the era of reionization based on simple analytical models. Our current model is not able to consider following important aspects.

(1) redshift evolution of parameters: In this work, we have assumed the escape fraction of photons and the shape of line profile did not change with redshift. These quantities depend on physical properties of interstellar medium that can change with redshift (e.g., Mao et al. 2007; Yajima et al. 2014).

(2) overlap of ionized bubbles: As ionized bubbles grow, they overlap each other, leading to formation of giant Hii bubbles. This increases the IGM transmission of photons with the frequency near the line center.

(3) surface brightness: Multiple scattering processes of photons in galaxies can make extended surface brightness distribution. Recent observations indicated star-forming galaxies had extended haloes (e.g., Steidel et al. 2011; Momose et al. 2014, 2016). Therefore, some fraction of flux could be lost in the observations of LAEs. Understanding the surface brightness distribution and the flux loss are quite important in liking between theoretical models and observations. However, the fraction of lost flux sensitively depends on the sensitivities of the observations, and it is uncertain now. In this work, we have assumed that all photons escaped from galaxies can be observed, if they are not scatted by IGM. This assumption can be reasonable when sizes of galaxies are not much larger than angular resolutions of observations. Recently Mas-Ribas et al. (2017) showed JWST will be able to reach the sensitivity for extended emission sources. Here we roughly estimate the surface brightness as . For detections of bright LAEs of , the source size should be smaller than , corresponding to at . This is comparable with the virial radius of halos with masses of . Therefore the flux loss due to the faint surface brightness is unlikely to be significant. In addition, Momose et al. (2014) suggested that typical scale lengths of halos of observed LAEs at were and do not change with redshift significantly. Note that, however, the surface brightness sensitively depends on complicated distributions of gas and velocity fields (e.g., Yajima et al. 2015). In reality, the surface brightness distribution should be considered in the estimation of the flux loss for extended sources, which is beyond the scope in this paper. We will study the detectability of distant LAEs with the calculation of surface brightness by combining cosmological hydrodynamics simulations and radiative transfer calculations in future work.

5 Summary

We present models of emitting galaxies (LAEs) with IGM transmission considered at the era of reionization. Based on halo merger trees and a simple star formation model, we estimate cosmic star formation and cosmic reionization history. Our model uses 5000 realizations of halo merger trees with the halo mass range from to at . As a result, our model reproduces the observed cosmic star formation densities, stellar mass densities, luminosity functions of Lyman-break galaxies with a tuning parameter, . Our modeled star formation history, with an escape fraction of ionizing photons , also provides a cosmic reionization history consistent with the Thomson scattering optical depth indicated by Plank (2016).

Based on the above parameters, we model LAEs and Hii bubbles using individual halo merger trees. Our models show the distribution function of the size of Hii bubble, and indicate giant Hii bubbles associated with bright LAEs at can be probed by future 21-cm observation using SKA. We find that flux is tightly related with stellar/halo mass, while there is a large dispersion in the relation between flux and SFR. By comparing our models with the observed luminosity function of LAEs at by Konno et al. (2014), we indicate that LAEs are likely to have the Hi column density of and the outflowing gas with . Using these parameters, we predict that future wide deep survey can detect LAEs at with and for the flux sensitivity of