The CHIME Fast Radio Burst Population Does Not Track the Star Formation History of the Universe

The redshift distribution of fast radio bursts (FRBs) is not well constrained. The association of the Galactic FRB 200428 with the young magnetar SGR 1935+2154 raises the working hypothesis that FRB sources track the star formation history of the universe. The discovery of FRB 20200120E in association with a globular cluster in the nearby galaxy M81, however, casts doubts on such an assumption. We apply the Monte Carlo method developed in a previous work to test different FRB redshift distribution models against the recently released first CHIME FRB catalog in terms of their distributions in specific fluence, external dispersion measure ($\rm DM_E$), and inferred isotropic energy. Our results clearly rule out the hypothesis that all FRBs track the star formation history of the universe. The hypothesis that all FRBs track the accumulated stars throughout history describes the data better but still cannot meet both the $\rm DM_E$ and the energy criteria. The data seem to be better modeled with either a redshift distribution model invoking a significant delay with respect to star formation or a hybrid model invoking both a dominant delayed population and a subdominant star formation population. We discuss the implications of this finding for FRB source models.


INTRODUCTION
The engines that power cosmological fast radio bursts (FRBs) are not well identified. Different types of engines may follow different redshift distributions (e.g. Fig.1 of Zhang et al. (2021)), so constraints on the FRB redshift distribution would offer clues to the origins of FRBs. The discovery of the Galactic FRB 200428 in association with an X-ray burst from the magnetar SGR 1935+2154 (CHIME/FRB Collaboration et al. 2020;Bochenek et al. 2020;Li et al. 2021;Mereghetti et al. 2020;Ridnaia et al. 2021;Tavani et al. 2021) suggests that at least some FRBs originate from young magnetars produced from deaths of massive stars. As a result, the working hypothesis that the majority of FRB sources follow the star formation history of the universe has been widely adopted in the community. Studies of the host galaxy properties and FRB position offsets from the hosts suggest that the FRB population is generally consistent with such a hypothesis (Li & Zhang 2020;Bhandari et al. 2020;Heintz et al. 2020;Bochenek et al. 2021;Mannings et al. 2021;Fong et al. 2021), even though the alternative hypothesis that FRBs follow a stellar population with a significant delay with respect to star formation (e.g. the distribution that track binary neutron star mergers) is not ruled out (Li & Zhang 2020). Recently, a repeating source FRB 2020120E was discovered to be located in a globular cluster of the nearby galaxy M81 (Bhardwaj et al. 2021;Kirsten et al. 2021;Nimmo et al. 2021), suggesting that at least some FRBs are associated with old stellar populations. It is desirable to know whether such a delayed population makes up a significant fraction of the entire FRB population.
So far, the FRB redshift distribution is poorly constrained. Only more than a dozen FRBs have direct redshift measurements (Tendulkar et al. 2017;Bannister et al. 2019;Ravi et al. 2019;Marcote et al. 2020;Macquart et al. 2020;Bhandari et al. 2021). These bursts were detected with a variety of radio telescopes that have very different instrumental selection effects, so this limited sample cannot give a reliable constraint on FRB redshift distribution. The dispersion measure (DM) of an FRB can be a rough proxy of its redshift, as is verified by the observational confirmation (Macquart et al. 2020) of the DM − z relation long suggested theoretically (Ioka 2003;Inoue 2004;Deng & Zhang 2014). In principle, based on the observed fluence distribution, DM distribution, and inferred energy distribution of FRBs, one can constrain the redshift distribution of FRBs (Zhang et al. 2021;James et al. 2022). However, the situation is so far inconclusive due to the small FRB samples and complicated observational selection effects. In particular, the small Parkes and ASKAP FRB samples are not inconsistent with either a model tracking the star formation history of the universe (Zhang et al. 2021;James et al. 2022) or models invoking a significant delay with respect to star formation (e.g. those related to binary neutron star mergers) (Zhang et al. 2021).
Recently, the CHIME/FRB Collaboration published their first FRB catalog (The CHIME/FRB Collaboration et al. 2021) reporting 536 FRBs including 62 bursts detected from 18 repeating sources. This uniform large sample provides an ideal resource to constrain the FRB redshift distribution. Chawla et al. (2021) performed a population study of the CHIME catalog. They focused on DM and scattering distributions and constrained the properties of the circumgalactic medium. They only assumed that the FRB rate evolves with redshift and tracks the star formation history of the universe without testing a range of redshift distribution models.
In this Letter, following the Monte Carlo method described in our earlier work (Zhang et al. 2021), we systematically investigate the consistency of various redshift distribution models with the CHIME FRB catalog. We rule out the hypothesis that the entire FRB population tracks the star formation history of the universe with high significance. We also reject the hypothesis that the FRB population tracks the accumulated stars throughout the history of the universe, despite appearing to match the data better. The data seem to instead approach a model that either requires the FRB sources to have a significant delay with respect to star formation or a hybrid model that includes both a dominant delayed population and a subdominant star formation population. Our method is reviewed in Section 2. The results are presented in Section 3, and conclusions are drawn in Section 4 with some discussion.

THE METHOD
The details of our Monte Carlo method have been described in Zhang et al. (2021). Here, we only outline the key ingredients of the method. Basically, for a uniform FRB sample detected by the same telescope (as is the case of the CHIME sample), the fluence, distance (and hence, DM E ), and energy distributions of the observed FRB population depend on three factors: (1) the intrinsic redshift distribution, (2) the intrinsic FRB isotropic energy (or luminosity) distribution, and (3) the telescope's sensitivity threshold and instrumental selection effects near the threshold.
The first factor is the focus of investigation of this paper. Interestingly, the other two factors are largely decoupled from the first factor and from each other, which allows us to treat them independently. The second factor (energy distribution) has been well constrained from observations. Independent studies regardless of the assumed redshift distribution (e.g. Luo et al. 2018Luo et al. , 2020Lu & Piro 2019;Lu et al. 2020;Zhang et al. 2021) have reached the consistent conclusion that the energy distribution of the entire FRB population is roughly a power law dN/dE ∝ E −α covering at least 8 orders of magnitude, with the index α ∼ (1.8 − 2.0). There might be a high-energy exponential cutoff (Luo et al. 2020;Lu et al. 2020) but the cutoff energy E c is not well constrained (Zhang et al. 2021) 1 . In principle, there could be a redshift evolution of FRB energy distribution. However, current data do not require such an evolution. Furthermore, the mechanism of FRB sources (e.g. magnetars) is likely related to their own physical properties rather than redshift. As a result, we have assumed a universal energy function for all FRBs throughout the universe. The third factor (instrumental effects) is difficult to characterize. The CHIME catalog data show that the telescope specific fluence cutoff is about 0.3 Jy ms, or log F min −0.5 (see panel (a) of Figs.2-5). However, due to many instrumental or human-related uncertainties of CHIME observations (e.g. unknown positions of most bursts that introduce large errors in the estimated fluences, nonuniform CHIME sensitivity on the sky due to large gaps between beams, and DM dependence and scattering time dependence of the detection efficiency The CHIME/FRB Collaboration et al. 2021), there is a "gray zone" in the F distribution within which the CHIME telescope has not reached full sensitivity to all sources. These effects are independent of the redshift distribution and may be corrected using an empirical model (see below) 2 .
By adopting an intrinsic energy (E) distribution and a redshift (z) distribution, one can simulate a large num-1 Even though the energy distribution was constrained from the observations using the telescopes (e.g. Parkes, ASKAP in ∼ GHz) with observing bands higher than that of CHIME (400 to 800 MHz), it is reasonable to assume that the shape of the E distribution in the CHIME band is similar to that constrained in the GHz band. Throughout the paper, our E is defined as the isotropic energy as observed in the CHIME frequency band.
2 Effectively, this is to assume that the efficiency for CHIME to detect FRBs is essentially independent of the DM value of the FRB, especially around the DM E distribution peak region ∼ (250 − 500) pc cm −2 . This is consistent with the analysis by the CHIME team, who stated that "the selection effects in DM are modest" and that "we appear to be detecting the full range of DMs represented in the population detectable at CHIME/FRB's sensitivity" (The CHIME/FRB Collaboration et al. 2021).
ber of mock FRBs. The specific fluence of each mock burst can be calculated based on its assigned E and z values. After screening them using a telescope sensitivity model, one can finally obtain a mock "observed" sample of FRBs. Based on the DM − z relation (Deng & Zhang 2014;Zhang 2018a;Pol et al. 2019;Cordes et al. 2021) , Ω m , Ω λ are cosmological parameters whose values are adopted from the latest Planck results (Planck Collaboration et al. 2016), G is the gravitational constant, m p is proton mass, and f IGM is the fraction of baryons in the IGM, which is adopted as 0.84), one can estimate the IGM portion of the dispersion measure (DM IGM ) of each mock FRB. Adopting a model for the host galaxy dispersion measure (DM host ), one can finally simulate the excess DM distribution of the mock sample, which can be compared with the DM E = DM IGM + DM host /(1 + z) data directly retrievable from the CHIME catalog (The CHIME/FRB Collaboration et al. 2021) 3 .
In order to test a certain z-distribution model, we make use of three observational criteria (see Figs.2-5): (1) the specific fluence (log F ν ) distribution, (2) the isotropic energy (log E) distribution, and (3) the excess dispersion measure (DM E ) distribution. The specific fluence is a convolution of the energy and redshift distributions and is insensitive to either distribution. 4 This is because the fluence distribution (also called log N − log F ν distribution) should follow a simple N ∝ F −3/2 ν distribution regardless of the energy function if the sources are uniformly distributed in a Euclidean space. 5 The non-Euclidean geometry of cosmological models may break the simple scaling, but only in the low fluence regime. However, in this low fluence regime, the instrumental selection effects become so important that any redshift-distribution-related features are removed. As a result, the criterion 1 is most easily satisfied among all three criteria for all models. For every pair of z and E distribution models, it is possible to find an empirical instrumental selection effect model to satisfy the log F ν criterion. However, many models that satisfy the log F ν distribution criterion could fail the log E and DM E criteria. We therefore come up with the following strategy: for each z-distribution model, we adjust the E distribution model and the empirical sensitivity model to make the F ν distribution of the mock "observed" sample not be rejected by the Kolmogorov-Smirnov (K-S) test (all the K-S test statistics in this paper are reported with 95% confidence) against the observed F ν distribution. We then go on to evaluate the log E and DM E distribution criteria. The model is ruled out if the same mock FRB sample fails both criteria.
For instrumental sensitivity threshold modeling in the "gray zone" of the lower end of the specific fluence distribution, we adopt log F min ν,th = −0.5 as the minimum threshold specific fluence (as shown by the data). We define a maximum threshold specific fluence, log F max ν,th , whose value is adjusted to match the observation, and define a detection efficiency parameter, η det , that depends on the ratio R = (log F ν,th − log F min ν,th )/(log F max ν,th − log F min ν,th ) for fluences between log F min ν,th and log F max ν,th , such that η det → 0 at log F min ν,th and η det → 1 at log F max ν,th . We model the dependence as η det = R n and find n = 3 can model the log F ν distribution and adopt this empirical function form and the typical value in our modeling.
Strictly speaking, a more rigorous but cumbersome method would be to vary all the free parameters (for the intrinsic energy distribution model and the empirical sensitivity model) for each redshift distribution model to compare against all three observational criteria. However, a salient feature is that both the fluence distribution and the instrumental effects near the fluence sensitivity threshold essentially do not depend on the redshift distribution model of FRBs. This allows us to adopt the abovementioned simpler approach to test all the zdistribution models. Nonetheless, for the star formation history model, which is the main model of interest, we also apply the more rigorous method to confirm the validity of the simpler method. We test the range (1.8, 2.0) inclusive with step size 0.1 for energy powerlaw index α, (41.5, 43.5) inclusive with step size 1 for log E c , (0.65, 0.85) inclusive with step size 0.1 for F max ν,th , and (1, 3) inclusive with step size 1 for the index n in the instrumental empirical function. We find a total of 15 models that are not rejected by the 1D K-S test for the F ν distribution. All 15 of those models, however, are rejected by the 1D K-S test for both the energy and DM E distributions. Thus, as expected, even with the best sets of parameters, the model fails to reproduce the data, so our main conclusion of the paper, i.e. FRBs do not track the star formation history of the universe, is robust. Since our purpose is not to find the best parameters to meet the data, and since there are more parameters in the delayed and hybrid models, we do not perform a multi-dimensional parameter search for those models.
In our simulations, we use the central value of DM host , 107 pc cm −3 , as constrained from the data (e.g. . In principle, both the DM host and DM IGM values for individual FRBs should have a distribution around the central values. However, since we are simulating a large sample of mock FRBs, the resulting distributions of various parameters by adopting central values should be similar to the more realistic case involving these distributions. In particular, our treatment of DM host is valid if its true distribution is normal. The true DM host distribution is, however, difficult to constrain from current observations and may have a more complicated form. For example, some FRBs have very large DM host values (Niu et al. 2021). Nonetheless, a larger DM host would imply a smaller DM IGM and a smaller z. If a larger DM host value than assumed is prevailing among FRBs, the true FRB z distribution would have a lower peak than what is inferred in the following discussion. This only strengthens the conclusion drawn in Section 3.

RESULTS
We investigate four families of z-distribution models in detail: 1. Star formation rate history (SFH) model: FRB sources follow the star formation history of the universe; 2. Accumulated model: FRB sources follow the number of total stars in the universe; 3. Delayed model: FRB sources follow stellar populations that have a significant delay with respect to star formation; 4. Hybrid model: a fraction of FRB sources follows star formation and another fraction follows a delayed population.
For each model, we first build a model of their redshiftdependent event-rate density dN/(dV dt). We then convert it to a redshift rate distribution using (e.g. Sun et al.

2015)
where and D L is the luminosity distance at z. Since we only care about the z distribution, we use the normalized probability distribution functions (PDFs) of dN/(dtdV ) and dN/(dt obs dz) for all of the models, which are presented in the upper and lower panels, respectively, of Figure 1. For the data, we take all nonrepeaters and only the first detected burst for each repeating FRB from the CHIME catalog.
For each model, we provide example parameters we use to get the respective fluence distributions not rejected by the K-S test in Table 1. We note that the SFH model requires a higher E c than the other models because the SFH model predicts on average higher redshift FRBs, which requires higher energies to get the same fluence distribution.

Star formation history (SFH) model
This is the most well-motivated model. Galactic magnetars are usually believed to be young neutron stars born from supernova explosions, some of which are found to be associated with young star clusters or supernova remnants (Kaspi & Beloborodov 2017). If magnetars are the sources of most FRBs in the universe (e.g. Li & Zhang 2020;Bhandari et al. 2020;Heintz et al. 2020;Bochenek et al. 2021), one would then naturally expect that FRBs follow the star formation history of the universe (Zhang et al. 2021;James et al. 2022;The CHIME/FRB Collaboration et al. 2021).
To test this model, we use the analytical threesegment empirical model of Yüksel et al. (2008). This model is consistent with the widely used two-segment empirical model of Madau & Dickinson (2014) but more precisely catches the SFH at high redshifts mapped by long gamma-ray burst observations. We assume that dN/(dtdV ) of FRBs is proportional to the volumetric star formation rate and derive its redshift distribution dN/(dt obs dV ) (purple curves in Fig.1). Even though this model was found consistent with the smaller Parkes and ASKAP FRB samples (Zhang et al. 2021;James et al. 2022), Figure 2 shows that it fails to account for the CHIME data. The model overpredicts the number of FRBs at relatively high redshifts, and  hence, high DM E (Fig.2c) 6 , which requires an E distribution that has a higher peak than observed (Fig.2b). The K-S tests for these two criteria show very clearly that the star formation history model is rejected to describe the data. In Fig.2d, we show the two-dimensional distribution of the data and simulated mock FRBs in 6 Note that even though the fluence and energy distributions are plotted in the logarithmic scale, the DM E distribution, which is a proxy of the z distribution, is plotted in the linear scale. This is because unlike Fν and E which have power-law distributions, none of the physical parameter scales with z as a power law. A more relevant parameter to delineate power-law behavior in cosmology is (1 + z), which is inversely proportional to the scale size of the universe. This parameter only changes by a factor of 2 from z = 0 to z = 1. For the redshift range we are interested in, it is more reasonable to present the DM E distribution in the linear space.
the DM E − log E space. A deficit of low DM E , low log E FRBs is clearly seen.
Since the SFH model is the most speculated and discussed model, we also perform two additional tests. First, even though both the 1D K-S tests for DM E and E criteria have rejected the model, we still perform a 2D K-S test for fluence and DM E to see whether the model has any chance of survival. We find that the null hypothesis that the data and star formation rate history simulation are from the same inherent distribution is rejected. Second, we also perform the more rigorous method as discussed in Section 2 to explore the compliance of the model with the data from multidimensional parameter space. Our results show that except for a small parameter space where K-S tests can pass the fluence distribution criterion but fail both DM E and E criteria, in the majority of the parameter space, all three criteria fail. These two additional tests strengthen our conclusion that the hypothesis that the CHIME FRB population track the star formation history of the universe is ruled out by the data with high significance. The main reason that this model fails to meet the data is that the observed dN/dDM E distribution drops significantly above ∼ 250 pc cm −2 , which corresponds to z ∼ 0.3. This is in sharp contrast to the prediction of the SFH model (see Fig. 1). This demands producing many more FRBs in the nearby universe than what the SFH model predicts, which requires various delayed models as discussed below.

Accumulated model
Next, we test a z-distribution model that tracks the accumulated stars throughout history (green curves in Fig.1). Hashimoto et al. (2020) advocated a model with nearly constant FRB rate during the past 10 Gyr of look-back time. They interpreted this as FRBs tracking the accumulated stellar population. We therefore test such an accumulated SFH model in detail. As shown in Fig.3, this model describes the data better than the SFH model. However, given a fluence distribution that is not rejected by the data, the resulting energy and DM E distributions are still both rejected by the K-S test when tested against the data. Thus, we still reject this model to describe the CHIME FRB population. Nonetheless, since this model is a better representation of the data, it suggests that a large fraction of FRBs likely are produced by a stellar population that is significantly delayed with respect to star formation.

Delayed models
We next explore a family of models that invoke a significant delay from star formation. As described in detail in Zhang et al. (2021), we first generate dN/(dtdV ) based on the SFH model and calculate the look-back time t L distribution of the sample. Next, we introduce a distribution model for the delay time τ , e.g. in the form of power-law, Gaussian, or lognormal functions (e.g. Virgili et al. 2011;Sun et al. 2015;Wanderman & Piran 2015). We then subtract the look-back time of the SFH model by the delay time τ for each mock FRB and finally obtain the look-back time distribution of the delayed population. The cases of negative t L are dropped out since they stand for future events. Finally, we convert the new t L distribution to the dN/(dtdV ) distribution of the new model, and simulate their redshift distribution dN/(dt obs dz) using Eq.(2). The same approach is then applied to test the three criteria in log F ν , log E, and DM E distributions.
Even though a family of delayed models was not rejected by the K-S test with the Parkes and ASKAP data as shown in Zhang et al. (2021), we find that all of the models previously tested (with a characteristic delay timescale of ∼ (2 − 3) Gyr and consistency with the short GRB data) are rejected by the K-S test with the CHIME data sample. Rather, the data require a model with a much longer delay. Figure 4 (also red curves in Fig.1) presents an example of such a lognormal delay model with a central value of 10 Gyr and a standard deviation of 0.8 dex. It reproduces the data much better than the SFH model, with both the fluence and energy criteria being not rejected by the K-S test. Even though the DM E criterion is still rejected by the K-S test, the K-S statistic is much closer to the critical value (to not reject the null hypothesis) than the SFH model. We note that there is a wide range of the inferred DM host in the FRB data, which would cause more complicated features in the modeled DM E distribution than our simple model. This could partially account for the DM E discrepancy between our model and the data. Since there are many parameters in the model and since our goal is to test the general trend of the models rather than identify model parameters, we do not make efforts to search for the best parameter set to satisfy the observational constraints. In any case, we can conclude that such a significantly delayed model offers a better description to the CHIME data than the SFH model.

Hybrid model
Finally, we test a family of models that includes the mixture of a young component that tracks star formation history and an old component that has a significant delay with respect to star formation (orange curves in Fig.1). Such a model is motivated by the discoveries of the Galactic FRB 200428 (CHIME/FRB Collaboration et al. 2020;Bochenek et al. 2020), which tracks star formation, and the M81 globular cluster FRB 20200120E (Bhardwaj et al. 2021;Kirsten et al. 2021;Nimmo et al. 2021), which tracks an old population. Since the SFH model fails badly, the proportion of the young population cannot be high. To compensate the high-DM events predicted from the young population, the old population in the hybrid model needs to have an even longer delay from star formation. We test a model with a lognormal delay distribution (central value = 13 Gyr, standard deviation = 0.8 dex) of the old population that is mixed with the SFH model of the young population. The proportions are 80% for the delayed component and 20% for the SFH component. The results are shown in Figure  5. We see similar results to that of the delayed model, with the fluence and energy criteria being not rejected, but DM E still rejected. It is clear that this model is also a much closer model to describe the data than the SFH Figure 2. A test of the SFH model against the three observational criteria. The full FRB sample is tested against the model predictions. (a) The log N>F ν − log Fν test, (b) the log E distribution test, (c) the DME distribution test, and (d) the 2D DME − log E distribution for illustrative purpose (not for testing). For panels (a), (b), and (c), the SFH simulations are scaled to the full data. For the energy distribution model, we adopt α = 1.9 and log Ec = 41.5.
model. Again, since there are even more parameters in this model than the delayed model and since our goal is to test the general trend of the models, we do not explore the best parameter set to fit the data. In any case, we deem that some hybrid models with a dominant delayed population component offer a better description to the CHIME data than the SFH model.

CONCLUSIONS AND DISCUSSION
Applying the Monte Carlo method developed in Zhang et al. (2021), we have systematically tested a variety of FRB redshift distribution models against the first CHIME FRB catalog data. We draw the robust conclusion that the CHIME FRB population do not track the star formation history of the universe. The hypothesis that FRBs track all the stars formed in the universe (the accumulated SFH) model, despite describing the data better, is also rejected. Instead, the models that invoke a significant delay with respect to star formation seem to be working toward the correct direction in describing the data. A model that FRBs track a population of sources that have a significant delay (∼ 10 Gyr) with respect to star formation is one possibility. A hybrid model of young and old components is another possibility, but the old population needs to be the dominant component. Even though it is difficult to pin down the parameters due to many unconstrained parameters, this general trend is quite robust. Our conclusion is in general consistent with several previous suggestions that FRBs do not track star formation (e.g. Cao et al. 2017;Hashimoto et al. 2020;Safarzadeh et al. 2020), Figure 3. Similar to Figure 2, but for a test of an accumulated star formation model. All the simulations are scaled to the data. For the energy distribution model, we adopt α = 1.9 and log Ec = 41.
even though an even longer characteristic delay time scale is needed to match the data better. Note that the reason that the delayed models are preferred is not because of their larger number of model parameters, but is rather because of the smaller peak value in the DM E distribution (peaking around 250 pc cm −3 ) of the CHIME sample compared with previous samples. Our results suggest that the "most conservative" scenario of FRB origin (Zhang 2020), i.e. magnetars can make them all, may have to be abandoned. Studies to interpret the M81 globular cluster FRB 20200120E Kremer et al. 2021;Lu et al. 2021) still invoke magnetars formed from other channels other than massive star deaths, e.g. binary white dwarf mergers, accretion induced collapses, and even binary neutron star mergers, to interpret FRB 20200120E. Since observationally, most known magnetars are formed from recent supernova explosions (Kaspi & Beloborodov 2017), one would expect that a small fraction of FRBs follow a delayed population if the magnetar hypothesis is correct. This is not what is inferred from the CHIME data, which requires that the delayed component is the dominant population. Challenges are raised to theorists regarding how to make FRBs with very old stars. One possibility is that FRBs are produced by reactivation of very old neutron stars across the universe, and old magnetars may fall into such a category (e.g. Wadiasingh et al. 2020;Beniamini et al. 2020). Other possibilities, e.g. FRBs being powered by interactions with old neutron stars (Mottez & Zarka 2014;Dai et al. 2016;Zhang 2017), may also deserve reinvestigation. (c) (d) Figure 5. Similar to Figure 2, but for a test of a hybrid model that includes a young component tracking the SFH (20% of the total population) and an old component (80% of the total population) with lognormal delay with a central value of 13 Gyr and standard deviation of 0.8 dex. All the simulations are scaled to the data. For the energy distribution model, we adopt α = 1.9 and log Ec = 41.