The arrival time and energy of FRBs traverse the time-energy bivariate space like a Brownian motion

The origin of fast radio bursts (FRBs), the brightest cosmic explosion in radio bands, remains unknown. We introduce here a novel method for a comprehensive analysis of active FRBs' behaviors in the time-energy domain. Using ``Pincus Index'' and ``Maximum Lyapunov Exponent'', we were able to quantify the randomness and chaoticity, respectively, of the bursting events and put FRBs in the context of common transient physical phenomena, such as pulsar, earthquakes, and solar flares. In the bivariate time-energy domain, repeated FRB bursts' behaviors deviate significantly (more random, less chaotic) from pulsars, earthquakes, and solar flares. The waiting times between FRB bursts and the corresponding energy changes exhibit no correlation and remain unpredictable, suggesting that the emission of FRBs does not exhibit the time and energy clustering observed in seismic events. The pronounced stochasticity may arise from a singular source with high entropy or the combination of diverse emission mechanisms/sites. Consequently, our methodology serves as a pragmatic tool for illustrating the congruities and distinctions among diverse physical processes.


Introduction
Fast radio bursts (FRBs) are intense pulses of radio emission that last just a few milliseconds.First discovered in 2007 [1], FRBs have since been observed by a variety of radio telescopes around the world.Despite their ubiquity, however, the origin of these mysterious signals remains unknown.FRBs have a wealth of observational parameters carrying information about the sources and the propagation paths, including arrival time, energy, duration, bandwidth, polarization, dispersion, scintillation, scattering, etc.Recent studies on the polarization of FRBs suggest that FRBs are located in complex magnetized environments [2,3,4,5].These studies demonstrate information about the propagation paths of FRBs.Time and energy are two observational parameters directly related to the radiation nature of FRBs.The observation of FRB 20121102A with the Five-hundredmeter Aperture Spherical radio Telescope (FAST) for the first time revealed a bimodal distribution of FRB energy, suggesting that there may be different radiation mechanisms for FRBs [6,7].
The study of time and energy sequences of repeating FRB bursts enables a deeper understanding of the origin of FRBs, which is difficult to achieve with non-repeating bursts because their only burst is just a point in time-energy phase space.Here, we use two active repeating FRBs as our analysis objects, which are the only two FRBs known to have associated compact persistent radio sources (PRSs) [8,9,10], namely FRB 20121102A and FRB 20190520B.Both have gone through highly active episodes over a substantial dynamic range of time scales (milliseconds to months), as detected by FAST [6,10].Due to the high sensitivity and high cadence coverage of FAST's observations of these two FRBs, we can obtain more complete event series than other telescopes.We report here a first systematic examination of the FRB behaviors in the time-energy bivariate domain.

Data
FRB 20121102A was the first FRB found to be repeating and to be precisely localized [11,8,12].In one extremely active epoch, FAST detected 1652 pulses from FRB 20121102A between August 29 and October 29, 2019 [6].Since then, we have been carrying out regular monitoring of FRB 20121102A once every one or two months.On August 17 and 23, 2020, FAST caught another 12 bursts.Since then, no more burst has been detected.FRB 20190520B is the first repeating FRB discovered by FAST during the Commensal Radio Astronomy FAST Survey (CRAFTS) [13], which has been accurately located with a PRS [10].Since its discovery on May 20th, 2019, more than 200 pulses have been detected by FAST and Parkes [5].Fig. 1 shows the observation coverage and the detected bursts of FRB 20121102A and FRB 20190520B.

Time domain analysis and results
We employed the Lomb-Scargle Periodogram (LSP) [14,15] to search for periods of FRB 20121102A and FRB 20190520B, yet no statistically significant period signals were identified within the range of 1 ms to 1000 s.To quantify the significance (signal-to-noise ratio, S/N) of signals at different periods, we compared the power differences between randomly generated time series and the time series of FRB bursts in the LSP.For FRB 20121102A, we generated a time series with the same number as the FRB bursts between the time of the first and last burst's arrival time, drawn from a uniform distribution.Subsequently, we applied the LSP to this randomly generated sequence.We repeated this process 100 times, resulting in a distribution of periodogram powers for random signals at different periods.We conducted the same procedure for FRB 20190520B.Within the range of 1 ms to 1000 s, no significant periodic signals for both FRBs surpassed 5 of the distribution of periodogram powers from the random time series (Fig. 2).Thus we have significantly constrained the possibility of periodic signal existence.Several effects, including variable emission altitudes/sites, may complicate the inference of a magnetar's spin period from the arrival times of FRB bursts [16].An alternative interpretation could be that the period of the magnetar powering FRB is much longer than that of Galactic magnetars, so there is not short-timescale periodicity in the FRB data [17].Regardless of the specific cause, the absence of a detected periodic signal still highlights the inherent randomness in the emission of FRBs.We now look into the more fundamental aspects of FRBs' time-domain behaviors through analysis of waiting times.The waiting time between two events is  =  +1 −   , where  +1 and   are the arrival times for the ( + 1) ℎ and  ℎ events, respectively.
Since the discovery of FRB 20200428, the origin of FRBs being from magnetars has gradually gained popularity [18,19].As a possible source model for FRBs, the production of FRBs by magnetars relies on some trigger mechanisms [20], including crust cracking at the neutron star surface [21,22,7], sudden magnetic reconnection events in the magnetosphere [23], or triggers from external events [24,25].Motivated by this, we compared FRBs with earthquake and solar flare with similar but not identical mechanisms.
The seismic data is from Southern California from the Southern California Earthquake Data Center [26], which contains information such as the occurrence time, latitude, longitude, and magnitude of the earthquake since 1932.
For earthquake, we selected events within a region of 2 • × 2 • as a continuous seismic sequence.The magnitudes of all earthquakes are converted into energy (erg) through an empirical relation [27].Solar flare data is from the Hinode Flare Catalogue [28].We simulated a 100-step Brownian motion.In mathematics, Brownian motion is described by the Wiener process, a stochastic process  () concerning time .According to the definition of Brownian motion, for time  and , the increments of Brownian motion  () −  () follow a normal distribution  (0,  − ).Therefore, our simulation proceeded as follows: firstly, we sampled the waiting time  for each step of the Brownian motion from the exponential distribution  () =  − , where  = 1, corresponding to a Poisson process with a event rate of 1. Subsequently, based on each step's waiting time, we sampled the step size from the normal distribution  (0, ).
The waiting time distributions of these phenomena are presented in Fig. 3a-e.We defined an average time as T = ∑   ∕ ∑   , where   and   are the length and event count of the observation session.All waiting times are normalized by corresponding T .For fitting the waiting time distribution, we employed the EMCEE package to perform maximum likelihood estimation on the fitting parameters.We define the likelihood function as where  is the waiting times,  (, ) is the model for waiting time distribution, and  is the free parameter(s) for the model. (, ) could be ) −1  −(∕)  where  represents the shape parameter  and the scale parameter .
where  represents the expected value  and standard deviation  of the variable's natural logarithm.
where  represents the event rate .
where  representing two event rates  1,2 and a proportionality factor .
The waiting time distribution for a Poisson process follows an exponential distribution.The waiting time distribution of the earthquake deviates from the exponential distribution and can be better described by a Weibull distribution.This means that earthquake events with shorter waiting times occur more frequently than expected from a simple Poisson process, i.e. clustering in time.The waiting time distribution  of solar flare can be fitted by a log-normal function, which means the event rate of solar flare changes randomly.
Regarding the waiting time of FRBs, the earlier assumption was a non-stationary Poisson process or a Weibull distribution [29,30], albeit influenced by the sparse sampling of early burst events.As depicted in Fig. 3, the bimodal waiting time distribution of FRBs evidently differs from a single Weibull distribution, which can be well applied to earthquake events.There are also diverse models used to describe the distribution of waiting times for FRBs, including LogNormal distribution [31], Weibull [32], time-dependent Poisson processes [33,34].
For comparative analysis, using FRB 20121102A as an example, we selected waiting times exceeding 1 second.Fitting was undertaken with exponential (corresponding to a Poisson process), Weibull, and lognormal functions.The goodness of fit was assessed using adjusted R-squared, which provides a refined assessment of the model's fitting effectiveness by accounting for the number of independent variables utilized.The parameter counts for exponential, Weibull, and lognormal are 1, 2, and 2, respectively.The calculated adjusted R-squared values were 0.920, 0.917, and 0.906 for these three functions.While all three functions appear adept at modeling FRB waiting times, we lean toward the simplest model, the exponential (Poisson process) function, for describing the distribution of FRB waiting times.The time-dependent Poisson process, which is the sum of multiple Poisson processes, entails a greater number of parameters and evidently excels in fitting residuals compared to models with fewer parameters.Nevertheless, under the simplest conditions, our results indicate the effectiveness of a single Poisson process.Consequently, in this context, we opt for employing two Poisson processes to fit the entire distribution of FRB waiting times.
We also explored the distribution of waiting times under different energy thresholds.Taking FRB 20121102A as an example, we set a series of energy thresholds, calculated the waiting time distribution for each, and employed an exponential function to fit the two peaks of waiting times.For exponential distributed waiting times,  () =  − d =  −  d log , the peak corresponds to the position where the derivative is zero  ′ () = 0 ⇒  = 1∕.Thus, the exponential distribution provides the typical waiting time corresponding to the peak.In Fig. 3g-h, we present the characteristic times corresponding to the two peaks of waiting time under different energy thresholds.The left peak (short timescale waiting times) is insensitive to the choice of energy threshold, implying an intrinsic correlation with the radiation mechanism.The right peak (long timescale waiting times) strongly depends on the energy threshold or the detected event number.As the energy threshold increases, it leads to a decrease in the detected bursts, resulting in longer waiting times.In panel H, two inflection points correspond to the two peaks in the energy distribution.After normalizing the right peak's waiting time using the detected number of bursts under different energy thresholds and observation duration, we can see that all the normalized right peak hovers around 1, indicating that the right peak is a stochastic process relying on sampling (event rate).When the waiting time of FRBs can be described by stochastic processes, it further disfavors the prospect of (quasi)-periodicity in FRBs.
The distinct characteristics of the waiting time distribution underscore that the emission of FRBs is not like seismic events.A recent study has proposed a semblance of aftershock characteristics between FRBs and earthquakes through a custom correlation function [35].However, the similarity in the "correlation function" between earthquakes and FRBs arises from the bias of waiting time distributions from one Poisson process.Here, we show that despite both earthquakes and FRBs deviating from a Poisson process in terms of waiting times, the manner of deviation differs.Therefore, it is hard to conclude that FRBs exhibit seismic aftershock characteristics.This assertion is further substantiated in our subsequent analysis of energy distribution.

Energy domain analysis and results
Despite the common perception of randomly occurring earthquakes, seismic events cluster in time and magnitude [36,37,38].The deviation of earthquakes from randomness origins from its nonlinear dynamical systems [36,39].Assuming in an event series, the subsequent event often possesses lower energy than its predecessor.After shuffling, the probability of subsequent energy being lower than the predecessor will decrease.Based on this fact, we adopted the analysis method used in [38] for earthquakes.
Consider conditional probability where Δ  =  +1 −   and Δ  =  +1 −   are the temporal and magnitude difference between subsequent events, (⋯) is the number of pairs of subsequent events satisfying the conditions specified in parentheses.
With fixed  0 and  0 , the quantity  (Δ *  <  0 |Δ  <  0 ) for independent reshuffled series are calculated, where Δ *  =  * +1 − *  is the magnitude difference between subsequent events of reshuffled series.We employed NUMPY.RANDOM .SHUFFLE to rearrange our event series, which is based on the Fisher-Yates shuffle algorithm to produce an unbiased permutation.The central limit theorem makes the quantity  (Δ *  <  0 |Δ  <  0 ) a Gaussian distribution with mean value ( 0 ,  0 ) and standard deviation ( 0 ,  0 ).The biases between ( 0 ,  0 ) and  (Δ  <  0 |Δ  <  0 ) for earthquake, FRB 20121102A and FRB 20190520B are 5.65, 0.79 and 0.40 (Fig. 4a-c), respectively.Generally, the detection of earthquakes is more complete than that of FRBs.Here, we conducted a test related to earthquake magnitude thresholds.For FRB 20121102A, the burst with the highest energy differs by a factor of 4000 from the one with the lowest energy, and for FRB 20190520B, the factor is 500.In the earthquake data, the highest magnitude is 7.5.By progressively increasing the energy threshold, we made the lowest-energy earthquake event reach 1/500 of the highest-energy event.Throughout this process, we calculated the S/N of the hypothesis "the energy of the subsequent event is lower than that of the preceding event" for different thresholds.In Fig. 4d, we can see that even when Figure 4: The probability (P) of the hypothesis "the energy of the subsequent event is lower than that of the preceding event".Panel a-c show this probability (P) of the original event sequences (dashed lines) for FRB 20121102A (blue), FRB 20190520B (red), and earthquakes (green).The histogram and the dotted lines show the distribution and median of this probability (P) of randomly shuffled event sequences.Panel d illustrates the variation in the S/N of the probability (P) under different earthquake magnitude thresholds.
compressing the dynamic range of earthquake energies to 1/500, the S/N still exceeds 2, far surpassing the two FRBs.This indicates that no clustering in the energy of FRB bursts.Simply put, there is no excess number of bursts following a bright one, while there are more after-quakes following a major event than pre-quakes.

Stochasticity and Chaos
In dynamical studies, chaos and stochasticity are two distinct concepts.Chaos is characterized by unpredictability that increases with time, whereas stochasticity's unpredictability remains stable over time.We use "Pincus Index" ( ) and "Maximum Lyapunov exponent" ( ), respectively, to quantify the stochasticity and chaos of event sequences.
The  is used to describe the degree of stochasticity based on the Max Approximate Entropy ( ) [40,41], by measuring the change in information entropy before and after shuffling a sequence. is zero for completely ordered systems and one for totally random systems.For an event series {  } =1,⋯, , the  [40,41] can be defined as where  = {  } =1,⋯,− = {{ 1 , ⋯ ,  1+ }, ⋯ , { − , ⋯ ,   }} is the reorganized sequence of the initial event series,  is the length of the initial event series, dist(  ,   ) is the distance between each pair of the reorganized sequence , and  is the distance threshold. quantifies the greatest difference in information between the segments of length  and  + 1, which needs to obtain a max entropy of the recombination sequence by changing different thresholds .
For making the  comparable between series, the  was defined as After calculating the  or called  initial for initial event series, we randomly reshuffled the series and calculated  for 100 times.The median value is  shuf f led , the standard deviation value is used to define the error of  .In the calculation of the  , the distance metric used is the Euclidean distance.To preserve the original sequence information as much as possible and minimize the introduction of artificial bias, we linearly mapped both the time and energy sequences to the range of 0 to 1 simply.The  represents the degree of dispersion of trajectories in phase space and is a numerical characteristic used to identify chaotic behavior in a nonlinear system [42]. less than 0 corresponds to a periodic motion or a stable system that is static in time-energy space. greater than 0 indicates the existence of chaos.We use NOLDS [43] method for  calculation.As the  is the maximum value of the whole spectrum of Lyapunov exponents, it is difficult to define an error.
The  and  values are shown in Fig. 5a.The pulsar data used here was a pulsar named J1840+2843 discovered within the CRAFTS project 1 .As expected, the most regular motion is a sine-wave, with small values in both, followed by pulsar with representative  ∼0.4 and  ∼0.08.Earthquakes are most chaotic with  ∼0.33.Solar flares are less chaotic, but more random than earthquakes.FRBs are even more random ( ∼0.84-0.97)and less chaotic ( ∼0.07-0.09)than solar flares, mimicking Brownian Motion ( ∼0.99) that is the most random in the bivariate time-energy space, among the systems considered here.
To ascertain the robustness of our calculations, we considered two types of tests.In the calculation of the  , only the parameter  remains as a free variable.We conducted tests by adopting  ∈ {2, 3, 4, 5}.As depicted in Fig. 6a, the deviation of the  did not exceed 0.1 for different values of .Even with slight variations in the  , the relative relationships among different physical phenomena remained unchanged when the same  was selected.Additionally, given that the detection of these physical phenomena cannot The top and middle panels present event series in time-energy space of these sources.The color changes from blue to red, implying increased stochasticity.( ,  ) value is (0.24, 0.06) for sine-wave, (0.40, 0.08) for pulsar, (0.72, 0.33) for Earthquake, (0.81, 0.21) for solar flare, (0.84, 0.09) for FRB 20121102A, (0.97, 0.07) for FRB 20190520B, and (0.99, 0.05) for Brownian Motion.
be exhaustive, we further examined the impact of sequence completeness on the results of the  and  calculations.By adjusting the energy threshold, we retained different proportions of events (ranging from 100% down to 50%) and subsequently computed the  and  for each case.From retaining 100% of the events to selecting only those with the top 50% energy levels, the  and  exhibited remarkable stability (Fig. 6b-c), confirming the robustness of our calculations.
Within the stochasticity-chaos phase space, FRBs are noticeably distinct from earthquakes, exhibiting lower levels of chaos and higher levels of stochasticity.Additionally, as previously analyzed, earthquakes exhibit time and energy clustering behavior (Figs.3-4).FRBs, however, do not show clustering in neither time nor energy.The random paths of FRB events in time-energy space favor origin models based on complex processes, such as in a thermodynamic system, the movement or diffusion of molecules is akin to the phenomenon of interest.

Conclusion
In summary, using  and  to quantify the stochasticity and chaos of event sequences may be an effective method for intuitively demonstrating the similarities and differences between various physical processes.Specifically, here we compare FRBs, pulsar, solar flares, earthquakes, and Brownian motion in the stochasticity-chaos phase space.The strong stochasticity, akin to Brownian motions, of active repeaters with PRS counterparts, along with the growing evidence of their multi-variate behaviors, such as their bimodal energy distribution [6,7,44,45], could be generated by a single source with high information entropy or the combination of multiple radiation mechanisms or emission sites.Either way, it is unlikely that active repeating FRBs originate from a stably spin compact objects within a clean environment.

Figure 1 :
Figure 1: The detected bursts with energy from FRB 20121102A and FRB 20190520B.In panel a, the blue dots and red dots show the energy and arrival time of bursts from the two FRBs.The blue and red bars are the observation session of the two FRBs with detection, the yellow bars are the observation session without detection.Panel b shows the kernel density estimate (KDE) energy distribution of the two FRBs.
Normalized right peak (long) waiting time

Figure 3 :
Figure 3: The normalized waiting time distribution of earthquake, solar flare, two FRBs, and Brownian motion.Colored lines in panel a-e are the KDEs of waiting times for these phenomena, which are normalized by T of 3435.1 s, 25810.4s, 88.4 s, 123.3 s, and 190.4 s, respectively.Gray regions in panel a-e indicate the fitting of waiting times with Weibull function, LogNormal function, and Exponential function (Poisson process).The gray solid lines denote the short characteristic waiting time of the two FRBs, which are ∼ 4 ms for FRB 20121102A and ∼ 36 ms for FRB 20190520B.Panel f shows the distribution of waiting times longer than 1 s for FRB 20121102A.The solid line, dashed line, and dash-dot line respectively denote the fitting results of the exponential function, Weibull function, and lognormal function.Panel g-i show the left peak, right peak and normalized right peak of waiting time under different energy thresholds.

Figure 6 :
Figure 6: The robustness test of  and  .Panel a shows the  values using different segments of length .Panel b-c show the  and  under different energy thresholds.