Lyman continuum leakage in faint star-forming galaxies at redshift z=3-3.5 probed by gamma-ray bursts

We present the observations of Lyman continuum (LyC) emission in the afterglow spectra of GRB 191004B at $z=3.5055$, together with those of the other two previously known LyC-emitting long gamma-ray bursts (LGRB) (GRB 050908 at $z=3.3467$, and GRB 060607A at $z=3.0749$), to determine their LyC escape fraction and compare their properties. From the afterglow spectrum of GRB 191004B we determine a neutral hydrogen column density at the LGRB redshift of $\log(N_{\rm HI}/cm^{-2})= 17.2 \pm 0.15$, and negligible extinction ($A_{\rm V}=0.03 \pm 0.02$ mag). The only metal absorption lines detected are CIV and SiIV. In contrast to GRB 050908 and GRB 060607A, the host galaxy of GRB 191004B displays significant Ly$\alpha$ emission. From its Ly$\alpha$ emission and the non-detection of Balmer emission lines we constrain its star-formation rate (SFR) to $1 \leq$ SFR $\leq 4.7$ M$_{\odot}\ yr^{-1}$. We fit the Ly$\alpha$ emission with a shell model and find parameters values consistent with the observed ones. The absolute LyC escape fractions we find for GRB 191004B, GRB 050908 and GRB 060607A are of $0.35^{+0.10}_{-0.11}$, $0.08^{+0.05}_{-0.04}$ and $0.20^{+0.05}_{-0.05}$, respectively. We compare the LyC escape fraction of LGRBs to the values of other LyC emitters found from the literature, showing that LGRB afterglows can be powerful tools to study LyC escape for faint high-redshift star-forming galaxies. Indeed we could push LyC leakage studies to much higher absolute magnitudes. The host galaxies of the three LGRB presented here have all $M_{\rm 1600}>-19.5$ mag, with the GRB 060607A host at $M_{\rm 1600}>-16$ mag. LGRB hosts may therefore be particularly suitable for exploring the ionizing escape fraction in galaxies that are too faint or distant for conventional techniques. Furthermore the time investment is very small compared to galaxy studies. [Abridged]


Introduction
Understanding the nature of the sources responsible for ionizing hydrogen in the intergalactic medium (IGM) remains one of the key challenges in studies of early structure formation. It has been established that active galactic nuclei by themselves provide enough ionizing radiation to keep the Universe reionized at redshifts z < 3 -4 (e.g., Cristiani et al. 2016). According to the predominant view, star-forming galaxies (SFGs) are the main contributors of ionizing radiation at earlier epochs (z 4; e.g., Fontanot et al. 2012;Robertson et al. 2013Robertson et al. , 2015. In order to quantify the contribution of SFGs to the reionization at z 4, we have to know: (i) the galaxy number density as a function of redshift and luminosity; (ii) their ioniz-Results based on observations carried out at ESO Observatory, Paranal, Chile, by the Stargate consortium under Program ID: 0104.D-0600, P.I.: N. Tanvir e-mail: jean-baptiste.vielfaure@obspm.fr Hubble fellow ing photon production efficiency (i.e., the number of hydrogenionizing photons relative to produced UV photons at ∼ 1500 Å, ( f 1500 / f 900 ) int ); (iii) the fraction of produced Lyman continuum (LyC) photons that can actually escape the local environment and ionize the IGM, called the escape fraction f esc . Especially the latter is poorly constrained from observations. Direct searches for LyC emission from galaxies is only possible up to z ∼ 4, beyond which the LyC is unobservable due to the increasing IGM opacity (e.g., Madau 1995). We therefore rely on lower redshift galaxies to be used as proxies for the higherredshift ones. The search for LyC emitters has been quite fruitful in recent years. A few LyC emitters have been found in the local Universe (Bergvall et al. 2006;Leitet et al. 2013;Borthakur et al. 2014). A population of green pea galaxies with strong LyC has been identified at z ∼ 0.3 (Izotov et al. 2016(Izotov et al. , 2018a. A number of LyC-emitting galaxies with typically M UV < −19.5 mag have also been securely identified at higher (z 3) redshifts (Vanzella et al. , 2018Shapley et al. 2016;Bian et al. 2017;Steidel et al. 2018;Fletcher et al. 2019). The growing number of discov-Article number, page 1 of 13 arXiv:2006.09377v1 [astro-ph.GA] 16 Jun 2020 A&A proofs: manuscript no. main ered LyC emitters allows us to look at their common properties that could be used both to make it easier to find new emitters at these redshifts as well as to identify likely emitters at z > 6. The emitters are typically (but not always) found to have a strong Lyα emission line (Verhamme et al. 2015(Verhamme et al. , 2017, a high ratio of [OIII]/[OII] nebular emission lines (Nakajima & Ouchi 2014;Nakajima et al. 2016) and a compact morphology (Izotov et al. 2016(Izotov et al. , 2018a. The compactness of the galaxy as a requirement for LyC photons to escape is also found in theoretical studies Cen (e.g., 2020), together with high star-formation rate surface densities.
Recently, important constraints on the LyC emission at high redshift come from an analysis of stacked spectra of galaxies at −19.5 > M UV > −22 mag at z ∼ 3 (Steidel et al. 2018). The measured mean escape fraction of 0.09 ± 0.01 is dominated by sub-L UV * galaxies, with negligible contribution from L > L UV * systems. These results are in agreement with the hypothesis that fainter galaxies have on average higher escape fractions and account for a large part of the ionizing-photon budget (e.g., Fontanot et al. 2014). Because the bulk of the high-redshift star formation occurred in very faint galaxies (e.g., Bouwens et al. 2015), the regime of LyC emission in those galaxies needs to be studied observationally (see, e.g., Japelj et al. 2017). Deep observations of galaxies in this faint regime (M UV −19 mag) in the LyC portion of their spectra are time-demanding and therefore sparse. Fletcher et al. (2019) identified four faint LyC emission candidates and Amorín et al. (2014) studied a faint lensed galaxy resulting in a poorly constrained limit on its LyC escape fraction.
A unique way to learn more about the escape fractions of faint galaxies is to make use of long gamma-ray bursts (LGRBs).
LGRBs mark the deaths of massive stars (Hjorth et al. 2003;Woosley & Heger 2006;Cano et al. 2017) and can be observed from the local to the z ∼ 8 Universe (Tanvir et al. 2009;Salvaterra et al. 2009).
LGRBs are generally found in sub-L* galaxies at any redshift (e.g., Tanvir et al. 2012;Vergani et al. 2015;Perley et al. 2016;Palmerio et al. 2019). The progenitors are associated with young star-forming regions that contribute hugely to the ionizing photon budget inside galaxies (e.g., Ramachandran et al. 2018). A burst of gamma-ray emission is typically followed by a longer-lived, multiwavelength afterglow (e.g., Kumar & Zhang 2015). Thanks to their brightness, the spectra of the UV/optical/NIR afterglows can reveal the detailed properties of the interstellar medium (ISM) along the lines of sight toward star-forming regions, such as hydrogen column density, dust properties, metallicity and kinematics (e.g., De Pasquale et al. 2003;Thöne et al. 2014;Heintz et al. 2018;Zafar et al. 2018;Bolmer et al. 2019). In addition, the afterglow brightness allows us to search directly for possible leaking LyC in their spectra. With LGRB afterglows we can therefore study the transparency of the LGRB star-forming regions to ionizing radiation Fynbo et al. 2009;Tanvir et al. 2019). Being based on afterglow spectroscopy, LGRB LyC studies do not suffer from galaxy magnitude selection and offer the advantage of a direct determination of the ( f 1500 / f 900 ) int ratio thanks to the intrinsically featureless (power-law) spectra of GRB afterglows.
LGRB afterglows are thus a very useful tool to search for LyC leakage in very faint galaxies at z ∼ 2 − 4.
In this paper we present the detection of LyC photons in the afterglow spectrum of GRB 191004B at z = 3.5055. We further investigate LyC emission in LGRB afterglow spectra, also adding to this case GRBs 050908 and 060607A (z = 3.3467 and z = 3.0749, respectively), the only two previously known GRBs whose spectra show significant LyC emission (Fynbo et al. 2009;Tanvir et al. 2019). This paper is organized as follows. In Sec-tion 2 we present the observations of GRB 191004B. We analyze them in Section 3. In Sections 4 and 5 we determine the escape fraction along the lines of sight to the three GRBs. Section 6 is focused on the Lyα emission modeling. All the results are discussed in Section 7, especially in the context of other known LyC-emitting galaxies. We provide our conclusions in Section 8.
All errors are reported at 1σ confidence and all magnitudes are reported in AB unless stated otherwise. We consider a flat ΛCDM cosmology with the cosmological parameters provided in Planck Collaboration et al. (2016): H 0 = 67.8 km s −1 Mpc −1 , Ω m = 0.308 and Ω Λ = 0.692. In the following we describe our ESO/VLT observations of the GRB 191004B afterglow and host galaxy obtained by the Stargate consortium under Program ID: 0104.D-0600 (P.I.: N. Tanvir). The log of observations is provided in Table 1. Notes. The columns indicate the instrument, the time ∆t after the GRB detection, the total exposure time T exp in each slit/filter, the slit/filter configuration in the case of spectral/imaging observations and the average DIMM seeing at the time of the observation.

VLT/X-shooter afterglow imaging
GRB 191004B was observed with the VLT/X-shooter acquisition camera in the g , r , z bands ∼ 7.2 hours after the burst. The observations were taken in a sequence of 3 × 40s, 3 × 30s and 3 × 60s with the g , r and z filters, respectively. The afterglow is clearly detected in all bands at the coordinates RA(J2000) = 03 h 16 m 49.14 s , Dec(J2000) = −39 • 38 04. 09, with a magnitude of g = 22.41 ± 0.03 mag, r = 21.23 ± 0.02 mag and z = 20.91 ± 0.04 mag. The photometry was calibrated with a single nearby star due to the small field of view. The coordinates of this star RA (J2000)= 03 h 16 m 54.42 s , Dec (J2000)= −39 d 37 41. 9 are taken from the Gaia catalog and the photometry from the Pan-STARRS catalog. . The corresponding afterglow position (error radius r e = 0. 6, as observed in the X-shooter r image ∼ 7.2 hours after the burst) is indicated by a black circle. The slit at PA= −75 • used for the X-shooter afterglow spectroscopy is also overplotted. Galaxy A corresponds to the intervening system at z = 2.137 identified in the VLT/X-shooter spectrum.

VLT/X-shooter afterglow spectroscopy
The first epoch spectrum of GRB 191004B was obtained with the VLT/X-shooter echelle spectrograph (Vernet et al. 2011) just after the g , r , z imaging, ∼ 7.2 hours after the Swift/BAT trigger. The observation consisted of a total integration time of 2400 s in each arm using the 1. 0/0. 9/0. 9JH slits at a position angle of PA= −75 • (see Fig. 1). The observation was carried out under good conditions with an average seeing of ∼ 0. 75 and an airmass of 1.1. We used a nodding scheme with a 5. 0 nodding throw. This observation covers the afterglow emission of the GRB across the wavelength range 300 -2100 nm. We reduced the X-shooter data using version 3.3.5 of the Xshooter data reduction pipeline (Modigliani et al. 2010). All observations of the afterglow and host galaxy, the telluric stars, and the spectrophotometric standards were reduced in the same way. Before processing the spectra through the pipeline, the cosmicray hits and bad pixels were removed following the method of van Dokkum (2001). Then, we subtracted the bias from all raw frames and divided them by the master flat field. We traced the echelle orders and calibrated the data in spatial and wavelength units using arc-line lamps. The flux calibration was done using spectrophotometric standards (Vernet et al. 2009) and a correction for mirror flexure was applied. Lastly, the sky-subtraction and the rectification and merging of the orders was done to obtain the final two-dimensional spectra. To optimally select the extraction region we chose the spatial extension of the afterglow continuum independently in each of the three (UVB, VIS, NIR) arms of the spectrograph. We corrected the 1D spectra for Galactic extinction using the Milky-Way extinction curve of Pei (1992) (assuming the ratio of total-to-selective extinction R V = 3.08) and the extinction map of Schlafly & Finkbeiner (2011). This final data set typically has a signal to noise ratio S/N > 3 per spectral bin and a spectral resolution of R ∼ 6700/11000/7800 in the UVB/VIS/NIR arms, respectively. The resolution was estimated from the average seeing reached during the observation corrected for airmass and wavelength dependence and by interpolating the nominal values available on the ESO instrument website 1 . We also compare these values to the resolutions obtained using the relation derived in Selsing et al. (2019) and we find a good agreement between both. We clearly detect the afterglow continuum from 3675 Å to the end of the X-shooter spectrum. We identify the Lyman-alpha (Lyα) emission and absorption line at ∼ 5474 Å and detect the metal absorption lines described in detail in Sect. 3.1.

Host-galaxy imaging and spectroscopy
We observed the field of GRB 191004B with the VLT/FORS2 (Appenzeller et al. 1998) camera in imaging mode on 2019 December 21, at 03:25:03 UT (∼ 77 days after the Swift/BAT trigger). In the 30 min R-band image (see Fig. 1), we do not detect the host galaxy. The magnitude limit is R ≥ 26.6 mag.
A further VLT/X-shooter spectrum was taken at the afterglow position, once the afterglow faded (∼ 108 days after the Swift/BAT trigger), using the blind offset technique to look for host galaxy emission lines. The average seeing and airmass during the observation were ∼ 0. 6 and 1.2, respectively. The observations were divided in two OBs for a total integration time of 5840/5440/6000 s in the UVB/VIS/NIR arm, respectively. We used the 1. 6/1. 5/0. 9 slits at a position angle of −40 • and a nodding scheme with a 5. 0 nodding throw. The data reduction was carried out in a similar fashion as described in Sect. 2.1. We corrected the 1D spectra for Galactic extinction as for the afterglow spectrum (Sect. 2.2). We clearly detect the Lyα emission line at the same wavelength as in the afterglow spectrum but we do not detect other nebular emission lines. Unlike the afterglow spectra, no continuum is detected in the host-galaxy spectra, so to optimally select the extraction region we chose the spatial extension of the Lyα emission line and applied this 1D extraction throughout the whole spectrum. We provide flux and upper limits of the lines in Sect. 3.1. The good seeing obtained during this observation allows us to reach resolutions of R ∼ 7800/12900/9200 in the UVB/VIS/NIR arms, respectively, calculated similarly as in Sect. 2.2.

HI and metal absorption lines
The X-shooter afterglow spectrum reveals Lyα, Lyβ, Lyγ, and Lyδ, Si iv and C iv absorption lines from the host galaxy gas (D'Elia et al. 2019) (see Fig. 2).
Firstly, we fit the saturated Si iv and C iv absorption with three components using VoigtFit (Krogager 2018). We determine lower limits on the C iv and Si iv column densities of log(N X /cm −2 ) > 15.6 and 14.1, respectively, and a redshift of z = 3.5028, 3.5039, and 3.5055 for the three components (corresponding to relative velocity from the reddest component of ∼ -105 and -180 km s −1 , respectively). The Doppler parameters are b = 25, 20, 30 km s −1 (from the red component to the blue one, respectively).
The high-ionization metal lines such as C iv and Si iv do not necessarily originate from the same place as the neutral gas (e.g., Heintz et al. 2018). However, since they are the only detected for the GRB system, we fixed these components and their redshift to fit the Lyman lines. The Lyman lines are saturated, therefore their fitting could not give a precise estimate of the hydrogen column density. Fixing the b parameters to those inferred from the corresponding metal line components, we find an acceptable fit for log(N HI /cm −2 ) = 17.2 ± 0.3. A fit letting the redshift of the components and the b parameters free to vary would give consistent N HI values. This is one of the lowest hydrogen column densities measured among the LGRB afterglow sample (Tanvir et al. 2019) and places the absorbing system on the boundary between being a Lyα forest absorber (log(N HI /cm −2 ) < 17; Rauch 1998) or a Lyman-limit system (LLS, 17 < log(N HI /cm −2 ) < 20.3; Péroux et al. 2003). The reddest component is the richest in neutral hydrogen, whereas the bluest component is the strongest for the high-ionization lines.
We do not detect the other common singly ionized metal lines such as Fe ii, Mg ii, Al ii or Si ii (Christensen et al. 2011;de Ugarte Postigo et al. 2012), that are typically detected in LGRB afterglow spectra. The part of the spectrum at the corresponding wavelengths of most of the low-ionization lines is quite noisy or contaminated by residual sky lines. We place approximate upper limits on their rest-frame equivalent widths of ∼ 0.4 Å. We can determine a 3σ upper limit on the C ii column density of log(N HI /cm −2 ) > 14.7. Low column densities of low-ionization metal transitions are also found for a part of the LLS sample detected through quasar spectroscopy (see, e.g., Fox et al. 2013;Cooper et al. 2015), as well as in the other few LGRBs having log(N HI /cm −2 ) < 19 (see, e.g., Thöne et al. 2011). In these systems a significant fraction of the gas is ionized, whereas usually LGRB hosts show higher HI absorption values (median of log(N HI /cm −2 ) = 21.6; Tanvir et al. 2019), typical of systems probing neutral gas. At such low HI column densities, even for approximately solar metallicity gas, the strongest low-ionization metal absorption features can have column densities below the detection threshold. De Cia et al. (2011) found very low values for the column densities of low-ionization metals detected in GRB 070125 afterglow spectra, at z = 1.5477. They conclude that the LGRB was likely located in the outskirts of a massive star-forming region inside a faint and small host galaxy (see also Updike et al. 2008;Cenko et al. 2008).
We detect Fe ii, Mg ii absorption lines and [O iii] emission corresponding to an intervening absorber at z = 2.137. This system is at 3. 8 from the afterglow position and corresponds to galaxy A in Fig. 1. Its emission trace in the 2D spectrum covers 1. 5. We do not detect the presence of another galaxy brighter than R=26.6 mag in a radius of 3. 0 around the position of the GRB 191004B afterglow, corresponding to 22 kpc at this redshift. Therefore, we can rule out contamination of the observed LyC emission from foreground interlopers.
We detect several Lyα absorbers in the Lyα forest from z = 2.8 to z = 3.45, with hydrogen column densities spanning log(N HI /cm −2 ) ∼ 15 − 17. The lack of LyC emission below ∼ 3650 Å is due to the Lyman limit absorption by some of these systems, and the dimming of the LyC emission around 3800 Å corresponds to the Lyman limit of a system at z = 3.19 and to the likely strong Lyα absorption of Galaxy A (see Sect. 3.1).
The host-galaxy spectrum does not show LyC emission above a 3σ flux density limit of 4.5 × 10 −18 erg s −1 cm −2 Å −1 . A rough binning of the spectrum in the 820 -910 Å range allows to place a more stringent 3σ upper limit on its LyC flux density to 6.0 × 10 −19 erg s −1 cm −2 Å −1 .

Line-of-sight extinction
We measured the host-galaxy extinction along the line of sight to the GRB directly from the afterglow spectrum corrected for Galactic extinction (see Sect, 2.2) following the procedure as described in Japelj et al. (2015) and Zafar et al. (2018). Due to the sparse X-ray afterglow observations around the epoch of the X-shooter spectrum, we did not include X-ray data into this fit. We ignored the spectral region blueward of the Lyα emission/absorption feature in the fit, since it is significantly affected by Lyα forest absorption. The edges of the VIS and NIR spectra and the regions affected by telluric absorption were also masked out. We cleaned the spectra of absorption lines by recursively fitting a polynomial to the continuum and removing everything that deviate from the fit by more than 3.5σ. Following the procedure of Japelj et al. (2015), we fitted the spectrum with a powerlaw function modulated by dust extinction assuming the Milky Way and Magellanic Clouds extinction curves from Pei (1992). We found that the spectrum can be well described by a powerlaw with small extinction. The best fit gives an optical spectral index β o = 0.49 ± 0.03, an extinction A V = 0.03 ± 0.02 mag and χ 2 /dof = 1.007 for the Small Magellanic Cloud extinction curve. We present the best fit (solid orange line) together with the corresponding intrinsic power-law corrected for extinction (dashed line) in Fig. 4.

Lyα emission
We detect the Lyα emission line in both the X-shooter afterglow and host-galaxy spectra (see Fig. 5). The emission line is asymmetric. We fit its profile with a skewed Gaussian parametrized as in Claeyssens et al. (2019); Matthee et al. (2019), where A is the amplitude of the Lyα line; λ 0 is its peak wavelength; a asym is the measure of asymmetry of the line and d the parameter that controls the line width. These parameters are left free during the fitting with broad uniform priors. In order to determine both a robust fit and its associated uncertainty, we use the Markov Chain Monte Carlo sampler emcee (Foreman-Mackey et al. 2013) to maximize the likelihood function using Eq.
(1). The best fit of the X-shooter Lyα line detected in the host spectrum provides a peak position of λ 0 = 5478.85 +0.13 −0.12 Å, an asymmetry of a asym = 0.33 +0.05 −0.07 and a width parameter of d = 0.6985 +0.110 −0.103 . We convert the latter into FWHM according to the analytic expression: The corresponding rest-frame FWHM of the line is FWHM 0 = 104 ± 30 km s −1 . The flux of the Lyα line corrected for Galactic extinction is F(Lyα) = (1.0 ± 0.1) × 10 −17 erg s −1 cm −2 . This corresponds to a Lyα luminosity of L Lyα = (1.18 ± 0.12) × 10 42 erg s −1 and a rest-frame equivalent width of EW 0 (Lyα) = 7.4 ± 2.6 Å. We convert this Lyα flux to a star formation rate of SFR(Lyα) ≈ 1 M yr −1 by using the recombination factor 8.7 from the case B of the theory of recombination (Brocklehurst 1971) and the relation from Kennicutt (1998), scaled to the Chabrier (2003) initial mass function. Nevertheless, we know that this SFR can be largely underestimated due to the complex scattering and destruction process of Lyα photons in neutral gas and dust. Therefore the SFR derived here only represents a lower limit.

Other emission lines
We do not detect other emission lines at the afterglow position in the afterglow and host-galaxy spectra. We determine 3σ upper limits of 0.3 × 10 −17 erg s −1 cm −2 for Hβ and [OII] λ3727 , and 1.0 × 10 −17 erg s −1 cm −2 Å −1 for [OIII] λ5007 . By assuming that the negligible dust extinction in the line of sight is also valid for the integrated host-galaxy dust content, the Hβ upper limit can be converted to an SFR < 4.7 M yr −1 considering the intrinsic Balmer decrement Hα/Hβ = 2.86 (Osterbrock 1989) and using the relation from Kennicutt (1998), scaled to the Chabrier (2003) initial mass function (IMF). This value is consistent with the lower limit derived from the Lyα line flux and allows us to put an additional constraint on the SFR. This implies an SFR of a few M yr −1 and an escape fraction of Lyα photons f esc (Lyα) > 0.13.

LyC escape fraction
The LyC escape fraction is the ratio of the observed flux below the Lyman limit, corrected for the Milky Way extinction and the IGM transmission, and the intrinsic flux. Usually in galaxy studies this last quantity is inferred from the galaxy SED (and therefore it is model dependent), and the IGM transmission used is an average value obtained from line-of-sight simulations.
Taking advantage of afterglow spectroscopy it is possible to determine f esc (i): directly using the N HI value derived by fitting the Lyα absorption line or the envelope of the residual , as commonly done in galaxy studies. The former case is a direct determination from the data. The fact that the afterglow continuum is a simple power-law is an advantage (in both cases) because the intrinsic f 900 is much less model dependent than in galaxy studies.

Direct determination
Following Siana et al. (2007), we define f esc as the fraction of emitted Lyman continuum photons that escapes into the IGM. This corresponds to the formula, where τ LL is the optical depth at wavelengths below the Lyman limit and is given by (see also Prochaska et al. 2010), where λ LL ≈ 912 Å is the Lyman limit wavelength. In principle we can directly use the value of N HI obtained by the fit of the HI absorption. In the case of GRB 191004B, the HI column density belongs to the flat part of the HI curve of growth, therefore it cannot be determined precisely. Considering the value and errors determined in Sect. 3.1, we find f esc = 0.43 +0.22 −0.24 . The large error range reflects the N HI uncertainties.
We can determine more precise N HI and f esc values in the following way. The observed afterglow spectrum blueward of the Lyman limit is the afterglow continuum (a power-law modulated by extinction, i.e. the solid line in Fig. 4), heavily absorbed by the neutral hydrogen in the host galaxy (Eq. (4)), and furthermore absorbed to a much lesser extent by Lyman alpha forest at the corresponding redshift. From the shape of the observed spectrum and the modeled afterglow continuum we can use Eq. (4) to directly determine the opacity in the host galaxy's line-of-sight by varying the N HI (see Fig. 6). Because the observed spectrum is partly absorbed by the Lyman alpha forest we don't aim to match the spectral continuum but instead its envelope. The best match is obtained with log(N HI /cm −2 ) = 17.20±0.15, which corresponds to f esc = 0.43 +0.12 −0.13 . The errors reflect the uncertainty on the absorption due to the forest. If we take into account the Fig. 7. Distribution of relative escape fraction ( f esc,rel (LyC)) for GRB 050908, GRB 060607A and GRB 191004B (from left to right) with and without uncertainty on the observed measure of ( f 1500 / f 900 ) obs from the afterglow spectrum (respectively gray and black distribution). The high values with f esc,rel (LyC) > 1 are non-physical and due to low IGM transmission, hence we restrict the statistics and axis range to 0-1 for f esc,rel (LyC). For GRB 191004B, 57% of the gray distribution is at f esc,rel (LyC) ≤ 1 while the fraction is 92% and 96% for GRB 050908 and GRB 060607A, respectively. The solid red and orange lines in each panel represent the median and average value of this f esc,rel (LyC) ≤ 1 part of the distribution. The dashed red lines represent the corresponding 16th and 84th percentiles. The f esc, abs are similar to the values quoted here for GRB 050908 due to the null extinction hypothesis. For GRB 191004B and GRB 060607A the corresponding absolute values are derived in Sect. 4.2 and 5. effect of dust extinction on the afterglow continuum, we obtain f esc, abs = 0.35 +0.10 −0.11 . We can also determine a lower limit to f esc and N HI from the observed continuum over the same wavelength range, as it corresponds to the flux of the escaping photons dimmed by the IGM opacity. We obtain log(N HI /cm −2 ) < 17.4, f esc > 0.25, and f esc, abs > 0.22 taking dust extinction into account.

Constraints from IGM simulations
In galaxy studies, the common way to estimate the fraction of escaping ionizing photons is by calculating the ratio between the observed fraction of escaping LyC photons (900 Å), corrected for IGM transmission, relative to the fraction of escaping nonionizing photons (1500 Å) (Steidel et al. 2001;Siana et al. 2007). The relative escape fraction can be expressed as: where ( f 1500 / f 900 ) int is the intrinsic flux density ratio, ( f 1500 / f 900 ) obs is the observed flux density ratio and T IGM 900 = e −τ IGM is the IGM transmission factor at 900 Å along the sightline, and τ IGM the line-of-sight opacity to the IGM for LyC photons. As shown in Siana et al. (2007), this formula is equivalent to Eq. (3).
The absolute escape fraction, f esc,abs is the ratio of the escaping to intrinsic ionizing flux density. Knowing the dust attenuation A 1500 , f esc,abs can be written as (e.g., Inoue et al. 2005;Siana et al. 2007), f esc,abs = f esc,rel 10 −0.4 A 1500 .
The ( f 1500 / f 900 ) int ratio is very difficult to constrain observationally and is usually estimated by using spectral synthesis models (e.g., Bruzual & Charlot 2003). Therefore, this value has generally significant uncertainties due to the assumptions on which these models rely, stellar population age, metallicity, starformation history and IMF. In our case, we have the possibility to directly measure the intrinsic flux density ratio as the intrinsic afterglow emission corresponds to a single power-law. This power-law is estimated as described in Sect. 3.3 and extrapolated to the ionizing UV domain. We then calculate the average intrinsic flux density in the same rest-frame range as for the observed flux densities ( f obs 1500 and f obs 900 , Sect. 3.2). We finally derive an intrinsic flux ratio of ( f 1500 / f 900 ) int = 0.43 from the best fit corrected for extinction (dashed orange curve in Figure 4).
To recover the T IGM 900 factor, we simulated a large number of sightlines using observational constraints on the properties of the IGM probed in intergalactic absorbers. We used the same Monte Carlo (MC) simulation as in Japelj et al. (2017) and described in Vanzella et al. (2015) to generate 10 000 line-of-sight transmissions from a source at z = 3.5. For each simulated sightline we calculate the average transmission over the same rest-frame range as used for f obs 1500 and f obs 900 , we then correct the observed ratio ( f 1500 / f 900 ) obs for this quantity and calculate f esc,rel (LyC). In order to take into account the uncertainty on ( f 1500 / f 900 ) obs , for each of the 10 000 sightlines we randomly select a value in a normal distribution defined by its measurement (center) and its error (width). The resulting distribution of f esc,rel is shown in grey in Fig. 7. In the case of GRB 191004B, the best fit of the afterglow provides an estimate of the extinction and hence the possibility to calculate the absolute escape fraction. In the following, we report the absolute values derived using Eq. (6) but for a comparison we also provide the relative ones, which do not rely on a dust-extinction model (see Sect. 3.3).
For the estimation of the escape fraction we proceed similarly to what is done in Shapley et al. (2016) in their study of the direct detection of LyC emission from a galaxy at z ∼ 3. We only keep the lines of sight resulting in an escape fraction of f esc ≤ 1, since the values higher than one are not physical. We find that 57% of the initial 10 000 sightlines are at f esc ≤ 1 and we measure that 95% of the distribution is at f esc, abs ≥ 0.31 ( f esc, rel ≥ 0.35). The median escape fraction is f esc, abs = 0.46 +0.23 −0.11 ( f esc, rel = 0.52 +0.26 −0.12 ), with an average value of f esc, abs = 0.50 ( f esc, rel = 0.56), see also Table 2.

>0.31
Notes. The columns correspond to the name of the GRB; z syst : the systemic redshift of the host galaxy; f HI esc,abs (rel) : the absolute (abs) and relative (rel) escape fraction determined from the estimation of the HI column density, see Sect. 4.1;f esc,abs (rel) : the absolute (abs) and relative (rel) lower limit of the escape fraction of ionizing photons. The relative value is given between brackets when a non-zero extinction is considered, otherwise (*) means same value as the absolute one (see Sect. 4.2); f lim esc,abs : the absolute lower limit of the escape fraction of ionizing photons (see Sect. 4.2).

LyC leakage among LGRB host galaxies
To date, in addition to GRB 191004B presented in this paper, there are only two more LGRBs for which direct LyC emission was detected, namely GRB 050908 and GRB 060607A. Their afterglow spectra have already been presented in Fynbo et al. (2009) and briefly discussed in Tanvir et al. (2019). Here, we summarize the main known features of these LGRB afterglows and host galaxies and derive their f esc (LyC), following the procedures of Sect. 4. We stress that, apart from their afterglow optical spectrum showing low HI column densities associated with the LGRB host galaxy, the general properties of these LGRBs (prompt emission, energetic, afterglow light curves) are similar to that of the general LGRB population (see, e.g., Kann et al. 2010).
GRB 050908 occurred at z = 3.3467. Its afterglow spectrum was observed with VLT/FORS1 (Fynbo et al. 2009), KeckII/Deimos and Gemini North/GMOS , and showed very strong high-and low-ionization absorption lines at the GRB redshift. In this paper we use the spectrum obtained by VLT/FORS. In addition to the absorptions associated with the GRB host galaxy, two strong Mg ii intervening absorbers are identified in the GRB afterglow spectrum at z = 0.6915 and 1.4288. The emission counterparts of these absorbers have been identified by Schulze et al. (2012) a few tens of kpc away from the GRB host. Deep HST-F775W images of the field allowed the detection of the GRB host galaxy with m AB = 27.67 mag (Hjorth et al. 2012;Blanchard et al. 2016). The VLT/FORS spectrum covers the ionizing domain at wavelengths > 818 Å rest-frame. A sufficient signal-to-noise ratio is only reached down to 880 Å. We measure the ionizing ( f obs 900 ) and non-ionizing UV ( f obs 1500 ) flux density from this spectrum and find a ( f 1500 / f 900 ) obs ratio of 6.9 ± 0.5. Considering a null extinction (A V = 0 mag), according to the results published in Vergani (2008), we fit the afterglow spectrum of GRB 050908 with a simple power-law and following the method described in Sect. 4.1, we derive ( f 1500 / f 900 ) int = 0.35 with f int 900 measured over 880-910 Å. Similarly to GRB 191004B we can derive the absolute escape fraction of ionizing photons which is equal to the relative one in this case, since no extinction correction is applied.
The distribution of f esc,rel (LyC) is presented in Fig. 7 and shows that 95% of the distribution is at f esc, abs ≥ 0.07 and 92% with f esc, abs ≤ 1. The median escape fraction for the latter 92% is 0.10 +0.08 −0.02 , with an average value of f esc, abs = 0.15. This estimation is fully consistent with the value of f esc = 0.08 +0.05 −0.04 derived by Tanvir et al. (2019) from the HI column density (log(N HI /cm −2 )= 17.6), following the method presented in Sect. 4.1.
GRB 060607A occurred in a galaxy at z = 3.0749. It has the lowest neutral-hydrogen column density observed among the LGRB afterglow sample (log(N HI /cm −2 ) = 16.95 ± 0.03; Tanvir et al. 2019) and, like GRB 191004B, the only metal absorption lines at the GRB redshift are those of C iv and Si iv. The VLT/UVES spectrum of its afterglow shows a very complex line of sight, rich in intervening absorbers Prochaska et al. 2008;Fynbo et al. 2009). There is evidence of LyC emission over ∼ 16 Å below the LL, likely absorbed blueward from an intervening system. It is not firmly possible to exclude that the apparent LyC emission is due to an interloper. Nonetheless, deep HST-F775W imaging observations are available for the GRB 060607A field (see Blanchard et al. 2016) and no galaxy is detected at m AB > 28.9 mag over a radius of ∼ 15 kpc from the afterglow position. From the UVES spectrum, we measure the ionizing ( f obs 900 ) flux density over the 16 Å where LyC emission is detected and we find a ( f 1500 / f 900 ) obs ratio of 5.0±0.6. The determination of dust extinction on this line of sight is not trivial as the afterglow experienced some rebrightening during the Xray observations and UVES spectroscopy. Following the same prescriptions as in Sect. 3.3, we find a quite high extinction with A V = 0.13±0.04 mag for the Small Magellanic Cloud extinction curve (see also Kann et al. 2010). Following the same method as for the previous cases, we derive ( f 1500 / f 900 ) int = 0.42. The distribution of f esc,rel (LyC) is presented in Fig. 7 and shows that 95% of the distribution is at f esc, abs ≥ 0.05 ( f esc, rel ≥ 0.09) and 96% with f esc,rel ≤ 1. The median escape fraction for the latter 96% is f esc, abs = 0.08 +0.03 −0.02 ( f esc, rel = 0.13 +0.04 −0.03 ) with an average value of f esc, abs = 0.08 ( f esc, rel = 0.14). The value derived from the HI column density, following the method presented in Sect. 4.1, is f esc = 0.45 ± 0.13 and f esc, abs = 0.20 ± 0.05, with log(N HI /cm −2 ) = 17.1 ± 0.15 (in agreement with the column density reported in the literature). The difference of the f esc values determined by the two methods could be due to a higher T IGM 900 average from simulations, as this line of sight is particularly rich in absorbing systems. We stress that both f esc, abs values are tentative due to the difficulties in determining the dust extinction.
The f esc results for the three GRBs are summarized in Table  2.

Statistical estimation of the ionizing escape fraction in
LGRB host galaxies The HI column density of log(N HI /cm −2 ) = 17.2 ± 0.15 derived for GRB 191004B corresponds to one of the lowest N HI measured for an LGRB host galaxy. Adding this case to the sample studied in Tanvir et al. (2019) we re-evaluate the average escape fraction derived from the compilation of the 141 column densities available to-date. To do so, we use the formula described in Sect. 4.1 (Eq. (3)) averaged over the full sample and we find f esc = 0.007 instead of 0.005 presented in Tanvir et al. (2019). Following again the method used in  and Tanvir et al. (2019) we performed a bootstrap resampling by allowing replacement of the 141 sightlines in the above sample. From the 10 6 N HI distributions simulated with this method and taking into account scatter produced by the uncertainty on each data point, we find a 98% c.l upper limit of f esc < 0.020. This new estimation increases somewhat the result found in Tanvir et al. (2019) but is still consistent with their f esc upper limit in the range 0.015 -0.02.

Modeling of the Lyα emission line
The non detection of nebular emission lines in the spectra of GRB 191004B prevents firm determination of the exact redshift of the host galaxy. We can, however, use the empirical correlation found by Verhamme et al. (2018) between the FWHM and the velocity shift of the peak of the Lyα line. We derive a predicted redshift of z predict = 3.5060 ± 0.0008. This value is fully consistent with the HI component seen in absorption at z = 3.5055. Due to its resonant nature, the Lyα radiation produced by active star-forming regions carries the signature of the physical properties of the neutral gas where the photons scatter before escaping. In that sense, the detection of the Lyα line gives us a unique chance to investigate if the properties probed by the GRB afterglow lines of sight are similar to the ones of the material probed by the Lyα emission. To do this comparison, we perform a shell-model fitting of the observed Lyα line using an improved model grid and process originally described in . The shell model is a commonly used model (first introduced by Ahn et al. 2003, also see Verhamme et al. 2006) that consists of a central source emitting a continuum and a Gaussian Lyα line defined by an intrinsic width and equivalent width (σ i , EW i ). This source is surrounded by a shell of neutral hydrogen and dust described by four parameters: a radial expansion velocity v exp ; a HI column density N HI ; an effective temperature of the gas T; an optical depth of dust close to Lyα wavelength τ d . We fit the X-shooter host-galaxy data leaving these parameters free in the grid range values described in . As the final fitting result is sensitive to the redshift of the emitting source, we also leave the systemic redshift as a free parameter with a Gaussian prior of z syst = 3.5060 ± 0.01, based on the estimation derived above. Note, however, that the prior on z syst is very wide, and thus, essentially a free parameter.
We compare the observed Lyα profile to the best-fit model in Fig. 8. We find that the best-fit profile describes the observation very well. Indeed, the normalized residuals between the observed and modeled Lyα line (brown dots in Fig. 8) provides a Shapiro-Wilk coefficient of 0.985 and a p-value = 0.316, therefore consistent with Gaussian noise. The best-fit results show a small blue component at -210 km s −1 , consistent with the spectrum. This blue component is only detected at a 2.5σ confidence level, with a flux of F(Lyα, blue) = (0.5 ± 0.2) × 10 −17 erg s −1 cm −2 . If real, the separation between the blue and red peaks is 323 km s −1 .
The best-fit parameters are consistent with the observed values, (i): the systemic redshift z syst = 3.5053 +0.0005 −0.0003 is in good agreement with the value of the reddest HI absorption component seen at z = 3.5055 and also with z predict = 3.5060 ± 0.0008; (ii): the low expansion velocity of the gas v exp = 40.0 +12.0 −11.0 km s −1 is consistent with the velocity of the reddest HI component which also is the richest in neutral hydrogen (cf. −33 km s −1 , EW i (Lyα) = 24 +13 −9 Å) imply a rest-frame EW(Hβ) = 1.0 +0.6 −0.4 Å which is fully consistent with its non detection in the X-shooter spectrum (EW(Hβ) < 3.0 Å, see Sect. 3.5).
We stress that, while it is well known that the shell-model can reproduce observed Lyα spectra extremely well (e.g., Gronke 2017), its clear oversimplification of realistic HI configurations in and around galaxies is a source of debate in the literature parameters (e.g., Gronke et al. 2017;Orlitová et al. 2018).

Comparison with Lyα-emitter LGRBs (LAE-LGRBs)
The Lyα line profile of GRB 191004B is typical of what is observed for known LAE-LGRBs. It shows a main peak redshifted by 113 ± 30 km s −1 which is just below the range of values (150-750 km s −1 ) determined for the known LAE-LGRB host-galaxies in Milvang-Jensen et al. (2012). The low-velocity shift of the line associated with LyC leakage is consistent with the theoretical predictions of Verhamme et al. (2015). Indeed, they found that the classical asymmetric redshifted profile of the Lyα line should be shifted by less than 150 km s −1 in the case of LyC leakage. This was also observationally confirmed in Verhamme et al. (2017) for low-redshift LyC emitters. The host galaxy has a Lyα luminosity (L Lyα = (1.18 ± 0.18) × 10 42 erg s −1 ) and a UV magnitude lower limit (M UV > −19.2 mag) found to be close to the median values of the LAE-LGRB host-galaxies sample, respectively, L Lyα = 1.7 +3.5 −0.9 × 10 42 erg s −1 (Fynbo et al. 2003;Jakobsson et al. 2005;Milvang-Jensen et al. 2012) and M UV = −19.9 +1.4 −1.0 mag (Tanvir et al. 2019). However, GRB 191004B is the only LAE-LGRB for which LyC leakage has been detected (over the eight LAE-LGRBs with a spectrum covering the spectral range beyond the Lyman limit) or inferred from the HI column density probed by the afterglow. It is in fact the only LAE-LGRB having log(N HI /cm −2 ) < 18, over the twenty-two LAE-LGRBs with an afterglow spectrum covering the Lyα absorption spectral range. Indeed, LGRBs usually probe dense HI column densities. Among the 140 LGRB lines of sight of the sample of Tanvir et al. (2019) only two have log(N HI /cm −2 ) < 18, namely GRB 050908 (log(N HI /cm −2 ) =17.6) and GRB 060607A (log(N HI /cm −2 ) =16.95).
Consistently, together with GRB 191004B, they also are the only cases showing nonzero afterglow continuum emission below the Lyman limit. Intriguingly they do not show any trace of Lyα emission above a 3σ flux upper limit of F(Lyα) = 6.4 × 10 −18 erg s −1 cm −2 and F(Lyα) = 7.3 × 10 −18 erg s −1 cm −2 , respectively . We should consider, however, that those are extremely faint galaxies. Even if no firm conclusion can be drawn from such a small sample, the fact that the highest escape fraction is obtained for GRB 191004B, which is also the only LAE, seems to support the idea that strong LyC emission is correlated with Lyα emission. Inversely, substantial Lyα emission is not necessarily associated with LyC leakage, as GRB 191004B is the only LyC leaker among eight LAE-LGRBs of comparison.
LyC and Lyα escape are expected to be somehow correlated (see, e.g., Verhamme et al. 2015Verhamme et al. , 2017Dijkstra et al. 2016;Kakiichi & Gronke 2019) since the escape of both radiations is facilitated by low-N HI sightlines. Nevertheless, there exist physical phenomena depending mainly on the density of the HI gas and the presence of low density paths in the medium, that can boost the Lyα against LyC escape, or inversely (see, e.g., Ji et al. 2020, Sect. 3.2). Indeed, recent studies such as the one of Bian & Fan (2020) show that strong Lyα emission is not necessarily a good indicator of LyC leakage (see also e.g., Guaita et al. 2016;Grazian et al. 2017). Observationally, the Lyα EW does seem to correlate with LyC escape as shown by Steidel et al. (2018). There are no deep enough observations to put significant Lyα EW limits for GRB 050908 and GRB 060607A hosts, but in general, galaxies with comparably high LyC f esc as GRB 191004B have higher Lyα EW 0 (e.g. Izotov et al. 2016Izotov et al. , 2018aVanzella et al. 2016). A relation between the separation of the observed Lyα emission blue and red peaks and the LyC escape fraction has also been found by e.g. Verhamme et al. (2015Verhamme et al. ( , 2017 and Izotov et al. (2018b) for low-redshift LyC leaking LAEs. The separation of ∼ 300 km s −1 found for GRB 191004B would correspond to f esc ∼ 10%, thus, lower than the value we determined. In our case, we are probing the LyC leakage along the GRB line of sight. Therefore, it is possible that the average LyC leakage of the galaxy and the N HI through which Lyα photons escape are different than those probed by the GRB afterglow. A higher average N HI value, as the best value of the shellmodel fitting (log(N HI /cm −2 ) = 17.7), would imply an f esc (LyC) of 6.9%, in better agreement with the Lyα EW 0 and peak separation relations. In this case, our results could support a medium with low density holes and the escape would not be homogeneous. This is a more realistic scenario, as the ISM in galaxies is not homogeneous. We would like to point out also that the re-lations above have not been tested for faint, high-redshift galaxies as those studied in this paper. Steidel et al. (2018) explore brighter galaxies, whereas Fletcher et al. (2019) study systems with Lyα EW 0 > 15 Å.

Comparison with known LyC emitters
We collected absolute escape fraction measurements from the literature, especially from galaxies lying at z 3. The comparison between our LGRB measurements and those from the literature is presented in Fig. 9 as a function of galaxy UV luminosity.
The data in the plot is not homogeneous, owing to the difficulty of finding LyC emitters and measuring their escape fractions. The IGM transmission T IGM 900 depends on the assumed model; Steidel et al. (2018) and Shapley et al. (2016) use a combined transmission of the IGM and circumgalactic medium (CGM), resulting in on average lower transmission than our model (and that of Vanzella et al. 2016). Fletcher et al. (2019) adopt an average transmission at z ∼ 3. If they used a similar probabilistic approach as in this work, the measurements from their study would have lower values, and they would be only a lower limit. The choice of ( f 1500 / f 900 ) int also varies based on different galaxy population synthesis studies. Finally, we overplot the sample of z ∼ 0.3 green pea LyC emitters. Their low-redshift origin makes the measurements less affected by the uncertainties in T IGM 900 and ( f 1500 / f 900 ) int . Fig. 9 shows that LGRBs can be powerful tools to probe LyC leakage from extremely faint high-redshift galaxies, making it possible to explore LyC leakage from galaxies beyond the current absolute magnitude limits of galaxy studies. The number of objects at M 1600 > −19.5 mag is similar to those of galaxy studies and the power to access even fainter galaxies is evident. LyC leakage investigations through LGRB afterglows do not suffer from galaxy apparent magnitude selection. Furthermore the results are obtained as a by-product of routinely performed observations for the GRB redshift determination, requiring a very small observing time (∼ 1 -3 hours).
We stress that there is a fundamental difference between our study, using LGRB afterglows, and the common ones using galaxy observations. In our case we look into the transparency of the galactic material to the ionizing radiation in one line of sight (toward the studied LGRB), while the other studies measure the averaged transparency of a whole galaxy facing us. Therefore through LGRB afterglow we can probe directly the lines of sight through which LyC photons escape.
The use of LGRBs to constrain the leakage of LyC photons s can be a powerful way to understand what sources reionized the Universe. Indeed, LGRBs can probe f esc not only at higher redshift -thus, requiring less extrapolation in terms of galactic properties -but also (UV) fainter galaxies. This is particularly noteworthy since those galaxies are thought to be the main driver of reionization because (i) models and analytical considerations show that their feedback is likely overruling gravity, thus, strong enough to clear channels through which LyC photons can escape (Paardekooper et al. 2015;Cen 2020), and (ii) the (increasing) steepness of the Lyα luminosity function hints toward a very steep LyC luminosity function Dressler et al. 2015;Dijkstra et al. 2016). A larger sample of LGRB afterglow spectra will allows us to test these theories.

Conclusions
We investigated LyC leakage in LGRB afterglow spectra through the observations of GRB 191004B at z = 3.5055 (presented here for the first time), GRB 050908 (z = 3.0749), and GRB 060607A (z = 3.3467). We found absolute escape fractions along their sightlines of 0.35 +0.10 −0.11 , 0.08 +0.05 −0.04 ,0.20 +0.05 −0.05 , respectively. By using simulations of the IGM opacity, we found similar results except for GRB 060607A likely because of its very rich line of sight.
Thanks to the fact that the intrinsic afterglow emission corresponds to a single power-law, the determination of the intrinsic flux density ratio between the fraction of escaping ionizing and non-ionizing photons does not suffer from the uncertainties due to the assumptions behind galaxy spectral-synthesis models. Furthermore afterglow-based studies do not suffer from galaxy magnitude selection effects as galaxy studies.
The results presented here are by-products of routinely performed observations to determine the redshift of GRBs obtained using only ∼ 1-3 hours each, compared to the many hours needed for galaxy-based studies to pre-select and observe such faint galaxies. They show that LGRBs can be powerful tools to study LyC leakage from faint, star-forming galaxies at high redshift. Indeed, the host galaxies of the LGRBs presented here all have M 1600 > −19.5 mag, with the extreme case of the host of GRB 060607A at M 1600 > −16 mag. Such faint galaxies are very common at very high redshift. Their role in reionization is still debated, but their global ionizing photon budget may contribute significantly to the reionization process.
LGRBs explode in young star-forming regions, major contributors of the LyC photons inside galaxies. LGRB afterglows probe lines of sight originating from these regions, therefore they shed light on the paths through which LyC photons escape. Uniquely LGRBs also have the potential to combine the study of these lines of sight (through afterglow spectroscopy) with those of the global properties of their hosting galaxies, through photometric and spectroscopic observations once the afterglow has disappeared. The galaxies presented in this paper are faint and at the limit of current instrumentation, but they will be perfect targets for E-ELT and JWST. Future observations will allow us to determine their characteristics and to define indirect indicators to be used to find similar galaxies at the reionization epoch.
The current limitation of the LGRB contribution to reionization studies is the limited number of LGRBs with afterglow detections, especially at high redshift. Future satellites such as Gamow Explorer (White 2020) and THESEUS 2 (pre-selected as M5 European Space Agency mission), will largely improve this situation. Indeed the latter will allow the detection and redshift estimate of about 150 GRBs at 3 < z < 4, during its 3-year nominal duration (Amati et al. 2018;Ghirlanda et al. 2015). Furthermore, its detection of ∼ 100 LGRBs at z > 6, combined with space and ground based observations will allow us to directly study the galaxies contributing to reionization.