Evidence for Extended Charging Periods Prior to Terrestrial Gamma Ray Flashes

We present an analysis of lightning interflash intervals in 219 terrestrial gamma ray flash (TGF) producing thunderstorms. Clustering was used to identify groups of lightning sferics, interpreted as individual thunderstorms, in combined World Wide Lightning Location Network and Earth Networks Total Lightning Network data. In these individual groups of sferics, analysis was done on the lightning flash frequency within ±10 min of the Fermi recorded TGF. We find that typical interflash intervals immediately prior to TGFs are 24% longer than mean interflash intervals in their individual producing storms, while the interflash intervals immediately following the TGFs are typically 8% shorter than normal. The significance of these results, tested using a numerical bootstrap method, was found to be highly significant for the pre‐TGF interval. These results could imply that a stronger electric field is necessary for the production of TGFs and may help to explain why some lightning strikes produce TGFs while others do not.


Introduction
Terrestrial gamma ray flashes (TGFs) are submillisecond, intense bursts of photons with energies often reaching tens of megaelectron volts . TGFs have been associated with thunderstorms since their discovery in 1994 (Cummer et al., 2005;Fishman et al., 1994. Observations of TGFs have been made primarily from spacecraft such as the Burst and Transient Source Experiment on board the Compton Gamma Ray Observatory (Fishman et al., 1994), the Reuven Ramaty High Energy Solar Spectroscopic Imager , the Gamma-Ray Burst Monitor on the Fermi Gamma-ray Space Telescope (Briggs et al., 2010), and by the Astrorivelatore Gamma a Immagini Leggero satellite (Marisaldi et al., 2010). However, TGFs have also been observed, in a few cases, from the ground (Bowers et al., 2017;Dwyer et al., 2004;Dwyer, Schaal, et al., 2012;Enoto et al., 2017;Hare et al., 2016) and by aircraft Smith et al., 2011).
In the past, the Relativistic Runaway Electron Avalanche (RREA) mechanism (Gurevich et al., 1992) acting on cosmic ray seed particles was used to explain TGFs (Lehtinen et al., 1996). However, Dwyer (2008) showed that this mechanism, on its own, was insufficient in explaining the production of TGFs because it would require unrealistic avalanche lateral sizes, unrealistic avalanche multiplication factors, or unrealistic seed fluxes to recreate the observed fluence of TGFs. Other viable explanations include the relativistic feedback mechanism (Dwyer, 2003) and the lightning leader tip mechanism (Carlson et al., 2010;Celestin & Pasko, 2011;Moss et al., 2006). The relativistic feedback model accounts for bremsstrahlung X-rays created by initial electron avalanches. These X-rays may create secondary avalanches either through Compton backscattering with other electrons or through pair production, in which backward propagating positrons may interact with other electrons, further seeding RREA (Dwyer, Schaal, et al., 2012). This model places a fundamental limit on the electric field strength, above which the feedback mechanism dominates and discharges the electric field. The lightning leader tip mechanism allows for low-energy electrons to be accelerated to thermal runaway electrons, due to extremely high local electric fields near lightning leader tips. This process is known as cold runaway (Gurevich, 1961). These runaway electrons can further act as seed particles for RREA in the ambient electric field, thus eliminating the need for external seeds (Dwyer, Smith, et al., 2012).
In this paper we examine the interflash intervals between lightning flashes in TGF-producing thunderstorms. We look specifically at how the flash intervals immediately prior to and after the TGF compare to typical interflash intervals within the same individual storms. We interpret these results in the context of the ambient electric field charge up time and the idea that longer interflash intervals (and thus longer charging times) result in stronger or more extended electric fields before the TGF associated lightning discharge. Higher electric fields may be an important factor in the production of TGFs as suggested by Smith et al. (2018) who observed two similar upward lightning flashes, one of which produced a TGF, while the other did not. The TGF-producing flash occurred during a gamma ray glow (a state of continuous, low-level RREA; see Tsuchiya et al., 2011, which showed that glows exhibit the expected spectra due to RREA), implying that the electric field was higher in that case than for the non-TGF-producing flash. This evidence suggests that an enhanced electric field may be required for the production of TGFs. The conclusion of this paper, that TGFs occur after relatively extended charging periods, is important because, to our knowledge, it is the first study to show evidence of an enhanced thunderstorm field driving TGFs in a statistically rigorous way. This is also the first observation that would allow one to guess, in the absence TGF signatures, which flashes in a storm are more likely to produce TGFs.

Data and Methods
In this analysis we use TGFs identified by Fermi Gamma-Ray Burst Monitor. Located Radio atmospheric signal pulses (hereafter referred to as "sferics"), supplied by the World Wide Lightning Location Network (WWLLN), within ±10 min of identified TGFs are available in Fermi's online TGF catalog (Roberts et al., 2018). In addition to data from WWLLN, lightning sferic data have also been provided by the Earth Networks Total Lightning Network (ENTLN).
There are 1,339 TGFs for which positive WWLLN associations have been identified. Of this number, 1,169 TGFs had both ENTLN and WWLLN data available. This subsample of TGFs was analyzed in this work. A time correction corresponding to a photon propagation time to Fermi (at an altitude of roughly 555 km) from a source altitude of 20 km and a location of the TGF associated sferic is first applied to all associated ENTLN sferics (Briggs et al., 2010). The WWLLN sferics used have already been time corrected into Fermi's observing frame, so this time correction is necessary to put both WWLLN and ENTLN sferics in the same time frame. This time correction ranges from 1.78 to 3.21 ms for the TGFs used in this study. The time corrected ENTLN sferics are then merged with available WWLLN sferics in chronological order. The merged data set is then spatially clustered using the Hierarchical Density Based Spatial Clustering of Applications with Noise (HDBSCAN) algorithm (McInnes et al., 2017). This is done to identify lightning belonging to individual thunderstorms. The HDBSCAN clustering algorithm requires neither a user-specified scale size nor a specified number of clusters. HDBSCAN works by performing the DBSCAN algorithm over a range of varying distance parameters and selecting the most persistent clusters, that is, the ones that remain stable as the length scale decreases (Campello et al., 2013). Thus HDBSCAN is better suited to this application than typical clustering algorithms such as DBSCAN and k means (Ester et al., 1996;Wu et al., 2008), as it allows for clusters of varying sizes and arbitrary shapes to be identified and the number of clusters does not have to be known or guessed beforehand.
HDBSCAN assigns a probability score (0.0 ≤ p ≤ 1.0) to each point in a given cluster. Points with lower scores are not as confidently identified as belonging to their specified cluster, in contrast with points having higher probability scores. For this work, we set this probability score threshold at p ≥ 0.2. This threshold was chosen after manually examining cluster assignments over a range of thresholds.Threshold values larger than 0.2 have the effect of including only the inner cores of clusters while excluding the apparent edge points of clusters. Threshold values smaller than 0.2, however, tend to include obvious outlier points which do not appear to be a part of the cluster they are assigned to. A threshold of 0.2 seems to be a good choice as it still excludes obvious outliers while retaining the inclusion of the main bodies of clusters. As an extra measure, this parameter was varied over a range of 0.0 to 0.99, and the results of this study do not strongly depend on this threshold value. In addition, HDBSCAN classifies some sferics, typically isolated sferics or those which have ambiguous cluster assignments, as noise. Noise classified sferics are excluded from this analysis. Figure 1 shows an example of this spatial clustering. Land has been colored beige for reference, while oceans have been left white. (left) Example of a "good" or usable cluster in orange (containing the TGF associated sferic). (middle) Example of a clustering that was unusable for this analysis due to its poor isolation. (right) Example of rejected cluster, based on its dispersiveness. TGF = terrestrial gamma ray flash.
In addition to spatial clustering, the sferics associated with the TGF-producing cluster were also temporally grouped using a Gaussian kernel density estimation with a bandwidth of 0.25 s. This bandwidth parameter controls the width of the gaussian kernel used. Too large of a bandwidth would result in the grouping together of multiple distinct flashes, while too small of a bandwith would result in each sferic being grouped into its own flash. A choice of bandwidth within the range of 0.2-0.7 produces consistent results and allows for individual sferics to be grouped into distinct lightning flashes, with total durations typically less than 1 s, a common value for typical lightning flashes (Dwyer & Uman, 2014). This flash grouping also eliminates any double counting from merging the ENTLN and WWLLN sferics by merging said double counts into the same flash.
After spatially and temporally grouping sferics for all 1,169 events, the full set was then visually reexamined and narrowed down to include only those events with accurate cluster assignments, as determined by two different researchers. The accuracy of cluster assignment was evaluated using the criteria of density, isolation, and distinctiveness. A cluster was considered "good" if it was distinguishable (i.e., not overlapping with other clusters), not clearly made up of several subclusters, and not so close to other clusters that individual sferics might have been included from those nearby clusters (i.e., at least separated by more than the WWLLN location accuracy or distinct, well-separated areas of sferic densities are clearly visible). In addition, we required the clusters to have at least 15 distinct lightning flashes and a total extent (farthest distance between any two points in a cluster) of less than 60 km. We required at least 15 distinct lightning flashes to ensure the statistical quality of the mean interflash intervals, as this statistic tends to become very imprecise with a small number of samples. This analysis was also performed on clusters with down to 10 flashes and produced similar results. The criterion of 60 km was chosen to exclude the use of clusters which were made up of several obvious subclusters, but which HDBSCAN failed to identify, and to exclude very dispersed clusters which were given questionable cluster assignments by HDBSCAN (see Figure 1 [right] for example). Typical thunderstorm diameters can reach up to roughly 100 km with convection regions reaching up to a few tens of kilometers in diameter (Houze, 2004;Klimov et al., 2018). This, considered with WWLLN's location accuracy, makes the 60-km threshold a conservative choice. This process of selecting clusters was entirely blind to the temporal characteristics of the clusters so as to avoid any bias in terms of the parameter being studied (interflash intervals). In total, 219 TGF-producing clusters met these standards and were subsequently included in this study. 475 TGF-producing clusters did not have enough recorded flashes to be statistically reliable, 162 clusters were too large or dispersed to identify a single storm, 175 clusters were excluded based on isolation or distinctiveness, 13 clusters were associated with double-pulsed TGFs which were not considered in this study to avoid any systematic differences that may exist between single and double-pulsed TGFs, and in 125 cases the clustering algorithm failed to attribute the TGF strike to any cluster. For the 219 selected TGF-producing clusters, time separations between flashes in the cluster were calculated. Specifically, differences between the terminating sferic in a given flash and the initiating sferic of the next flash were calculated for all flashes. Calculating the time difference in this way results in the minimum possible interflash interval, as opposed to calculating the time differences from the start (or end) of every flash, for example. For all 219 events, we kept track of the time differences between the TGF-producing flash and the flashes immediately prior to and immediately after the TGF, referred to here as the "pre-TGF" and "post-TGF" flashes, respectively. In addition, we also kept track of the mean interflash interval of all flashes in the 20-min data sample of a given cluster, referred to here as the typical interflash interval. This analysis has also been done using the median interflash intervals, with very similar results. Figure 2 shows these interflash intervals for an example TGF-producing storm, with all important parameters labeled.
Flash rates in individual thunderstorms can vary dramatically depending on the specific environment in which they develop (Ushiro et al., 2001;Williams et al., 2005). For this reason, we compare the pre-TGF and post-TGF intervals to the typical interflash interval of their individual producing storms. This is done by dividing the pre-TGF and post-TGF intervals by the typical interflash interval of their individual storms. The result is a ratio which is then used to characterize the intervals. Comparing the pre-TGF, post-TGF, and typical interval distributions of all storms directly to each other, without dividing first, would fail to consider this flash rate variation as storms with higher charging rates will necessarily have shorter interflash intervals than storms with lower flash rates. It is impossible to determine if the pre-TGF or post-TGF intervals are significantly different by comparing storms in this way. Figure 3 shows the results of the above analysis. We find a mean ratio value of 1.24 for pre-TGF intervals and 0.92 for the post-TGF intervals. These results imply that the pre-flash pre-TGF charging period is typically 24% longer than the average charging rate in TGF producing storms, while the post-TGF charging rate is typically 8% shorter than normal. We also note a distinct difference between the pre-TGF and post-TGF ratio histograms. In particular, we see an apparent lack of TGFs where the pre-TGF interval is shorter than normal (a ratio of less than 1) as compared to the post-TGF intervals.

Results
To test the significance of the 24% and 8% deviations from the typical interflash intervals in both the pre-TGF and post-TGF distributions, we have used a numerical bootstrap method. For each of the 219 chosen clusters, we randomly select a flash to be the "TGF." We then calculate the corresponding pre-TGF and post-TGF intervals and take the ratio to the typical interflash interval for that storm. We then calculate the mean value of these ratios for both the pre-TGF and post-TGF intervals over the entire 219 TGFs. This is then repeated 5 million times for each of the TGFs to get sufficient random sampling. The results of this bootstrapping, shown in Figure 4, indicate that the 24% longer than normal interval immediately prior to the TGF seen in Figure 3 is extremely unlikely to have occurred by chance, as only twice in 5 million random trials was this statistic larger than 1.24. Conversely, the 8% shorter intervals seen immediately following the TGF are not highly significant, as the random scrambling produced a value lower than 0.92, roughly 6% of the time.
This significance test was also performed on 317 non-TGF-producing clusters which satisfied the above criteria (and which were within 500 km of the TGF-producing storms to reduce any systematic geographical variations). Two hundred nineteen storms were randomly selected from this sample of 317, and the numerical bootstrap method was performed on that subset of storms in the same manner as above. From each storm, a random pre-TGF and post-TGF interval was selected and divided by the mean interflash interval of that storm. After the 219 interval ratios were selected, the mean of that distribution was calculated. This process was repeated 5 million times, producing a distribution similar to that in Figure 4. To check the dependence of this statistic on storm selection, this entire process was repeated for 1,000 different sets of 219 storms, chosen from the 317 that met the aforementioned criteria. Because each of the 1,000 trials were composed of a different set of storms, the number of mean ratios larger than 1.24 varied in each trial. However, of all 1,000 different sets of storms, the maximum number of times the 1.24 statistic was exceeded was only 71 out of 5 million random trials (0.001%). For that same set of storms, values lower than 0.92 were found roughly 10% of the time. This bootstrapping represents the worst case scenario and shows that our finding of a 24% longer charging period prior to TGFs is very significant.
The offset between the two curves in Figure 4 is a real effect caused by the overall decline of the flash rates in these storms. TGFs are known to occur while the flash rate of a storm is declining (Smith et al., 2010), and so post-TGF intervals in random samples will tend to be slightly longer than pre-TGF intervals.

Conclusions
In this work we have examined lightning interflash intervals in 219 TGF producing thunderstorms. We find that the flash interval immediately prior to the TGF is longer than the typical interflash interval by an average of 24%. Conversely, we find that the interval immediately after the TGF is shorter than normal by a typical value of 8%. The significance of these results was tested using a numerical bootstrap method. The extended pre-TGF interval was found to be very significant, as only 0.00004% of random trials produced a value larger than 1.24. However, in the case of the post-TGF interval, values with means shorter than the actual value were found, by random, in 6% of the trials. The post-TGF interval is thus consistent with being a "typical" flash.
There are several influences which could produce such extended pre-TGF intervals. For example, lightning flash rates are known to be strongly correlated to factors such as updraft velocity and content, cloud base height, and concentrations of cloud condensation nuclei (MacGorman et al., 2010;Williams et al., 2005;Yoshida et al., 2017). It is possible that any of these factors could affect flash rates by affecting charging rates; however, these factors likely all vary on time scales much longer than typical interflash intervals. Furthermore, it is difficult to imagine why variations in these factors would influence the production of TGFs.
Longer interflash intervals prior to the TGF would imply that the electric field had more time to charge up before a discharge, resulting in stronger or more spatially extended electric fields. As a simple model, if one imagines the charge structure in a thunderstorm to be similar to that of a capacitor, for which the volume remains roughly constant over the span of tens of seconds, a shorter interflash interval would limit the field to a lower value, while a longer interval would allow the field to grow to relatively stronger values before lightning discharges. The assumption of a roughly constant charge structure over tens of seconds is supported by the results of Pawar and Kamra (2004) who demonstrate the large spatial and temporal scales of the lower positive charge center in a thunderstorm. This mental picture is also supported by the recovery curves of the electric field between lightning flashes as demonstrated in Pawar and Kamra (2002; see their Figure 2), which show a consistent charging profile between flashes and that the field at the time of each flash depends on the interflash interval.
It is also well established that low-flash-rate storms produce more energetic lightning flashes with higher peak currents (Chronis et al., 2015;Hutchins et al., 2013), meaning that the electric field must either be enhanced or strong over a larger volume than for less energetic flashes. It is also not unexpected that higher flash rates limit the growth of the background electric field (Chronis & Koshak, 2017). An enhanced electric field makes it more likely that a leader channel could push the field well above the feedback threshold over a significant volume (Dwyer, 2008), resulting in more favorable conditions for the production of TGFs. A stronger field could explain why some lightning flashes produce TGFs, while similar flashes do not .Finally, in terms of the electric field, the fact that the post-TGF interval is consistent with typical flashes would imply that the TGF returns the ambient electric field to a "normal" state, as opposed to overdischarging or underdischarging it.