On the non-Poissonian repetition pattern of FRB121102

TL;DR

This paper analyzes the irregular repetition pattern of FRB121102, showing that a Weibull distribution better models burst intervals than Poissonian statistics, with implications for detection strategies.

Contribution

It introduces a Weibull distribution model for FRB121102 burst intervals, capturing burst clustering and improving understanding of its repetition pattern.

Findings

Burst intervals are better described by Weibull distribution than Poisson.

Estimated mean repetition rate is approximately 5.7 per day.

Clustered burst behavior affects detection probabilities and observation planning.

Abstract

Rhe Fast Radio Burst FRB121102 has been observed to repeat in an irregular fashion. Using published timing data of the observed bursts, we show that Poissonian statistics are not a good description of this random process. As an alternative we suggest to describe the intervals between bursts with a Weibull distribution with a shape parameter smaller than one, which allows for the clustered nature of the bursts. We quantify the amount of clustering using the parameters of the Weibull distribution and discuss the consequences that it has for the detection probabilities of future observations and for the optimization of observing strategies. \new{Allowing for this generalization, we find a mean repetition rate of r=5.7^{+3.0}_{-2.0} per day and index k=0.34^{+0.06}_{-0.05} for a correlation function \xi(t)=(t/t_0)^{k-1}.

Figures

Tables

Table 1. Table 1: Fraction of the simulations with independent intervals ( second column ) and those with finite observational separations ( third column ) for which the parameter values assumed in the simulations lie within the 68 %-, 95 %-, and 99 %-confidence regions derived from the simulated data.

confidence level	independent	finite
68 %	59 %	63 %
95 %	93 %	93 %
99 %	98 %	98 %

Equations

$$\mathcal{P}(\delta |r)=r{\mathrm{e}}^{-\delta r},$$

$$\mathcal{W}(\delta |k,r)=k{\delta }^{-1}{\left[\delta r\Gamma \left(1+1/k\right)\right]}^k{\mathrm{e}}^{-{\left[\delta r\Gamma \left(1+1/k\right)\right]}^k},$$

$$\Gamma (x)={\int }_0^{\infty }\mathrm{d}t{t}^{x-1}{\mathrm{e}}^{-t}$$

$$r^{-1}={⟨\delta ⟩}_{(\delta |k,r)}={\int }_0^{\infty }\mathrm{d}\delta \delta \mathcal{W}(\delta |k,r).$$

$$\mathcal{P}(N,t_1,\dots ,t_N|k,r)$$

$$={\int }_{t_N}^{\infty }\mathrm{d}{t}_{N+1}P(N|t_1,\dots ,t_N,t_{N+1})\mathcal{P}(t_1,\dots ,t_N,t_{N+1}|k,r).$$

$$P(N|t_1,\dots ,t_N,t_{N+1})$$

$$\mathcal{P}(t_1,\dots ,t_N,t_{N+1}|k,r)$$

$$=\mathcal{P}(t_1|k,r)\prod _{i=1}^{N-1}\mathcal{W}(t_{i+1}-t_i|,k,r).$$

$$P(\text{obs. start during }{\delta }_{01}|{\delta }_{01})$$

$$\mathcal{P}({\delta }_{01}|k,r,\text{obs. start during }{\delta }_{01})$$

$$\propto P(\text{obs. start during }{\delta }_{01}|{\delta }_{01})\mathcal{P}({\delta }_{01}|k,r)$$

$$\propto {\delta }_{01}\mathcal{W}({\delta }_{01}|k,r).$$

$$\mathcal{P}({\delta }_{01}|k,r,\text{obs. start during }{\delta }_{01})=\frac{{\delta }_{01}}{\bar{\delta }}\mathcal{W}({\delta }_{01}|k,r),$$

$$\bar{\delta }={\int }_0^{\infty }\mathrm{d}\delta \delta \mathcal{W}(\delta |k,r)=\frac{1}{r}.$$

$$\mathcal{P}(t_1|{\delta }_{01})=\frac{1}{{\delta }_{01}}\theta ({\delta }_{01}-t_1).$$

$$\mathcal{P}(t_1|k,r)$$

$$=\frac{1}{\bar{\delta }}{\int }_{t_1}^{\infty }\mathrm{d}{\delta }_{01}\mathcal{W}({\delta }_{01}|k,r)$$

$$=\frac{1}{\bar{\delta }}\mathrm{CCDF}(t_1|k,r),$$

$$\mathrm{CCDF}(\delta |k,r)$$

$$\mathcal{P}(N,t_1,\dots ,t_N|k,r)=$$

$$\prod _{i=1}^{N-1}\mathcal{W}(t_{i+1}-t_i|k,r).$$

$$\mathcal{P}(N=1,t_1|k,r)=$$

$$P(N=0|k,r)$$

$$=\frac{1}{\bar{\delta }}{\int }_{\Delta }^{\infty }\mathrm{d}{t}_1\mathrm{CCDF}(t_1|k,r).$$

$$P(N=0|k,r)$$

$${\Gamma }_{\mathrm{i}}(x,z)={\int }_z^{\infty }\mathrm{d}t{t}^{x-1}{\mathrm{e}}^{-t}$$

$$\mathcal{P}(k,r)\propto k^{-1}{r}^{-1}$$

$$\mathcal{P}(\log k,\log r)=\mathrm{const}.$$

$$\mathcal{P}(k,r|N,t_1,\dots ,t_N)$$

$$r_{\alpha }^{(\mathrm{Poiss})}=-\frac{\ln (\alpha )}{\Delta },$$

$$\bar{N}=r\mathrm{d}t.$$

$$P(\nu )={⟨{\left|n(\nu )\right|}^2⟩}_{(N|k,r)}.$$

$$\xi(t)=\left\{\begin{array}[]{ll}k\,\left[\,\Gamma\left(1+1/k\right)\right]^{k}(rt)^{k-1}&t\lesssim 10/r\\ \frac{\Gamma(1+2/k)/\Gamma(1+1/k)^{2}-2}{rt}&t\gg 10/r\end{array}\right..$$

Full text

††pubyear: 2016††pagerange: On the non-Poissonian repetition pattern of FRB121102–References

1 Introduction

Fast Radio Bursts (FRBs) are bright, short-duration (ms), pulses detected in the radio waveband. They are characterized by high dispersion measures ($10^{2}$–$10^{3}\,\mathrm{pc}\,\mathrm{cm}^{-3}$) and their origin is not yet understood. Several bursts have been observed to arrive from a single location and with the same dispersion measure, collectively referred to in the following as FRB121102. None of the other FRBs have been observed to emit repeated pulses. However, it is possible that these have merely escaped detection.

The bursts of FRB121102 seem to happen in an irregular fashion, but appear to be clustered to some degree. This clustering may contain hints to the emission mechanism. Additionally, it has important consequences for the strategy that will have to be followed in order to detect the maximum number of bursts in the future and for the maximum repeat rate that one may infer from non-detections of repetitions for other FRBs, as discussed by Connor et al. (2016).

In this paper, we suggest to model the distribution of intervals between successive bursts of FRB121102. This allows us to consider simple deviations from Poissonian statistics and to quantify the clustering. Any clustering parameter estimation explicitily or implicitly requires a Bayesian prior to quantify error bars on the inferred parameters. To this end, we discuss the Weibull distribution in the next section. We introduce the data we use in Sect. 3, derive the statistical formalism in Sect. 4, and discuss our results and their consequences in Sect. 5. A brief summary of our findings and their consequences is given in Sect. 6.

2 Beyond Poissonian statistics

For a Poissonian point process, a sufficient statistic is the total number of events. Here we consider deviations from purely Poissonian statistics and make use of the additional information given by the distribution of events.

If the bursts of FRB121102 were distributed in time according to a Poisson point process with a constant expectation value, then the distribution of intervals $\delta$ between subsequent bursts would be exponential, i.e.,

[TABLE]

where $r$ is the constant rate of bursts. One possible generalization of this distribution of intervals is the Weibull distribution,

[TABLE]

which has an additional shape parameter $k$. Here,

[TABLE]

is the gamma function. We choose a parameterization in which the parameter $r$ can again be interpreted as a rate, since

[TABLE]

For $k=1$ the Weibull distribution reduces to the Poissonian case. For values of $k$ different from 1, the Weibull distribution describes a certain degree of clustering. Specifically, for $k<1$, small intervals are favored compared to the Poissonian case and thus the presence of a burst makes an additional burst in the near future more likely.

With this parameterization, we can ask whether the existing observations of bursts and observations with an absence of bursts are consistent with purely Poissonian statistics and, if not, we can quantify how strong the clustering of bursts is.

3 Data

For this analysis, we use the bursts listed in Table 1 of Spitler et al. (2016) and Table 3 of Scholz et al. (2016), as well as the starting and stopping times of the observations listed in Table 2 of Scholz et al. (2016). To get all times onto equal footing, we convert the starting and stopping times to barycentric modified Julian dates and apply a frequency-dependent dispersion correction. These corrections have already been applied to the published burst arrival times and we apply them here to the observation intervals to obtain the correct time differences. The lengths of the resulting observational intervals and the location of each burst within these intervals are shown in Fig. 1.

Further follow-up observations of FRB121102 have been conducted and more bursts have been detected (e.g., Chatterjee et al., 2017). However, the exact duration of each observation is not published and thus we do not include these in our analysis. While more data will constrain the parameters more tightly, the qualitative conclusions of this study will not change once the necessary data from newer observations are included.

4 Formalism

Consider first a single observation of a finite duration $\Delta$. The information from this observation consists of the number of observed bursts, $N$, and the times at which these bursts happened, $t_{1},t_{2},\dots,t_{N}$. For simplicity, we measure times with respect to the start time of the observation. Figure 2 illustrates the variables we use.

We calculate the likelihood for $N$ and $t_{1},\dots,t_{N}$ by marginalizing over the time of the next burst, $t_{N+1}$,111We use $\mathcal{P}$ to denote probability densities and $P$ to denote probabilities.

[TABLE]

The first ingredient is the probability that $N$ bursts are observed, given the times of the $N+1$ first bursts. This is independent of the parameters of the Weibull distribution and trivially given by

[TABLE]

where $\theta(\cdot)$ is the Heaviside step function, which is 1 for a positive argument and 0 otherwise.

The second ingredient, the probability for the arrival times of the $N+1$ first bursts, given the parameters of the Weibull distribution, can be split into a product of probabilities for the inter-burst intervals,

[TABLE]

The first factor in the last equation, the probability density for the time of the first observed burst, is slightly more complicated than the others. Since there is no previously observed burst, this is not a pure Weibull distribution. However, we will relate the two distributions in the following way: Consider the interval $\delta_{01}=t_{1}-t_{0}$ between the last unobserved burst and the first observed burst. The probability distribution for the length of this interval is characterized by the Weibull distribution, which describes the distribution for interval lengths in general, and the fact that the observation started within this interval. The probability for the latter fact is proportional to the length of the interval so that we can write

[TABLE]

and thus

[TABLE]

To ensure the correct normalization, the full probability distribution has to be

[TABLE]

where

[TABLE]

Further, for symmetry reasons, we have

[TABLE]

Thus, we can write the first factor in Eq. (4) as

[TABLE]

where we have defined the complementary cumulative distribution function

[TABLE]

in the last line.

Plugging everything into Eq. (4), we are left with a product of $N-1$ Weibull distributions, a cumulative Weibull distribution due to the time difference $\Delta-t_{N}$ at the end of the observation, and another cumulative Weibull distribution due to the time elapsed at the beginning of the observation before the first burst. As argued above, this last factor has an additional factor $1/\bar{\delta}$. In total, we obtain

[TABLE]

This expression is valid for any distribution of intervals. To obtain the specific expression for the case of the Weibull distribution, we plug in Eqs. (2), (11), and (14).

Strictly speaking, Eq. (4) is only valid for $N>1$. For $N=1$ it simplifies trivially to

[TABLE]

For the case $N=0$, i.e., an observation without detected bursts, we follow a similar argument as above and marginalize over the time of the first burst after the start of the observations,

[TABLE]

Plugging in Eq. (14) yields the result for the Weibull distribution,

[TABLE]

where

[TABLE]

is the incomplete gamma function.

This completes the derivation of the likelihood for a single finite-duration observation. In the limit that successive observations are widely spaced when compared to the spacing of the bursts, the individual observations become independent and the likelihood turns into a simple product of the likelihoods of the individual observations. In the following, we follow this approximation. In Appendix A we show that the error introduced by this simplification is small.

To infer the parameters $k$ and $r$ we need to multiply the likelihood with a prior and calculate the posterior distribution for the two parameters. Since both parameters are strictly positive we work with their logarithms. Specifically, we choose independent Jeffreys priors for $k$ and $r$, i.e.,

[TABLE]

so that

[TABLE]

We then calculate the posterior as

[TABLE]

5 Results and discussion

We make use of the data shown in Fig. 1 and insert them into Eqs. (4), (16), and (18). Figure 3 shows the resulting two-dimensional posterior distribution for $k$ and $r$, as well as the one-dimensional marginal posteriors for the individual parameters.

5.1 Non-Poissonian nature of FRB121102

From Fig. 3 we immediately see that the Poissonian case of $k=1$ is strongly disfavored by the data. The posterior mean value for the shape parameter is $\left<k\right>_{(k|N,t)}=0.34^{+0.06}_{-0.05}$ and the posterior mean for the rate is $\left<r\right>_{(r|N,t)}=5.7^{+3.0}_{-2.0}/\mathrm{day}$ The uncertainties given here are the 68 % confidence intervals. The estimated rate is similar to the Poissonian rate, but the uncertainty interval is about twice as wide. Since the value of $k$ is less than 1, there is indeed an excess of short inter-burst intervals. A consequence of this is that, if a burst is observed, the probability that another burst happens shortly after is increased.

Although we have not exhausted all possible distributions to model the data, it can be quantitatively shown that the Weibull distribution describes the data much better than Poissonian case. We examine this by Kolmogorov-Smirnov (KS) test on the distribution of burst intervals. For Poissonian case, the distribution of intervals disagrees with the exponential best fit at level $\alpha=1.03\times 10^{-6}$. If we assume a Weibull distribution, we have to consider 3 kinds of intervals ($N=2$, $N=1$ and $N=0$ as discussed in §4) separately. For each of these, the distribution depends on the Weibull parameters $k,r$, and additionally the distribution of observation durations (see Appendix B for more details). These three distributions disagree with the Weibull model at levels $\alpha=1.46,1.18\times 10^{-4},0.37$ respectively. It shows that Weibull model performs much better than Poissonian. Note that, although much better than Poissonian in $\alpha$ the Weibull distribution still does not describe the data well, especially in the distribution of “$N=1$”-kind intervals. The reason is obvious – from Fig. 1 we see many clustered bursts to accumulate in short observations, and even very likely to be in the beginning or ending of one observation.

To illustrate the effect this has on the number of bursts during an observation, we consider the case of a Weibull distribution with parameters fixed to the posterior-mean values. We numerically calculate the probability to see $N$ bursts during an observation that lasts ten times the mean burst separation, i.e., $\Delta=10/r$. For the posterior-mean repetition rate of FRB121102, this is a 42-hour observation. The result is shown in Fig. 4. For comparison, we also plot the Poissonian probability. While the mean number of bursts is the same for both cases, the Weibull case has a much higher probability of yielding no bursts at all. The excess of short intervals makes it likely that either no burst is seen (if the observation happens to fall within a long inter-burst interval) or several bursts are seen, since the presence of one burst makes the presence of other bursts in the temporal vicinity more likely.

An important consequence is that the inference that can be done on the mean repetition rate from observing an FRB, or a random location in the sky, for some time $\Delta$ and not seeing any bursts is much weaker than in the Poissonian case. A possible line of reasoning is to equate the probability of seeing no bursts in such an observation with a probability threshold $\alpha$, to solve for the threshold rate $r_{\alpha}$, and to conclude that rates above this threshold value are disfavored by the observational evidence. For the Poissonian case, this leads to a threshold rate of

[TABLE]

whereas the Weibull result can be obtained by numerically solving Eq. (18) for $r$. Figure 5 shows the ratio of the threshold rate in the Weibull case, $r^{\mathrm{(Weib)}}$, and the threshold rate in the Poissonian case, as a function of the shape parameter $k$. For the posterior-mean value, we find a ratio of 13 for $\alpha=5\,\%$, meaning that rates that are consistent with the data of such an observation can be 13 times as high when allowing for $k\neq 1$ than when assuming Poissonian statistics. Choosing a lower probability threshold $\alpha$ makes this ratio even larger.

Thus, it is possible within the Weibull model that an FRB that has been observed not to repeat for some time does indeed repeat with a fairly high average rate, but with a value of $k$ that is significantly less than 1. This opens the door to the possibility that all observed FRBs are in fact the same class of objects, even though only one object has been observed to burst repeatedly.

5.2 Clustering

To further quantify the degree of clustering, we consider the power spectrum of the burst density for the Weibull distribution with the parameters fixed to the posterior-mean values. If we discretize the time axis into bins of width $\mathrm{d}t$, the mean number of bursts during a bin is

[TABLE]

We can calculate a density of bursts by dividing by the bin width, $n(t)=N(t)/\mathrm{d}t$, and use its Fourier-space version $n(\nu)$ to define the power spectrum,

[TABLE]

For a Poissonian point process, the power spectrum is flat, i.e., independent of $\nu$. We choose the normalization of the Fourier transform in such a way that the value of the power spectrum, except for the zero frequency, is the mean number of bursts per day. In Fig. 6 we show the numerical result for the power spectrum of Weibull-distributed events after subtracting the Poissonian power, i.e., we show the excess of power due to the clustering. Two extremes are of interest: For frequencies much larger than the mean rate $r$, the clustering becomes unimportant and the Poissonian power spectrum is recovered. The clustering power declines as a power law in frequency. For frequencies much smaller than the mean rate the power spectrum becomes flat as well, but with an increase in power by a factor $(\Gamma(1+2/k)/\Gamma(1+1/k)^{2}-1)$. We derive this result in Appendix C. At intermediate frequencies, the power spectrum smoothly transitions between these two regimes and it is this range that describes the non-trivial clustering.

A common parametrization of temporal clustering is through the Fourier Transform of the power spectrum, called the two-point correlation function, $\xi(t)\equiv(2\pi)^{-1}\int_{0}^{\infty}P(\nu)\,\exp(it\nu)\,\mathrm{d}\nu$, which describes the fractional excess over Poissonian at time lag $t$. The asymptotic forms of the correlation function at short and long lags are

[TABLE]

5.3 Observational strategies

The fact that the bursts appear clustered should inform future observational campaigns that aim to detect more bursts. In the Poissonian case, the probability of detecting a burst depends only on the total duration of the observation. In the presence of clustering, however, an observational strategy that spreads the same observing time over a larger period by introducing gaps has a higher chance of detecting a burst.

To illustrate this effect, we consider again the example of an observation whose total duration is such that the expected number of bursts is 10. However, we now split this time into ten equally long observations with equally long gaps between them. Figure 7 shows the numerically calculated probability of detecting at least one burst in such an observation, as a function of the length of each observing gap. Evidently, spreading the observations out increases the odds of a successful detection. Here we have again fixed the parameters to the posterior-mean values and the Poissonian probability of detecting at least one burst is 99.995 %.

Calculations such as the one shown here can be done for arbitrary observing strategies and may therefore be useful for their optimization.

6 Summary

We have used the Weibull distribution as a generalization from Poisson statistics. This is not the only possible extension. Examples for other extensions could be a gamma distribution of intervals or a Poissonian model with time-varying rate, which could for example be modeled as a log-normal random field. We use the Weibull distribution as a simple, mathematically tractable, model that allows us to demonstrate that the available information on burst arrival times and the observation timing strongly disfavors Poissonian statistics and favors instead a clustered distribution. The most important consequences are:

• •

For a continuous observation, a non-detection of bursts is much more likely than in the Poissonian case.

• •

Consequently, upper limits on repetition rates derived from non-detections become much looser and rate estimates in general less certain.

• •

A spread-out observation interspersed with gaps is better suited for detecting further bursts than continuous observations with the same total observation time.

Appendix A Effect of the finite time between observations

In our analysis we treat each of the observational intervals shown in Fig. 1 as independent. In reality, however, they are separated by finite gaps and the correlations in the burst pattern extend, at least in principle, across these gaps and thus correlate the observations across observational intervals. Here we investigate the severity of this effect.

We run two series of simulations of Weibull-distributed events, with parameters fixed to the posterior-mean values. In the first series, we simulate events in the intervals shown in Fig. 1 assuming that they are completely independent, i.e., infinitely separated. In the second series, we simulate the intervals at their actual separation. We then apply the analysis described in Sect. 4 to each simulated data set. If the finite separation of the observations had a great effect, we would expect to see posterior distributions that are different for the two simulation series and are shifted away from the simulated parameters in the case of the realistic simulation. However, neither appears to be the case for the two times 240 simulations that we have run.

We check the accuracy of our model by computing the fraction of simulations for which the true parameter values lie within each contour. The result of this test is shown in Table 1 and confirms that the size of our uncertainty region is accurate.

As a null test, we have also run a similar set of simulations where we set $k=1$, i.e., the Poissonian case. In this case, we do find a posterior that is centered on $k=1$, thus confirming that the non-Poissonian structure we infer is not an artifact of our analysis.

Appendix B Probability density functions of intervals

We consider a number of events (e.g. FRB bursts here) detected by an independent series of observations with finite times. The duration of each observation $T_{\rm obs}$ obeys $f_{\rm obs}$. In the case of $T_{\rm obs}\rightarrow\infty$, we can observe the intrinsic statistics of intervals between two events – we denote it $f^{\infty}_{2}(\Delta t)$. Realistically, when $T_{\rm obs}$ is finite, we observe 3 kinds of intervals – combinations between observation boundaries and events. We define $t_{0}$ as time of interval between two boundaries, $t_{1}$ as time of interval between one event and one boundary, and $t_{2}$ as time of interval between two events. We study the distribution of $\{t_{0},t_{1},t_{2}\}$ given $f^{\infty}_{2}$ and $f_{\rm obs}$.

One can easily proof that the PDF of $t_{0}$ can be written as $f_{0}\propto P_{0}(t)f_{\rm obs}(t)$, where $P_{0}$ is the probability of no detections within $T_{i}$. For discrete cases, $f_{\rm obs}=N^{-1}\sum_{i=1}^{N}\delta(t-T_{i})$, where $\delta(t)$ is the Dirac delta function. It is straightforward to integrate and get the CCDF of $t_{0}$,

[TABLE]

where $m$ is the smallest index s.t. $T_{m}>t$.

To produce $t_{1}=t$, we need condition(1): Given an event, the probability that the next/previous event happens at least $t$ after/before it222condition(1)$\,\sim\,$condition(1*′*): Given a boundary, the next/previous event happens $t$ away from it. $P_{1}(t)=1-F^{\infty}_{2}$ 333In Eq.(27), $P_{0}(t)=1-F^{\infty}_{1}(t)$., and condition(2): the probability that $T_{i}>t$, which is $1-F_{\rm obs}(t)$. Thus,

[TABLE]

To produce $t_{2}=t$, we need the probability of $\Delta t=t$ selected from $f^{\infty}_{2}(\Delta t)$, the probability of $T_{i}>t$, and the probability that $t$ can be put in $T_{i}$, which is $T_{i}-t$. So

[TABLE]

where again $m$ is the smallest index s.t. $T_{m}>t$.

Recall that for Weibull distribution,

[TABLE]

Appendix C Power spectrum of a Weibull distribution

For a Weibull distribution of inter-burst intervals, the power spectrum on temporal scales much smaller than the mean burst separation and on scales much larger than this separation is flat. Here we derive the ratio of power between these two regimes.

Consider the time elapsed until $N$ successive events have occurred, $\Delta_{N}$. For large $N$ this will correspond to large time-scales and thus small frequencies. Since the inter-event intervals are considered to be independent, the central limit theorem holds and the distribution of $\Delta_{N}$ is Gaussian for large $N$. The mean and variance of this Gaussian are straightforwardly calculated from the mean and variance of an individual interval. The mean is

[TABLE]

both for the Poissonian and the Weibull case. The variance is

[TABLE]

for the Poissonian case and

[TABLE]

for the Weibull case.

Now consider a fixed time-scale $\Delta$. For large $\Delta$, the number of bursts that occur during this time, $N_{\Delta}$ will again be approximately Gaussian distributed with a mean

[TABLE]

and a comparatively small variance, which is proportional to the variance of the duration of $\bar{N}_{\Delta}$ events, given above. We thus obtain as the ratio of the large-scale power in the Weibull case and the Poissonian power the factor

[TABLE]

At high frequencies, the Poisson noise subtracted power spectrum scales as $\propto\nu^{-k}$. The numerical result shown in Fig. 6 confirms the validity of this argument.

Note that the line of reasoning followed here is applicable more generally than just to the case of a Weibull distribution. For any distribution of independent intervals for which the central limit theorem has validity, i.e., any distribution with finite non-zero moments, a similar result can be derived.

Acknowledgements

We thank Christina Peters and Liam Connor for helpful discussions and comments on the manuscript. This research has made use of NASA’s Astrophysics Data System. The figures were produced using the matplotlib library (Hunter, 2007). Some of the results have been derived using the NIFTY package (Selig et al., 2013). We acknowledge NSERC support.