The Origin of Low-redshift Event Rate Excess as Revealed by the Low-luminosity Gamma-Ray Bursts

The relation between the event rate of long gamma-ray bursts and the star formation rate is still controversial, especially at the low-redshift end. Dong et al. confirmed that the gamma-ray burst rate always exceeds the star formation rate at a low redshift of z < 1 in spite of the sample completeness. However, the reason for low-redshift excess is still unclear. Since low-luminosity bursts are at smaller redshifts generally, we choose three Swift long burst samples and classify them into low- and high-luminosity bursts in order to check whether the low-redshift excess is existent and if the excess is biased by the sample size and completeness. To degenerate the redshift evolution from luminosity, we adopt the nonparametric method to study the event rate of the two types of long bursts in each sample. It is found that the high-luminosity burst rate is consistent with the star formation rate within the whole redshift range, while the event rate of low-luminosity bursts exceeds the star formation rate at a low redshift of z < 1. Consequently, we conclude that the low-redshift excess is contributed by the low-luminosity bursts with possibly new origins unconnected with the star formation, which is also independent of the sample size and the sample completeness.


Introduction
Gamma-ray bursts (GRBs) can be traditionally classified into long and short classes according to the duration (T 90 ), with a dividing line of ∼2 s (e.g., Kouveliotou et al. 1993;Zhang & Choi 2008;Wang et al. 2020) or ∼1 s (Bromberg et al. 2013;Gruber et al. 2014;Zhang et al. 2020b).The exact dividing line between long and short GRBs depends on the sample and the energy band.For example, the T 90 values of Swift GRBs with good duration measurement (Gruber et al. 2014) and wellmeasured spectrum (Zhang et al. 2020b;Deng et al. 2022) are bimodally distributed with a boundary of T 90 ∼ 1 s.In any case, it proves that the apparent middle class in some T 90 -selected GRB samples is inexistent or artificial (e.g., Zhang et al. 2016;Tarnopolski 2019aTarnopolski , 2019b;;Zhang et al. 2020b;Deng et al. 2022;Tarnopolski 2022).Long GRBs (lGRBs) with T 90 > 2 s are believed to come from the core collapse of massive stars (e.g., Woosley 1993;Paczyński 1998;Hjorth et al. 2003;Stanek et al. 2003;Woosley & Bloom 2006;Zhang 2018), so that the long GRB (lGRB) rate was frequently compared with the star formation rate (SFR) in a low-metallicity environment (Wang & Dai 2009;Elías & Martínez 2018;Palmerio & Daigne 2021).On the other hand, short GRBs (sGRBs) with T 90 < 2 s are usually interpreted as the mergers of compact objects such as double neutron stars (DNS) or neutron starblack hole (NS-BH) binaries (e.g., Eichler et al. 1989;Paczyński 1991;Narayan et al. 1992).lGRBs can also be divided into high-and low-luminosity classes (see Zhang 2018, for a detail).They can even be grouped into classical GRBs, X-ray-rich GRBs, and X-ray flashes (XRFs; Zhang et al. 2020b) by considering their peak energy, or into more interesting subclasses by considering whether they have extended emissions or not (Zhang et al. 2020a).Considering the joint observational features, Zhang et al. (2009) proposed to classify GRBs into two subgroups of long/soft type I and short/hard type II in physics.
Since most lGRBs are believed to have a massive star corecollapse origin, their redshift distribution is expected to trace the SFR, which has been revealed by various SFR indicators (namely Madau et al. 1998;Daigne et al. 2006;Hopkins & Beacom 2006;Grieco et al. 2012).Over the past decades, the event rates derived with different methods are still controversial and even conflicting with the star formation history.Interestingly, it was found that the lGRB rate exhibits an excess at either the low-redshift end (e.g., Petrosian et al. 2015;Yu et al. 2015;Zhang & Wang 2018;Lloyd-Ronning et al. 2019) or the high-reshift end (e.g., Kistler et al. 2008;Yüksel et al. 2008;Virgili et al. 2011;Lan et al. 2019) when it is compared with the SFR.Since the current high-z SFR is measured with a relatively large uncertainty (see, e.g., Dong et al. 2022, hereafter D22), the discrepancy between the lGRB formation rate and the SFR at the high-redshift end is not prominent yet.On the contrary, whether the lGRB rate traces the SFR in the low-z domain is still an important issue that has been deeply investigated recently (e.g., Petrosian et al. 2015;Yu et al. 2015;Pescalli et al. 2016;Dainotti et al. 2021;Palmerio & Daigne 2021;Dong et al. 2022;Nugent et al. 2022).Similarly, some authors claimed that the formation rate of sGRBs also exceeds the SFR at low redshift (Zhang & Wang 2018;Dainotti et al. 2021).Recently, Nugent et al. (2022) analyzed the stellar population of 69 sGRB host galaxies and concluded that the event rate of sGRBs traces a combination of recent star formation and stellar mass, but is not singularly dependent on either property.It is interesting to note that some authors argued that low-luminosity GRBs are a distinctive class of bursts (e.g., Sazonov et al. 2004;Cobb et al. 2006;Liang et al. 2007).However, the relation between their event rate and the low-redshift lGRB rate has not been carefully examined yet.Since most low-luminosity GRBs occur at a relatively low redshift, it is quite possible that they could be the key factor leading to the low-redshift excess of the lGRB rate.Here we will study this issue in-depth.
People have employed various statistical methods to assess the GRB formation rate and compare it with SFR.There are mainly two approaches frequently applied for this purpose.The first one is the traditional parametric estimate, namely the direct fit (Gehrels et al. 2004;Guetta & Piran 2005;Daigne et al. 2006;Le & Dermer 2007;Liang et al. 2007;Yüksel et al. 2008;Wanderman & Piran 2010, 2015;Elías & Martínez 2018) or the maximum likelihood method (Lan et al. 2019) on the stellar population analysis (Ghirlanda & Salvaterra 2022;Nugent et al. 2022), in which the dependence of luminosity on redshift was neglected and the luminosity distribution should be assumed in advance.The second one is the nonparametric τstatistic method (Lynden-Bell 1971;Efron & Petrosian 1992;Petrosian et al. 2015) that was initially proposed to investigate the luminosity function of quasars and its evolution based on a flux-limited sample.Noticeably, the evolution of luminosity with redshift can be represented by a factor of (1 + z) δ .It has now been popularly adopted to calculate the luminosity function and event rate of both lGRBs and sGRBs (Wu et al. 2012;Petrosian et al. 2015;Yu et al. 2015;Pescalli et al. 2016;Tsvetkova et al. 2017;Zhang & Wang 2018;Lloyd-Ronning et al. 2019;Dainotti et al. 2021;Dong et al. 2022).It should be emphasized that the derived event rate could be somewhat biased by the statistical method, the instrumental effect, and the sample selection.By contrast, the nonparametric method has a special advantage for a truncated sample since it does not depend on any preassumptions of the luminosity property and can give a point-by-point description of the cumulative distributions.Therefore, we will use the nonparametric c − method to investigate different samples of Swift lGRBs in this study.
To reduce the instrumental bias, GRBs detected by a single satellite are usually utilized as the first choice.As a result, most previous works have paid particular attention to the Swift/BAT (Burst Alert Telescope) GRBs.However, note that such GRB samples collected from the literature are generally incomplete according to the selection standards defined in Jakobsson et al. (2006) and Salvaterra et al. (2012).The completeness of the GRB sample was especially emphasized in Pescalli et al. (2016, hereafter P16), but our recent work showed that it does not affect the resulting event rate too much as long as the selected GRB sample has a large size and a sufficient redshift coverage (Dong et al. 2022).It strongly indicates that the low-redshift excess of the GRB rate is not due to the sample incompleteness.Furthermore, even for the same Swift lGRB sample, different methods may lead to inconsistent results.For example, P16 used the nonparametric method to study a complete sample of 81 Swift lGRBs and found no low-redshift excess, which, however, may be caused by the neglect of the differential comoving volume in their calculations.Using the same complete lGRB sample in P16, Lan et al. (2019) employed the maximum likelihood method to examine the luminosity function and event rate and concluded that the GRB event rate is roughly consistent with the SFRs at z < 2, but shows a discrepancy at a high redshift of z > 2. However, D22 also applied the nonparametric method to the P16 samples and confirmed the similar low-redshift (z < 1) excess and the highredshift consistency in both complete and noncomplete Swift lGRB samples.Meanwhile, they concluded that the resulting event rate is basically threshold-independent provided that the number of bursts above the flux limit is sufficient enough.However, whether the event rate of low-luminosity GRBs differs from that of high-luminosity ones and how these two types of GRBs trace star formation are still unknown yet.It motivates us to utilize a larger Swift GRB sample to investigate such an important issue in this study.Furthermore, we will compare the event rates derived from complete samples with those from incomplete samples.
This Letter is arranged as follows.Section 2 describes our sample selection, data reduction, and method.The detailed results are presented in Section 3. We summarize our findings in Section 4 and end with a discussion in Section 5.

Data and Method
To obtain the less-biased luminosity function and event rates of low-and high-luminosity lGRBs and their dependence on the redshift, we need to collect those lGRBs with both redshift measurement and well-constrained spectra so that the bolometric luminosity can be precisely calculated.First, we select 365 lGRBs with redshift and peak photon flux in the Swift/ BAT catalog5 from 2005 January 26 to 2022 May 31.Note that the big lGRB sample (called the BL sample thereinafter) does not suffer from any selection effects except the Malmquist bias (Dainotti & Amati 2018).The method used to convert the peak photon flux into peak energy flux (F p ) can be found in Zhang et al. (2018).Owing to the narrow energy range of the Swift/ BAT detector, a single power-law spectrum of Φ(E) ∝ E − α (see also Tang et al. 2019) has been applied to calculate the bolometric luminosity , where F p is in units of erg cm −2 s −1 , K c is the K-correction factor, and denotes the luminosity distance at a redshift of z.Here, a flat ΛCDM universe with Ω m = 0.27, Ω Λ = 0.73, and H 0 = 70 km s −1 Mpc −1 will be used (Zhang et al. 2021).Second, we also take the sample of 127 Swift lGRBs with wellmeasured spectra from Yu et al. (2015;hereafter Y15) and the complete sample of 79 lGRBs in D22 in order to check whether the resulting event rates of lGRBs are dependent on the spectral quality, the sample size, and the sample completeness.
Subsequently, we divide the above three samples of lGRBs into low-luminosity and high-luminosity ones according to the differential luminosity distributions as shown in Figure 1.In this figure, the broken luminosity is marked by a vertical line in the left panels, and the luminosity versus the redshift of low-and highluminosity lGRBs is displayed in the right panels.To put the luminosity limit on the plots, we need to set a suitable flux limit in advance.For the samples of Y15 and D22, we have taken the threshold as . For our newly built sample BL, we plot the differential flux distributions of the two kinds of lGRBs in Figure 2 the threshold in order to reserve sufficient low-luminosity lGRBs.
It is worth noting that we have determined the boundaries between low-and high-luminosities via the following broken power-law function which has been widely adopted in previous studies (e.g., Lloyd-Ronning et al. 2002;Guetta & Piran 2005;Petrosian et al. 2015).Note that the broken power-law function presents a better description for the luminosity function in Figure 1, as compared with the single power-law function.
As a result, we get the broken luminosity as L b ∼ 4.92 × 10 51 erg s −1 , 3.4 × 10 52 erg s −1 , and 1.09 × 10 52 erg s −1 for the Y15, BL, and D22 samples, respectively, which are larger than the former artificial reference value of L b ∼ 10 50 erg s −1 by about 1-2 orders of magnitude.In practice, there are no strict criteria to distinguish between low-and high-luminosity lGRBs previously.For example, Cobb et al. (2006) suggested that supernovae (SNe) SN 1998bw/GRB 980425 and SN 2003lw/GRB 031203 are subluminous in gamma rays, but they did not provide a standard for the classification (see also Liang et al. 2007).The fractions of low-luminosity lGRBs on the basis of our criterion given by Equation (1) are ∼22% for sample Y15, ∼44% for sample BL, and ∼41% for sample D22.For the sample BL, we caution that the percentages of low-redshift (0 < z < 1) lGRBs are ∼36% and ∼6.2% for the low-and high-luminosity lGRBs, respectively.The fraction of low-luminosity lGRBs goes up to ∼88% in the redshift range from 0 to 1.This demonstrates that most low-luminosity  lGRBs are located at high-redshift instead of in the nearby Universe and will be absolutely dominant at the low-redshift end.
Following our recent work, we utilize the Lynden Bell's c − method and the τ-statistic to determine the dependence of the peak luminosity on the redshift, which is assumed to take the form of L ∝ (1 + z) k (see Dong et al. 2022

GRB Rate versus SFR
Using the cumulative redshift distribution function f(z), one can calculate the lGRB rate, ρ(z), as done in D22.The normalized cumulative redshift distributions of samples Y15, BL, and D22 are shown in Figure 4, where we find that the low-luminosity lGRBs are relatively closer to us than those high-luminosity ones in each sample.In comparison, the cumulative redshift distributions of low-luminosity lGRBs increase faster in the smaller redshift domain.Interestingly, the median redshifts of both low-and high-luminosity lGRBs in sample Y15 are evidently larger than those corresponding medians of sample D22, but smaller than the medians of sample BL.More interestingly, we find that the median redshift of high-luminosity lGRBs in each sample is about 2 times larger than that of low-luminosity lGRBs, for which the physical reason is unknown yet.However, we can speculate that weak GRBs at high redshifts are most probably lowluminosity bursts not detected by us due to the Malmquist bias effect (Dainotti et al. 2021).They need more effective ways to detect in the future (Jiang et al. 2021;Tanvir et al. 2021;Thomas et al. 2022).Now, we use the cumulative redshift distribution to compute the event rates of the two classes of lGRBs in each sample.The results are compared with the SFR in Figure 5.It can be obviously found that the dependence of the event rate on redshift for the low-and high-luminosity lGRBs is very different no matter whether the lGRB sample is complete or not, which is consistent with our findings in D22.In other words, the event rates of high-luminosity lGRBs match the SFR well in all redshift ranges and the low-luminosity lGRBs exhibit significant excess at low-redshift of z < 1, particularly for the complete sample D22.Even for the smaller and incomplete Y15 sample, the low-redshift excess is still significant.On the other hand, Figure 5(b) indicates that the event rate of a small fraction of low-luminosity lGRBs with  higher redshift of z > 1 is somewhat in agreement with the SFR.This is evidently different from the previous viewpoint that the low-redshift excess was biased by the incompleteness of GRB samples (e.g., Pescalli et al. 2016;Lan et al. 2019).Strictly speaking, the sample completeness can influence the estimation of GRB rate slightly, as displayed in Figure 5, but it would not change the evolutionary trend of GRBs in nature and its importance may be overestimated in some previous studies.Our results indicate that the high-luminosity lGRBs are produced by the core collapse of massive stars, while most low-luminosity lGRBs could originate from other processes, including compact star mergers, which are unrelated to the SFR.Note that the limited GRB samples may bias the estimation of the GRB rate, especially in the low-redshift region.To test the reliability of our results, we have doublechecked the data used in Figure 5, ensuring that the number of data points in each redshift bin is no less than 5.In this way, we still find that there is excess at low redshift in both the complete and incomplete lGRB samples based on Poisson statistics.The excess is mainly contributed by the low-luminosity lGRBs as well.

Conclusions
Based on the analysis of GRB rates for both complete and incomplete lGRB samples, we draw the following conclusions: (1) The relations between peak luminosity and redshift of lGRBs are described by L ∝ (1 + z) k , with = (2) After dividing the lGRBs into low-and high-luminosity ones, we find that the event rates of high-luminosity lGRBs match the SFR well in all redshift ranges and the low-luminosity lGRBs show significant excess at lowredshift of z < 1.Therefore, we conclude that the lowredshift excess of the GRB rate results from the lowluminosity lGRBs.
(3) We find that the lGRB rate always exceeds SFR at low redshift of z < 1 no matter whether the lGRB samples are complete or not.That is to say, the low-redshift excess is intrinsically real and it cannot be attributed to the sample incompleteness.
Consequently, we argue that the high-and low-luminosity lGRBs might have different physical origins.High-luminosity bursts tracing the SFR should be produced by the death of massive stars and the low-luminosity bursts may come from other processes that are unrelated to the SFR.The results might shed new light on the investigation of GRB progenitors and central engines.

Discussion
Currently, the origin of the low-luminosity GRBs is still an open question.For example, Liang et al. (2007) and Guetta & Della Valle (2007) claimed the existence of two physically distinct classes of lGRBs but cannot exclude a single population.Lan et al. (2021) showed that the low-luminosity lGRBs are not the straightforward extension of high-luminosity ones through the luminosity distribution of their sample.Pais et al. (2023) mentioned that a cocoon breakout is most likely the origin of low-luminosity GRBs.Broad absorption lines produced in such a process have been observed in several SNe associated with low-luminosity lGRBs (e.g., Bufano et al. 2012;Xu et al. 2013;Izzo et al. 2019).Guetta & Della Valle (2007) investigated the GRB classification by considering their probability of being associated with an SN Ib/c.They found that low-luminosity GRBs are more likely to be associated with an SN Ib/c than high-luminosity GRBs.However, it is still too difficult to infer the nature of these two populations simply from their difference in the probability.Grieco et al. (2012) compared the Type Ib/c event rate with the GRB rate.It is found that the local ratio of GRB/SNe Ib/c is <3 × 10 −3 in spirals, and is <2 × 10 −2 in irregulars, which means that only a tiny fraction of SNe Ib/c sources end up exploding as lGRBs.Another possible origin involving the cocoon shock breakout from the active galactic nucleus disk was also proposed (Zhu et al. 2021). Particularly, MichałowskI et al. (2018) argued that the low-luminosity GRB 111005A represents an independent class of GRBs differing from the typical core-collapse events and inferred that the GRB rate may be biased toward low SFRs.
The discrepancy between the GRB rate and the SFR has been tried to be reconciled by several groups by considering the cosmic metallicity evolution (Salvaterra et al. 2012;Lloyd-Ronning et al. 2019) or an evolving luminosity function of GRBs (Lloyd-Ronning et al. 2002;Pescalli et al. 2016;Lan et al. 2021).Soderberg et al. (2006) derived the intrinsic rate of subenergetic GRBs as ∼230 Gpc −3 yr −1 , which is about 10 times more abundant than typical bright GRBs.The lowluminosity lGRBs are likely a unique population because of the high local event rate and small beaming factor (see, e.g., Guetta & Della Valle 2007;Liang et al. 2007;Zhu et al. 2021).However, the real reason for the event rate exceeding the SFR at low redshift remains unknown.The most important finding in this study is that the low-redshift excess of GRB rate is induced by low-luminosity events rather than faint lGRBs.Moreover, the low-luminosity lGRBs at low redshift may not be related to the death of massive stars but should be most probably connected to the merger of compact stars with a relatively low mass.
The reason that some previous works found no excess at low redshift may be that their samples mainly consisted of highluminosity lGRBs (Salvaterra et al. 2012;Wu et al. 2012;Pescalli et al. 2016;Ghirlanda & Salvaterra 2022) since they have utilized the complete GRB samples, such as the BAT6 sample (Salvaterra et al. 2012), which is dominated by bright Swift lGRBs.Dong et al. (2022) also applied the same nonparametric method to a complete sample, including weak bursts and found that low-redshift excess exists.Simultaneously, they verified that the GRB rate exceeds the SFR at low redshift even for the noncomplete GRB samples when the fraction of low-redshift GRBs is sufficient enough.Interestingly, we find in this study that the low-redshift excess is caused by the low-luminosity lGRBs instead of the highluminosity ones, which is independent of the completeness of GRB samples.
An interesting issue is that the distributions of both luminosity and redshift can be influenced by the detailed geometry of the jet.It is quite clear that a realistic GRB jet could have some kind of structure (Guo et al. 2020).Therefore, how different jet models affect the luminosity and redshift distributions is an important issue.Some authors have tried to solve the problem by assuming a particular jet profile in advance (e.g., Nakar et al. 2004;Guo et al. 2020;Banerjee & Guetta 2022).For example, Nakar et al. (2004) discussed the influence of a universal structured jet model on the redshift distribution.They also calculated the theoretical GRB rate as a function of redshift and viewing angle.It is found that the structured jet model is incompatible with the currently available observational data because of the redshift selection effects.More recently, Banerjee & Guetta (2022) proposed that the scatter of the luminosity in the comoving frame can be influenced by the angular emission profile of GRB jets.They also found that the mean luminosity is insensitive to the on-axis luminosity of the jet.But still, they did not find the actual effect of a particular structured jet model on the distributions of luminosity and redshift.In practice, it is usually difficult to discriminate whether the GRB ejecta is structured or isotropic through observations.It is even more difficult to distinguish different types of structured jets.In some cases, our line of sight may deviate from the jet axis (Xu et al. 2023), which will lead to more difficulties in deriving the luminosity function and event rate of GRBs.We hope that the rapid increase in the number of GRBs with redshift measured will bring more opportunities to us in the future.
, from which we can get the smallest flux as

Figure 1 .
Figure1.The differential luminosity distributions of samples Y15, BL, and D22, shown in Panels (a), (c), and (e), respectively.The solid thick lines denote the best fits with a broken power law of Equation (1) and the vertical thin lines mark the boundaries between the low-and high-luminosity events.Panels (b), (d), and (f) plot the luminosity vs. reshift for the Y15, BL, and D22 samples, correspondingly, together with the estimated threshold (dashed line).The observed low-and highluminosity lGRBs are shown with filled squares and circles, respectively.

Figure 2 .
Figure 2. The differential flux distributions of low (filled squares) and high (filled circles) luminosity lGRBs.The two vertical dashed lines represent the smallest fluxes of them.The two solid lines denote the homogeneous distribution of N ∝ F −3/2 in the Euclidian space.
for details).By varying the k value until τ reaches zero, we obtain = in the D22 sample for low-and high-luminosity lGRBs, correspondingly.The results are illustrated in Figure3.It can be seen that the power-law indices of k l and k h are significantly different in each sample.Interestingly, we notice that k l is larger than k h in the two samples of BL and D22, but the case is opposite for sample Y15, which may be caused by the lower proportion of low-luminosity lGRBs in the Y15 sample.Although the derived k values are sample-dependent, they are roughly coincident with previous estimates.For instance,Wu et al. (2012) derived the value of k as - 95 lGRBs detected by multiple satellites including Swift, Fermi/GBM, Konus-wind, and so on.Tsvetkova et al. (2017) got a smaller value of k ∼ 1.7 for a sample of 150 Konus-wind lGRBs with reliable redshift measurements.It indicates that the k values are less biased by the instrumental effect as compared with the sample selection effect.

Figure 3 .
Figure 3.The τ-statistic results of low-luminosity (thin solid line) and high-luminosity (thick solid line) GRBs for samples Y15 (panel (a)), BL (panel (b)), and D22 (panel (c)), respectively.The vertical lines show the expected k values when τ is equal to zero.

Figure 4 .
Figure 4. Cumulative redshift distributions of low (thin solid line) and high (thick solid line) luminosity lGRBs for the three samples.All data have been normalized to the maximum redshift.Panels (a), (b), and (c) correspond to the Y15, BL, and D22 samples, respectively.

Figure 5 .
Figure 5.Comparison between the SFR (gray points) and the event rates of low (red thin ladder lines) and high (blue thick ladder lines) luminosity lGRBs for samples Y15 in panel (a), BL in panel (b), and D22 in panel (c), respectively.The observed SFR data together with the best-fitting line (cyan solid lines) have been taken from D22.The light shaded regions indicate the 3σ range of SFR.