Implications from ASKAP Fast Radio Burst Statistics

TL;DR

This paper analyzes ASKAP fast radio burst data to understand their energy distribution and rate, using simplified models and Bayesian analysis, providing constraints on their energetics and occurrence rates.

Contribution

It introduces a combined analytical and Bayesian framework to infer FRB energetics and rates from fluence and DM distributions, applicable to future larger samples.

Findings

Power-law index gamma ~ 1.7 fits the data

Maximum energy E_max ~ 10^{33} erg/Hz constrained

Estimated local volumetric rate log N > 10^{32} erg/Hz ~ 2.6 Gpc^{-3} yr^{-1}

Abstract

Although there has recently been tremendous progress in studies of fast radio bursts (FRBs), the nature of their progenitors remains a mystery. We study the fluence and dispersion measure (DM) distributions of the ASKAP sample to better understand their energetics and statistics. We first consider a simplified model of a power-law volumetric rate per unit isotropic energy dN/dE ~ E^{-gamma} with a maximum energy E_max in a uniform Euclidean Universe. This provides analytic insights for what can be learnt from these distributions. We find that the observed cumulative DM distribution scales as N(>DM) ~ DM^{5-2*gamma} (for gamma > 1) until a maximum value DM_max above which bursts near E_max fall below the fluence threshold of a given telescope. Comparing this model with the observed fluence and DM distributions, we find a reasonable fit for gamma ~ 1.7 and E_max ~ 10^{33} erg/Hz. We then carry out a full Bayesian analysis based on a Schechter rate function with cosmological factor. We find roughly consistent results with our analytical approach, although with large errors on the inferred parameters due to the small sample size. The power-law index and the maximum energy are constrained to be gamma = 1.6 +/- 0.3 and log(E_max) [erg/Hz] = 34.1 +1.1 -0.7 (68% confidence), respectively. From the survey exposure time, we further infer a cumulative local volumetric rate of log N(E > 10^{32} erg/Hz) [Gpc^{-3} yr^{-1}] = 2.6 +/- 0.4 (68% confidence). The methods presented here will be useful for the much larger FRB samples expected in the near future to study their distributions, energetics, and rates.

Figures

Equations

$$dN/dE=A{E}^{-\gamma },$$

$$N(<D)={\int }_0^{D}4\pi D^{2}dD{\int }_{4\pi D^2{F}_{\mathrm{th}}}^{E_{\max }}\frac{dN}{dE}dE$$

$$=\frac{A{E}_{\max }^{1-\gamma }}{\gamma -1}{\int }_0^{D}4\pi D^{2}dD\left[{\left(\frac{4\pi D^2{F}_{\mathrm{th}}}{E_{\max }}\right)}^{1-\gamma }-1\right].$$

$$N(<D)\propto \{\begin{array}{l}{D}^3,\mathrm{for}\gamma <1,\\ D^{5-2\gamma },\mathrm{for}1<\gamma <5/2,\end{array}$$

$$D_{\max }=(E_{\max }/4\pi F_{\mathrm{th}}{)}^{1/2}.$$

$$N(<D)$$

$$N_{\mathrm{tot}}=2A{E}_{\max }^{1-\gamma }4\pi D_{\max }^3/[3(5-2\gamma )].$$

$$N(>F)={\int }_0^{\sqrt{\frac{E_{\max }}{4\pi F}}}4\pi D^{2}dD{\int }_{4\pi D^{2}F}^{E_{\max }}\frac{dN}{dE}dE.$$

$$N(>F)=N_{\mathrm{tot}}{\left(F/F_{\mathrm{th}}\right)}^{-3/2}.$$

$$D\simeq 6.7\mathrm{Gpc}\left(\mathrm{DM}/855\mathrm{pc}\mathrm{c}{\mathrm{m}}^{-3}\right),$$

$$E_{\max }\simeq 2\times 10^{33}\mathrm{erg}\mathrm{H}{\mathrm{z}}^{-1}\frac{F_{\mathrm{th}}}{50\mathrm{Jy}\mathrm{ms}}{\left(\frac{\mathrm{D}{\mathrm{M}}_{\max }}{750\mathrm{pc}\mathrm{c}{\mathrm{m}}^{-3}}\right)}^2.$$

$$A{E}_{\max }^{1-\gamma }\simeq 90\mathrm{Gp}{\mathrm{c}}^{-3}\mathrm{y}{\mathrm{r}}^{-1}.$$

$$N(>E)\simeq \frac{A{E}^{1-\gamma }}{\gamma -1}\simeq 1.1\times 10^3{E}_{32}^{-0.7}\mathrm{Gp}{\mathrm{c}}^{-3}\mathrm{y}{\mathrm{r}}^{-1}.$$

$$\frac{dN}{dE}=\frac{{\varphi }_0(1+z{)}^{\beta }}{E_{\max }}{\left(\frac{E}{E_{\max }}\right)}^{-\gamma }\exp \left(\frac{-E}{E_{\max }}\right),$$

$$N(\le z,\ge F)={\varphi }_0{T}_{\mathrm{eff}}{\int }_0^{z}dz\frac{dV}{dz}(1+z{)}^{\beta -1}{\int }_{x_{\min }}^{\infty }{x}^{-\gamma }{\mathrm{e}}^{-x}dx,$$

$$f(\mathbf{p})\propto L(\mathbf{D}|\mathbf{p})f_0(\mathbf{p}),$$

$$L(\mathbf{D}|\mathbf{p})=\underset{1\le i\le N_{\mathrm{tot}}}{\min }\left[P(k_1^{(i)},{\lambda }_1^{(i)}),P(k_2^{(i)},{\lambda }_2^{(i)})\right],$$

$$S/N\simeq \frac{F}{w_{\mathrm{obs}}}\frac{\sqrt{2B{w}_{\mathrm{obs}}}}{S_{\mathrm{sys}}},$$

$$t_{\mathrm{DM}}=3.6\mathrm{ms}\frac{\mathrm{D}{\mathrm{M}}_{\mathrm{tot}}}{10^{3}pcc{m}^{-3}}\frac{\Delta \nu }{1\mathrm{MHz}}{\left(\frac{\nu }{1.32\mathrm{GHz}}\right)}^{-3},$$

Full text

††slugcomment: Submitted for publication in The Astrophysical Journal

1. Introduction

Fast radio bursts (FRBs, Lorimer et al., 2007; Thornton et al., 2013) are millisecond radio pulses with large dispersion measures (DMs) strongly suggesting an extragalactic origin. This has been directly confirmed by the repeater FRB 121102 (Spitler et al., 2016a), which has been localized to a $z=0.19$ galaxy (Chatterjee et al., 2017; Tendulkar et al., 2017). FRBs promise to provide amazing probes of the baryonic distribution across cosmological distances, but before their full potential can be reached, a fundamental understanding of their sources is required (for an overview of the observations and potential progenitor models, see Platts et al., 2018; Petroff et al., 2019a)111http://frbtheorycat.org. A chief issue is the connection between repeating FRBs and those seemingly one-off bursts. Do all FRBs have the same source as the repeaters? What is the volumetric rate of FRBs as a function of their energies? And how do their nature and rate impact their utility as cosmological probes?

These questions have made discovering FRBs and measuring their physical properties some of the leading scientific goals of many current and future telescopes, such as Parkes (Thornton et al., 2013; Champion et al., 2016; Bhandari et al., 2018), Arecibo (Spitler et al., 2016a), UTMOST (Bailes et al., 2017; Caleb et al., 2017), ASKAP (Bannister et al., 2017; Shannon et al., 2018), CHIME (CHIME/FRB Collaboration et al., 2018), FAST (Li et al., 2013), and Apertif (Maan & van Leeuwen, 2017). This has led to a rapidly growing but highly heterogeneous sample of FRBs. In Figure 1, we summarize the DMs and fluences of currently published FRBs (Petroff et al., 2016)222http://frbcat.org. The DMs are obtained by subtracting from the total measured values the contributions of the interstellar medium (Cordes & Lazio, 2002) and halo ($\mathrm{DM}_{\rm halo}=30\rm\,pc\,cm^{-3}$) of the Milky Way.

Many studies have been carried out to statistically constrain the volumetric rate and luminosity/energy distribution function of the growing sample of FRBs, which may provide important clues about their progenitors (Katz, 2016; Lu & Kumar, 2016; Li et al., 2017; Nicholl et al., 2017; Macquart & Ekers, 2018; Luo et al., 2018; Bhattacharya et al., 2019). Such work has mainly focused on the Parkes sample. Although Parkes accounts for almost half of the bursts, there are many selection effects that make it fluence-incomplete (Keane & Petroff, 2015; Patel et al., 2018), e.g., different backend instruments (AFB vs. BPSR) and galactic-latitude rate dependence (Petroff et al., 2014; Burke-Spolaor & Bannister, 2014). Likely the most serious issue is that the fluences of Parkes bursts are often lower limits (by assuming they occurred at the beam center) and are uncertain by up to an order of magnitude due to poor localization by single-beam detection. Three Parkes bursts (010724, 110214, 150807, shown in Figure 1 with errorbars) had multi-beam detections and their fluences included the correction of antenna attenuation (Ravi et al., 2016; Petroff et al., 2019b; Ravi, 2019).

Recently, the Commensal Real-time ASKAP Fast Transient (CRAFT) survey, targeting the brightest portion of the FRB population, has published 23 bursts — the first and largest well-controlled sample to date. These bursts have well-measured fluences (uncertainty333One exception is FRB 170110 with a fluence of $200^{+500}_{-100}\rm\,Jy\,ms$ (90% confidence). This large uncertainty is due to detection in a corner beam with poor localization (Shannon et al., 2018). Nevertheless, it is well above the ASKAP fluence threshold and is hence included in our analysis. of $\sim$20%) because the overlapping beam arrangement allows the full focal plane to be uniformly sampled (Shannon et al., 2018). The survey had high galactic latitude pointings $|b|=50\pm 5$ deg, which eliminates potential biases due to varying sky temperature and Galactic DM contribution. For the above reasons, we focus on the relatively uniform ASKAP sample (Bannister et al., 2017; Shannon et al., 2018; Macquart et al., 2019). which can be shown to be roughly complete above a threshold fluence of \sim$$50\,{\rm Jy\,ms}.

This paper is organized as follows. We start with a simplified model to explain the observed DM distribution and the fluence distribution in §2. Then, we discuss the inferred model parameters and implications in §3. A numerical Bayesian analysis with a more sophisticated model is presented in §4 to compare with our simpler analytic approach. Possible sample completeness issues are discussed in §5. A few other potential caveats to keep in mind for future work are discussed in §6. A summary is provided in §7.

2. DM and fluence Distribution

In our sample, bursts with similar DM (or distance) show a large spread in fluences, which motivates us to consider a power-law volumetric rate of FRB events per unit (isotropic) energy

[TABLE]

with a maximum energy of $E_{\rm max}$ above which there are no FRBs.

Law et al. (2017) found a power-law index of $\gamma\sim 1.7$ from the statistics of the first repeater FRB 121102 (see also Wang & Yu, 2017). This is similar to the energy distribution of high-energy bursts from magnetars (Turolla et al., 2015), which perhaps lends support for models which argue for magnetars as the progenitors of FRBs (Popov & Postnov, 2010; Lyubarsky, 2014; Kulkarni et al., 2014; Pen & Connor, 2015; Katz, 2016; Kumar et al., 2017; Metzger et al., 2017; Beloborodov, 2017). Whether or not this is the case, such indices are natural consequences of self-organized critical phenomena (Katz, 1986; Bak et al., 1987; Aschwanden et al., 2016). More recently, Gourdji et al. (2019) reported 41 bursts from two hours of Arecibo observations444For comparison, earlier Arecibo searches at 1.4 GHz found 16 bursts in about 30 hours of observations (Spitler et al., 2016b; Scholz et al., 2016, 2017). , most of which are faint ($E\sim 10^{29}\rm\,erg\,Hz^{-1}$), narrow-band (\sim$$200\,MHz), and low signal-to-noise ones found by careful visual selection against RFIs. The authors obtained a much steeper power-law index $\gamma\sim 2.8$ from this sample. It is possible that the inferred $\gamma$ is affected by the selection biases in different analyses and the non-Poissonian occurrence pattern (Oppermann et al., 2018). In this paper, we focus on the energy distribution of all FRBs (including repeaters and non-repeaters) and do not impose a prior on $\gamma$ based on the knowledge of FRB 121102.

Setting $F_{\rm th}$ as the threshold fluence for detection, the total number of events observed below a distance $D$ is

[TABLE]

Since bursts in our sample have low redshifts, we ignored the cosmological effects and changes of the rate normalization $A$ with distance (these will be included in the numerical analysis in §4).

For sufficiently nearby bursts such that $4\pi D^{2}F_{\rm th}\ll E_{\rm max}$, the scaling of the cumulative number is

[TABLE]

where the first case comes from the fact that most bursts are near $E_{\rm max}$ (like standard candles) and the second case reflects the effect of a wide distribution of burst energies.

For sufficiently large $D$, the rise of $N(<D)$ becomes shallower, because the maximum energy $E_{\rm max}$ limits the number of high-redshift bursts. The maximum distance a burst may be detected is given by

[TABLE]

The integral in Equation (2) gives

[TABLE]

where the normalization $N_{\rm tot}=N(<D_{\rm max})$ is

[TABLE]

Since at low redshifts we roughly have $D\propto z\propto{\rm DM}$, Equation (5) can easily be generalized as a function of $z$ or ${\rm DM}$. In Figure 2, we compare the expression in Equation (5) to our ASKAP sample. We assume that the ASKAP sample is roughly complete above fluence $F_{\rm th}=50\rm\,Jy\,ms$ (see §5), which encompasses 20 out of 23 bursts in the ASKAP sample (thus we take $N_{\rm tot}=20$). An energy distribution power-law index of $\gamma=1.7$ gives $N(<{\rm DM})\propto{\rm DM}^{1.6}$ in low-DM end.

In the future when the energy distributions of individual repeaters are better measured, we can compare the power-law index $\gamma$ obtained from the non-repeating sample with that of the repeaters. This provides a constraint on whether they belong to the same population. We also include a comparison to $N(<{\rm DM})\propto{\rm DM}^{3}$ (short-dashed blue line), which is appropriate if most bursts have characteristic energy. This is inconsistent with the ASKAP distribution, as mentioned by Li et al. (2019).

The total number of events above fluence $F$ is

[TABLE]

Evaluating this integral and substituting the expression for $D_{\rm max}$ from Equation (4) and $A$ from Equation (6), we obtain

[TABLE]

Thus we expect the fluence distribution to be insensitive to $\gamma$ (assuming $\gamma<5/2$) and basically given by the value expected from a characteristic burst energy and Euclidean space. This is because, for Euclidean space, the number of bursts above a given fluence $F$ is dominated by those at the distance \sim$$\sqrt{E_{\rm max}/4\pi F} where the brightest bursts are detectable (Macquart & Ekers, 2018).

In Figure 3, we compare our fluence distribution given by Equation (8) with the ASKAP sample. A model with $N_{\rm tot}=20$ and $F_{\rm th}=50\,{\rm Jy\,ms}$ provides a good fit. The sample is incomplete below the sensitivity threshold $F_{\rm th}$ (see §5).

3. Constraints on the Energy Distribution

Combining $N_{\rm tot}$ with ${\rm DM}_{\rm max}$ and $F_{\rm th}$ estimated in the previous discussions allows us to constrain the energy distribution function of the FRB population. To relate luminosity distance and DM we use

[TABLE]

which is scaled to match the DM given by Zhang (2018) based on the latest Planck $\Lambda$CDM cosmology (Planck Collaboration et al., 2016). From the comparisons given in Figures 2 and 3, we can estimate the maximum isotropic FRB energy by using Equation (4),

[TABLE]

This in turn should be an important constraint on any emission model for FRBs (e.g., Lu & Kumar, 2019).

On the other hand, the total exposure of the ASKAP survey was $\Omega T\simeq 5.1\times 10^{5}\rm\,deg^{2}\,hr$ (Shannon et al., 2018), which corresponds to an effective all-sky survey time $T_{\rm eff}\simeq 12.4\,\rm hr$. Since most bursts in our sample have moderate redshifts $z\sim 0.5$, we replace the Euclidean volume $4\pi D_{\rm max}^{3}/3$ in Equation (6) with the comoving volume given by the Planck cosmology. For our fiducial values $N_{\rm tot}=20$, ${\rm DM}_{\rm max}\simeq 750\rm\,pc\,cm^{-3}$, and $\gamma\simeq 1.7$, the volumetric rate normalization is estimated to be

[TABLE]

Taking $E_{\rm max}\simeq 2\times 10^{33}\rm\,erg\,Hz^{-1}$, we can further estimate the volumetric rate of FRBs above energy $E=10^{32}E_{32}\rm\,erg\,Hz^{-1}$,

[TABLE]

This rate density is useful for comparison against potential FRB progenitors, although a critical unknown is how often each source may repeat at a given energy. Non-repeating models based on rare events, such as long gamma-ray bursts or binary neutron star mergers (Totani, 2013; Zhang, 2014; Wang et al., 2016), are inconsistent with the high rate of low-energy FRBs with $E\lesssim 10^{31}\rm\,erg\,Hz^{-1}$.

4. Bayesian Analysis of Full Parameter Space

With a rough analytical understanding of the ASKAP FRB statistics in hand, we next consider a Schechter-like model

[TABLE]

where $\phi_{0}$ (in $\rm Gpc^{-3}\,yr^{-1}$) is the volumetric rate normalization at redshift $z=0$, $(1+z)^{\beta}$ represents possible evolution of FRB rate with redshift, $\gamma$ is the energy power-law index, and $E_{\rm max}$ is the maximum energy beyond which the rate cuts off exponentially. Hereafter, we adopt the latest Planck cosmology and the numerical integral relation $z(\mathrm{DM})$ in Zhang (2018).

Equation (13) predicts the distribution of redshifts (or DMs) and fluences for the observed bursts within the effective all-sky exposure time $T_{\rm eff}$,

[TABLE]

where $dV/dz$ is the differential comoving volume, $x=E/E_{\rm max}$, and $x_{\rm min}=4\pi D_{\rm L}^{2}(z)F(1+z)^{0.5}/E_{\rm max}$ with $D_{\rm L}(z)$ as the luminosity distance. Note that the $(1+z)^{0.5}$ factor in $x_{\rm min}$ is the k-correction assuming an average intrinsic spectrum555For $E_{\nu}=E_{\nu_{0}}(\nu/\nu_{0})^{q}$ with arbitrary $q$ (where $\nu_{0}$ is the observer’s reference frequency), we have $x_{\rm min}=4\pi D_{\rm L}^{2}(z)F_{\nu_{0}}(1+z)^{-q-1}/E_{\nu_{0},\rm max}$. The cosmological rate evolution parameter $\beta$ is degenerate with the spectral index $q$ due to k-correction. Variation in spectral index $\Delta q$ roughly corresponds to a linear shift in $\beta$ by $\Delta\beta\simeq(1-\gamma)\Delta q$. We note that the average FRB spectrum across a wide frequency range has not been well measured so far, so the cosmological evolution of FRB rate is not meaningfully constrained. The other three parameters are not affected by this degeneracy. of $E_{\nu}\propto\nu^{-1.5}$ (Macquart et al., 2019).

We carry out a Bayesian analysis of the 4-dimensional parameter space $\boldsymbol{p}=(\mathrm{log}\phi_{0}$, $\beta$, $\gamma$, $\mathrm{log}E_{\rm max})$ by comparing the model prediction with the observed distribution as seen in Figure 1. The probability density function (PDF) of the parameter vector is given by (besides a normalization factor)

[TABLE]

where $L(\boldsymbol{D}|\boldsymbol{p})$ is the likelihood for the data $\boldsymbol{D}$ to occur under a given parameter vector $\boldsymbol{p}$, and we take a flat prior $f_{0}(\boldsymbol{p})$ for sufficiently wide ranges of parameters: $\mathrm{log}\phi_{0}\,\mathrm{[Gpc^{-3}\,yr^{-1}]}\in(-2,\ 3.5)$, $\beta\in(-5,\ 5)$, $\gamma\in(1.1,\ 2.4)$, $\mathrm{log}E_{\rm max}\,\mathrm{[erg\,Hz^{-1}]}\in(31,\ 36)$. The likelihood function measures the “goodness of fit” of how well the model at $\boldsymbol{p}$ fits the data $\boldsymbol{D}$. Since the data is a 2-dimensional distribution $N_{\rm obs}(\leq z,\geq F)$, information will be lost if we were to perform two marginalized 1-dimensional KS (or $\chi^{2}$) tests (as shown in Figures 2 and 3). In the spirit of the 2-dimensional KS test (Peacock, 1983; Fasano & Franceschini, 1987), we take the following likelihood function

[TABLE]

where $\lambda_{1}^{(i)}=N(>z^{(i)},\geq F^{(i)})$ and $\lambda_{2}^{(i)}=N(\leq z^{(i)},\geq F^{(i)})$ are the expected number of detections in the first and second quadrants, $k_{1}^{(i)}=N_{\rm obs}(>z^{(i)},\geq F^{(i)})$ and $k_{2}^{(i)}=N_{\rm obs}(\leq z^{(i)},\geq F^{(i)})$ are the actual number of detections, $i$ is the index of the observed bursts in our sample, and $P(k,\lambda)=\lambda^{k}\mathrm{e}^{-\lambda}/k!$ is the Poisson PDF. The cumulative numbers of observables ($z$ or $F$) are not independent of each other (Jauncey, 1967) and are hence only used once in the likelihood function in Equation (16).

The number of detections $N(>z^{(i)},\geq F^{(i)})$ is calculated by integrating Equation (14) from the redshift of the $i$-th burst $z^{(i)}$ to a maximum value of $z_{\rm max}=2$. The precise value of $z_{\rm max}$ is not important, as long as it is much larger than the maximum redshift $\max(z^{(i)})\sim 0.7$ in our sample. We do not make use of the cumulative numbers in the third and fourth quadrants $N(\leq z^{(i)},<F^{(i)})$ and $N(>z^{(i)},<F^{(i)})$, because these numbers are subjected to uncertainties of the telescope threshold fluence $F_{\rm th}$ (see §5 for a discussion).

We then generate $10^{6}$ Markov-Chain Monte Carlo (MCMC) samples of the parameter vector $\boldsymbol{p}$ based on the PDF $f(\boldsymbol{p})$ in Equation (15), as shown in Figure 4. The four parameters are constrained at 68% ($1\sigma$) confidence level to be $\mathrm{log}\phi_{0}\,\mathrm{[Gpc^{-3}\,yr^{-1}]}\simeq 1.1^{+0.8}_{-1.2}$, $\beta\simeq 0.8^{+2.6}_{-2.9}$, $\gamma\simeq 1.6\pm 0.3$, $\mathrm{log}E_{\rm max}\,\mathrm{[erg\,Hz^{-1}]}\simeq 34.1^{+1.1}_{-0.7}$. The cumulative rate above certain energy $E\ll E_{\rm max}$ at $z=0$, given by $N(>E)=\phi_{0}(E/E_{\rm max})^{1-\gamma}/(\gamma-1)$, is constrained to be $\mathrm{log}N(>10^{32}\rm\,erg\,Hz^{-1})\,\mathrm{[Gpc^{-3}\,yr^{-1}]}\simeq 2.6\pm 0.4$ (68% confidence). These constraints, although weak due to small-number statistics, broadly agree with our simple analytical results in §2 and §3. Our model will provide tighter constraints on the FRB rate function when applied to the much larger samples expected in the near future.

Finally, we compare the energy distribution function given by the standard “$1/V_{\rm max}$” estimator (Schmidt, 1968) with our Bayesian results in Figure 5. We calculate the weighted sum $\sum_{i}(V_{\mathrm{max},i}T_{\rm eff})^{-1}$ for each energy bin666For the small ASKAP sample, the result of the Schmidt estimator is sensitive to binning because of Poisson error. As we group the bursts into logarithmic energy bins, a smooth distribution is obtained when the number of bins is $\leq 5$. ($E$, $E+\Delta E$), where $V_{\mathrm{max},i}$ is the maximum observable comoving volume for event $i$ in this bin and $T_{\rm eff}$ is the effective all-sky survey time. For the Schmidt estimator, we assume the FRB rate to be independent of redshift, since the redshift evolution is not strongly constrained by the small ASKAP sample. We also show the energy distribution functions (Equation 13) evaluated at redshift $z=0.5$ (where most bursts are located) for $10^{3}$ randomly selected MCMC samples of the parameter vector $\boldsymbol{p}$ from our Bayesian analysis. The reasonable agreement between these two independent methods means that our understanding of the ASKAP FRB statistics is physical.

5. Sample Completeness

It is possible that the ASKAP survey has missed some bursts with fluence above our selection threshold $F_{\rm th}=50\rm\,Jy\,ms$, which would cause the sample to be biased. The missing bursts may affect the overall DM and fluence distribution and change the conclusions in the previous sections (which are drawn from a potentially incomplete sample). In this section, we show that the number of missing bursts above $50\rm\,Jy\,ms$ is small ($\lesssim$3) and that the differential incompleteness is a weak function of DM (i.e., the missing bursts do not preferentially have large or small DMs).

The signal-to-noise ratio (S/N) for a burst with fluence $F$ and duration $w_{\rm obs}$ can be estimated by (Shannon et al., 2018)

[TABLE]

where $B=336\,$MHz is the observing bandwidth, $S_{\rm sys}\simeq 2000\,$Jy is the system equivalent flux density of a single beam, and $w_{\rm obs}=(t_{\rm DM}^{2}+t_{\rm samp}^{2}+t_{\rm arr}^{2})^{1/2}$ is the total duration from a convolution of the DM smearing $t_{\rm DM}$, the sampling time $t_{\rm samp}=1.26\,$ms, and the burst width at arrival $t_{\rm arr}$ (which is given by the emission width, redshift, and scattering broadening). We have ignored the residual time dispersion $t_{\rm\delta DM}$ due to de-dispersion with a slightly incorrect DM (Cordes & McLaughlin, 2003). The DM smearing is given by

[TABLE]

where $\rm DM_{\rm tot}$ is the total DM of a burst, $\Delta\nu=1\,$MHz is the the spectral resolution, and $\nu=1.32\,$GHz is the central frequency of the ASKAP survey. From Equations (17) and (18), we see that finite temporal and spectral resolution always decreases the S/N by broadening the observed duration $w_{\rm obs}$.

We caution that Equation (17) only provides an approximate estimate of the S/N, because $S_{\rm sys}$ varies with telescope operational status and the burst’s localization (see James et al., 2019b, for detailed discussions), and some narrow-spectrum bursts only partially fill the bandwidth $B$. Additionally, the burst spectrum may not peak at $\nu=1.32\,$GHz, so the DM smearing may deviate from Equation (18) by up to \sim$$30\%. Nevertheless, the discussion in this section stays qualitatively true even under these variations.

In Figure 6, we show the distribution of the observed widths and fluences, following Keane & Petroff (2015). The observed widths cluster around $\sim$3 ms, because of the sampling time $t_{\rm samp}\sim 1\,$ms and the DM smearing $t_{\rm DM}\sim 2\,$ms (such that many shorter bursts with $t_{\rm arr}\lesssim 1\,$ms have been broadened). The longest duration burst has $w_{\rm obs,max}=7.3\,$ms (not strongly affected by DM smearing) and fluence $F=95\rm\,Jy\,ms$. We assume that the total duration $w_{\rm obs}$ is not correlated with the fluence (no such correlation has been reported), so the number of bursts in any $(w_{\rm obs},w_{\rm obs}+\Delta w_{\rm obs})$ bin goes as $N(>F)\propto F^{-3/2}$ as shown in Figure 3. Thus, a conservative estimate of the number of missing bursts in the blue-shaded region is given by $2\times(50^{-1.5}-95^{-1.5})/95^{-1.5}\simeq 3.2$, where we have assumed that all bursts of width near $\sim$7.3 ms and fluence below 95 Jy ms have been missed. Therefore, only a small fraction (about $10\%$) of the bursts above our fluence treshold $F_{\rm th}=50\rm\,Jy\,ms$ may have been missed in the ASKAP sample. This may be counter-intuitive given the relatively coarse spectral resolution of $\Delta\nu=1\,$MHz. It is the high fluence threshold (only selecting those extremely bright bursts) that pushes the incompleteness region to very long durations where FRBs are sufficiently rare.

Figure 7 shows the total DM and S/N of the main beam of detection for each burst. We also show a few curves of the S/N expected for bursts of a given fluence $F$ (in $\rm\,Jy\,ms$) and arrival width $t_{\rm arr}$ (in $\rm ms$). The variation $\rm{S/N}(DM)$ in each curve is only due to DM smearing $t_{\rm DM}$. We find that (i) long-duration faint bursts, represented by the $(F,\ t_{\rm arr})=(50,5)$ curve, are likely missed in our sample; and (ii) for long-duration bursts, represented by the two curves with $t_{\rm arr}=5$, the S/N depends weakly on $\rm DM_{\rm tot}$, because the total duration $w_{\rm obs}$ is dominated by the arrival width $t_{\rm arr}$ ($\gg t_{\rm DM}$). Thus, even if some long-duration faint bursts have been missed in our sample, the differential incompleteness777In the future, when a better understanding of the arrival durations ($t_{\rm arr}$, including scattering broadening) is available, it is straightforward to include a bias correction factor $f_{\rm bias}(\rm DM)$ for bursts of different DMs in our Bayesian analysis. is a weak function of DM (i.e., the missing bursts do not preferentially have large or small DMs). The abrupt cutoff in the number of bursts with $\rm DM_{\rm tot}\gtrsim 750\rm\,pc\,cm^{-3}$ at $\rm S/N>9.5$ and $F>50\rm\,Jy\,ms$ is not caused by loss of S/N due to excessive DM smearing. Instead, it can only be due to the lack of these events arriving at the telescope, which is explained by a cutoff in the volumetric rate above a maximum burst energy $E_{\rm max}$ in our model.

Recently, Connor (2019) studied how the finite time and frequency resolution affects the distributions of the observed burst duration, DM, and flux. The author argued that a large fraction of FRBs with short durations $t_{\rm arr}\ll 1\,$ms and high DMs may have been missed in current surveys, including ASKAP. Our approach is different in that we aim to have a fluence-complete sample above a threshold S/N whereas Connor (2019) focused on the observed bursts above a threshold detection flux (see their Equation 16). We know that $\mathrm{S/N}\propto Fw_{\rm obs}^{-1/2}\propto Sw_{\rm obs}^{1/2}$ (Equation 17), where $S=F/w_{\rm obs}=S_{0}t_{\rm arr}/w_{\rm obs}$ is the detection flux and $S_{0}$ is the “intrinsic” flux arriving at the telescope. Thus, fluence completeness depends crucially on the number of longest-duration bursts whereas flux completeness relies on shortest-duration ones. If there exists a large number of short bursts with $t_{\rm arr}\ll 1\,$ms (possibly down to nano-second timescale), as suggested by the Parkes sample (Ravi, 2019) and giant pulses of Galactic pulsars, it is extremely difficult to achieve a flux($S_{0}$)-complete sample. This is the motivation for why we choose to study the energy distribution function of FRBs but not the luminosity function.

6. Discussion

We note a few caveats and highlight some implications of our study.

(1) Instead of taking $F_{\rm th}=50\rm\,Jy\,ms$, we have also tested a more conservative choice of fluence threshold $F_{\rm th}=60\rm\,Jy\,ms$ (which includes 15 out of the 23 bursts). According to §5, this smaller sample is better fluence-complete with less bias against long-duration faint bursts. We find that our model still provides a good fit to the DM and fluence distributions, although the statistical constraints on the model parameters are slightly worse.

(2) We have assumed that the DM excess (beyond the Milky Way’s contribution) is largely due to the intergalactic medium and hence can be used as a distance indicator. However, as we know from FRB 121102 (Tendulkar et al., 2017), the host galaxy and the circum-burst medium may contribute a significant fraction of DM. Unfortunately, the properties of FRB host galaxies are still highly uncertain. We also tried subtracting a constant host-galaxy contribution of $30\rm\,pc\,cm^{-3}$. The resulting distribution $N(<{\rm DM})$ becomes slightly shallower on the low-DM end, and the Bayesian analysis gives similar parameter constraints. In the future, with a larger FRB sample, it is possible to constrain the (averaged) host-galaxy and circum-burst DM (e.g., from a supernova remnant, Connor et al., 2016; Piro, 2016; Piro & Burke-Spolaor, 2017; Piro & Gaensler, 2018) by studying the deviation of the DM distribution from a power-law on the low-DM end. On the other hand, FRBs at $z\gtrsim 0.8$ may randomly have their sight lines intersecting with a few galactic haloes (McQuinn, 2014; Prochaska & Zheng, 2019). This effect causes stochastic deviations of the DM distribution from our model.

(3) The scaling in Equation (8) breaks down at sufficiently small fluences. This is because (i) the spacetime is no-longer Euclidean at high redshift (the co-moving volume increases much slower than $z^{3}$), and (ii) $N(>F)$ is no longer dominated by bursts with $E\sim E_{\rm max}$ when the majority of them in the Universe have already been counted at higher fluences. For instance, at $F\lesssim E_{\rm max}/4\pi D^{2}_{\rm L}(z=5)\simeq 1\rm\,Jy\,ms$, the cumulative number $N(>F)$ should be dominated by bursts with $E\ll E_{\rm max}$, and hence the fluence distribution becomes $N(>F)\propto F^{1-\gamma}$. This means that a sensitive telescope with fluence threshold $\ll 1\rm\,Jy\,ms$ can directly measure the energy distribution slope $\gamma$ by source counting. This may explain the difference in logN-logF slopes given by the Parkes and ASKAP samples (James et al., 2019a), since the former has a much lower fluence threshold. In the future, with a much larger (well-controlled) sample down to a lower fluence threshold, it is possible to constrain the redshift evolution of FRB rate (Macquart & Ekers, 2018), which will in turn constrain the progenitor models for FRBs.

(4) The first repeater has been localized to a metal poor star-forming dwarf galaxy, suggesting that there may be some relation to a long gamma-ray bursts and superluminous supernovae (Metzger et al., 2017). As more FRBs are localized, an important constraint to the progenitor model will be how much variety is seen in the hosts.

If, for example, some FRBs are found to have distinct hosts from others, one might conclude that they have different progenitors. If the distributions of the various FRB populations continue to obey similar energy distributions and maximum energies as we describe here, it argues that the underlying source is the same even if the situation that generated the source is different. For instance, in the magnetar picture, one could imagine magnetars generated both from young massive stars and old stellar environments by neutron star or white dwarf mergers. FRBs from these two populations could vary in a number of ways, from the host galaxies, to the local DM, to the rotation measure and presence of a persistent radio source. But by looking at the burst statistics as we describe here, we may understand if the underlying source is the same.

7. Summary

This work uses the ASKAP sample, which we show to be reasonably complete above $50\rm\,Jy\,ms$, to study the energetics and cosmological rate of the whole FRB population.

We find that the ASKAP sample is well described with a power-law energy distribution with slope $\gamma=-d{\rm log}N/d{\rm log}E\sim 1.7$, independent of whether they individually repeat themselves. This is because the observed DM distribution scales as $N(>{\rm DM})\propto\mathrm{DM}^{5-2\gamma}$ (for $\gamma>1$). The abrupt cutoff in the number of bursts with $\rm DM\gtrsim 750\rm\,pc\,cm^{-3}$ suggests that the FRB population has a maximum specific energy of $E_{\rm max}\sim\mathrm{a\,few}\times 10^{33}\rm\,erg\,Hz^{-1}$ above which the volumetric rate rapidly drops. For a spectral width of \sim$$1\,GHz, this implies a maximum isotropic energy of order $10^{42}\rm\,erg$ or isotropic luminosity of order $10^{45}\rm\,erg\,s^{-1}$ for millisecond duration, which is a factor of a few to ten higher than that found by Luo et al. (2018) from the Parkes sample.

The existence of a maximum energy $E_{\rm max}$ causes the observed fluence distribution to be $N(>F)\propto F^{-3/2}$ at high fluences, because the number of bursts is dominated by those at distances $D\sim\sqrt{E_{\rm max}/4\pi F}$. At sufficiently low fluences, the distribution approaches $N(>F)\propto F^{1-\gamma}$, reflecting the intrinsic energy distribution. The transition between these two regimes is sensitive to the evolution of FRB rate at high redshift.

Besides simple analytical insights, we also present a numerical Bayesian analysis of the ASKAP sample by comparing the observed distributions with the predictions from a Schechter-like model for the FRB rate (Equation 13). The results are shown in Figure 4. The energy distribution power-law index and the maximum energy are constrained to be $\gamma\simeq 1.6\pm 0.3$ and $\mathrm{log}E_{\rm max}\,\mathrm{[erg\,Hz^{-1}]}\simeq 34.1^{+1.1}_{-0.7}$ (68% confidence), respectively. From the survey exposure time, we further infer a cumulative local volumetric rate of $\mathrm{log}N(E>10^{32}\rm\,erg\,Hz^{-1})\,\mathrm{[Gpc^{-3}\,yr^{-1}]}\simeq 2.6\pm 0.4$ (68% confidence). Finally, we compare the energy distribution function given by the standard “$1/V_{\rm max}$” estimator with that from our Bayesian approach. These two independent methods give consistent results as shown in Figure 5, which means that our understanding of the ASKAP FRB statistics is physical. Our model will give tighter constraints on the statistical properties of FRB rate when applied to the much larger samples to be collected in the near future.

We thank the two anonymous referees for carefully reading the manuscript and providing valuable criticisms and suggestions. We thank Jonathan Katz, Liam Connor, and Ue-Li Pen for useful comments. WL is supported by the David and Ellen Lee Fellowship at Caltech.