This study provides a preliminary, regional assessment of the viability of using spatiotemporal lightning patterns to classify storms into single- versus multi- and supercell storm modes. Total lightning flashes (intracloud and cloud-to-ground flashes) occurring during the afternoon and evening of the period May–August 2012 within an area of the central Appalachian Mountains region were grouped based on their spatial and temporal characteristics using single-linkage clustering. The resulting discrete thunderstorm clusters were characterized in terms of duration, motion, areal extent, and shape. These values were used to formulate four individual attribute scores representing the similarity to the expected values for a typical single-cell thunderstorm. The four scores were then combined into a storm index (SI) using relative weights determined through the analytic hierarchy process (AHP) performed on input from operational forecasters. Of the study days, 89 (72.4%) possessed appreciable lightning, of which 36 (40%) possessed a defined minimum amount of lightning activity required for further analysis. These 36 storm days were divided into two tiers according to the distribution of median daily SI values. The tier containing the 24 storm days (66.7%) with the largest median SI values possessed statistically significant smaller values of 0–6-km wind shear [13.8 knots (kt; 1 kt = 0.51 m s−1)] versus the 12 days in the lower tier of SI values (26.5 kt). This consistency between the total lightning-based classification scheme and increased vertical wind shear associated with lightning-defined multi- and supercells, also evident in synoptic atmospheric composites, lends credibility to the procedure.
It has long been known that the form of convective storms varies as characterized by longevity, degree of severity, mode of propagation, intensity of rainfall (Weisman and Klemp 1982), areal extent (e.g., Austin and Houze 1972), and, therefore, their radar signature (e.g., Smith et al. 2012). Taken together, these traits have given rise to conceptual models describing a set of storm modes that includes three broad categories: single-, multi-, and supercell thunderstorms (Weisman and Klemp 1982). However, for either of these storm modes there can be substantial variation in their characteristic traits. In this respect, explicitly categorizing and defining storm mode is a difficult task, prompting the suggestion of a continuous spectrum of storm modes (Weisman and Klemp 1982, 1984).
The least organized convective mode, single-cell ordinary thunderstorms, is generally defined as a storm forming in an environment characterized by weak vertical wind shear [<10 m s−1; Markowski and Richardson (2010)]. Although known to form in atmospheres possessing only a few hundred joules per kilogram (J kg−1) of convective available potential energy (CAPE), single-cell thunderstorms capable of producing severe weather, also termed pulse storms, often require CAPE in excess of 2000 J kg−1 (Markowski and Richardson 2010). This disorganized storm mode is prevalent during the warm season across much of the eastern United States.
Multicell convection, a more organized form of ordinary convection, is composed of ordinary cells whose outflows initiate new ordinary cells with the assistance of atmospheric wind shear (Markowski and Richardson 2010). While the individual cells are ordinary, the multicell cluster can last for hours and produce hazardous weather. Supercell thunderstorms, the most organized convective mode, are a nonordinary form of convection in which strong atmospheric wind shear enables an individual cell to last for hours and obtain a mesoscale circulation not present in ordinary cells (Markowski and Richardson 2010).
Discernment of the convective mode is a complicated task given the commonalities and overlap in the different characteristics of storm structures (Weisman and Klemp 1982, 1984). Furthermore, many of the characteristics associated with the storms themselves, or the atmospheric conditions that support their development, are qualitatively defined (i.e., “low shear” and “high instability” that reference the typical pulse storm environment). While general operational values exist for these thresholds (Thompson et al. 2003), cases where the observed storm mode violates the informal expectations of what a given environment would support further complicate the classification process (Csirmaz et al. 2013).
Currently, most storm mode classification efforts utilize radar-based observations to categorize an event (e.g., Gallus et al. 2008; Thompson et al. 2007; Smith et al. 2012). However, lightning-based classification could provide enhanced automation and objectivity to a currently labor-intensive process. This study provides a preliminary assessment of the viability of using spatiotemporal total lightning patterns to differentiate single-cell storms from more organized convective modes. For the purposes of this analysis, non-single-cell storm modes will be broadly referred to as the multi-/supercell mode. While it may be argued that there exist more than three distinct storm modes, the taxonomy utilized within the study is not meant to imply otherwise; it is simply a method to distinguish single-cell storms from non-single-cell modes.
Additionally, while multi- and supercell thunderstorms are by no means equivalent modes, the spatiotemporal traits that would differentiate a single- from a multicell storm would also differentiate a single- from a supercell storm. Thus, this preliminary study will make no attempt to discern multi- from supercell thunderstorms, instead combining them into one multi-/supercell category. Following the establishment of the Geostationary Lightning Mapper (GLM) aboard the Geostationary Operational Environmental Satellite-R (GOES-R) scheduled for launch in early 2016, lightning data are likely to assume more prominent roles in weather forecasting and in the research and study that underlies improved severe weather forecasting. A greater ubiquity of lightning data in the future, and the research need for stratification of storms by mode, leads to the question of whether archived lightning data can be used in research to identify storm mode.
In this preliminary study, the classification of storms as single- versus multi-/supercell storm modes during a summer season for a region of the eastern United States was predicated on the spatiotemporal distribution of total lightning, which is the sum of intracloud (IC) and cloud-to-ground (CG) flashes, detected by the Earth Networks Total Lightning Network (ENTLN). The data analyzed in this study were gathered during the summer of 2012 over southwestern Virginia, southeastern West Virginia, and northwestern North Carolina, covering an area roughly outlining the Blacksburg, Virginia, National Weather Service Forecast Office (NWSFO) County Warning Area (CWA; Fig. 1a). The approximately 75 630 km2 area of the study domain includes the mountainous terrain of the central Appalachians to the west and lower-elevation piedmont to the east (Fig. 1b). During summer, this region of the United States is frequently characterized by warm, moist air masses that give rise to disorganized convection during the peak of daytime heating, but it is also susceptible to more dynamically forced, organized convection associated with midlatitude synoptic cyclones.
Regional lightning and atmospheric data
The ENTLN operates a network of 700 lightning detection sensors worldwide (S. Heckman, Earth Networks, 2013, personal communication). The Earth Networks Lightning Sensor (ENLS) is capable of detecting electromagnetic frequencies between 1 Hz and 12 MHz (Liu and Heckman 2011). This wide range of frequencies allows the network to locate and classify the low-frequency waveforms emitted by CG flashes as well as the high-frequency emissions of IC flashes. Waveforms detected by the ENLS are classified as either IC or CG and, subsequently, are grouped into flashes if they occur within 700 ms and 10 km of each other (Liu and Heckman 2011). As implied by this data processing sequence, lightning flashes, which are the subject of this study, are composed of one or more strokes. In the preprocessed data obtained from the ENTLN, all detected strokes had already been combined into flashes.
In order for a flash to be formally recorded by the ENTLN, between five and eight sensors must detect the flash and agree upon the time and distance from the sensor (S. Heckman 2013, personal communication), minimizing the concern that erroneous flashes are recorded in the database. Comparisons of the ENTLN to regional Lightning Mapping Arrays (LMAs; Carey et al. 2011) and the Lightning Imaging Sensor (LIS) satellite (Sloop et al. 2014), have found reasonable consistency between the data sources. The detection efficiency of the ENTLN approaches 99% for CG flashes in the eastern United States while the detection efficiency for IC flashes over the study area is roughly 70% (Liu and Heckman 2011). Following an upgrade of the detection algorithm, comparisons to rocket-triggered flashes in Florida determined a median location error of 687 m (Mallick et al. 2013).
To characterize the atmospheric environment in which the lightning flashes occurred, daily 1200 and 0000 UTC, corresponding to 0800 and 2000 local time (LT), radiosonde data were obtained for the Blacksburg, Virginia (KRNK), launch site, located at approximately 37.2°N, 80.4°W. At an elevation of 654 m, the KRNK sounding site sits in the heart of the Appalachian Mountains near the center of the domain area (Fig. 1a). Sounding data were retrieved from the University of Wyoming upper-air data distribution website (http://weather.uwyo.edu/upperair/sounding.html). Wind shear over the 0–6-km layer [given in knots (kt; 1 kt = 0.51 m s−1)] was calculated by finding the magnitude of the difference of the wind vectors at the surface and an interpolated vector 6 km above the surface. A linear interpolation was performed using the wind vectors immediately below and above the 6-km level. Several studies have found 0–6-km wind shear to provide a meaningful measure of storm organization (e.g., Thompson et al. 2003).
While the 0000 UTC sounding is advantageous as a direct observation, in some cases the sounding profile may have been influenced/contaminated by convection occurring prior to this evening balloon launch. However, precipitation data obtained from the National Climatic Data Center’s (NCDC) Climate Data Online (CDO) portal (www.ncdc.noaa.gov/cdo-web/) revealed that precipitation was recorded at KRNK within 1.5 h of the radiosonde launch on only 10.7% of the days considered. These days were allowed to remain within the dataset because determining whether the sampled atmosphere had indeed been contaminated would be very difficult. Days with sounding malfunctions were excluded from the statistical comparisons.
a. Identification of lightning clusters
The 1 692 338 ENTLN flashes recorded between 1800 and 0300 UTC (1400 and 2300 LT) each day of the period from 1 May through 31 August 2012 were isolated for further analysis. Because the major forcing mechanism for a true single-cell storm is generally buoyancy produced by daytime surface heating via insolation (Markowski and Richardson 2010), any flashes occurring outside of this time frame were excluded from consideration. The summertime afternoon–evening period of the diurnal cycle represents a period during which all storm modes were possible. If a broader seasonal or diurnal range were considered, the possibility of biasing the sample toward more organized storm modes would increase. Thus, the potential differentiating capability of the lightning data was judged to be maximized during the selected time frame. Further, the 9-h daily temporal domain was established through consultation with forecasters at the Blacksburg NWSFO.
Flashes were organized by day of occurrence and then grouped into discrete lightning clusters using the single-linkage agglomerative hierarchical clustering technique as outlined by Gong and Richman (1995). Clustering techniques, including single linkage, are a common method of analysis within the atmospheric sciences (e.g., Fovell and Fovell 1993; Kalkstein et al. 1987; Marzban et al. 2009). This specific method of clustering lightning flashes began with each flash as a separate cluster. Clusters were then iteratively combined by merging the two clusters with the smallest calculated distance separating their nearest two elements. The iterative merges continued until reaching a desired number of final clusters. The strength of the single-linkage method lies in its flexibility to form clusters of unknown and often irregular shape. Additionally, this clustering procedure, along with all agglomerative hierarchical cluster techniques, was capable of classifying any number of final clusters between N (the number of flashes in the sample) and 1. While other hierarchical clustering techniques, such as complete and average linkage, also possess this trait, they did not perform as accurately on a pilot day used to test alternative clustering techniques by comparison to radar imagery.
For this analysis, the distance measure formulated for the separating distance between each flash was proportional to the magnitude of the three-dimensional vector connecting the flashes in two-dimensional space (x, y) and time t. The x and y coordinates of a flash were chosen to be the universal transverse Mercator (UTM) easting and northing, respectively. For computation purposes the easting and northing of each flash were divided by three orders of magnitude, essentially representing the coordinates in kilometers instead of meters. For t, the original 24-h time stamp was converted to thousands of seconds since midnight (kiloseconds; ks). The coordinate was then linearly scaled by a factor of 5.556 km ks−1 in order to standardize the relative contributions of space and time to the distance vector. This value was determined by calculating the ratio of space to time for a generic, discrete thunderstorm as defined by observations from Byers and Braham (1949, hereafter BB49) and Rauber et al. (2008)—a storm lasting 1 h and attaining a 20-km cross section at its widest (20 km/3.6 ks = 5.556 km ks−1). While a 20-km diameter is large for a single-cell storm, it is an intermediate value between unorganized and organized convective modes (Rauber et al. 2008). A smaller cross section might have biased the clustering process to identify more “single cell” clusters than truly existed.
No evidence of this time–space standardization could be identified in the previous literature within the atmospheric sciences. However, statisticians have long debated the most appropriate method for standardizing variables across different sets of units (e.g., Everitt et al. 2001; Milligan and Cooper 1988). Consequently, standardization is a complex challenge and many different methods have been suggested. The standardization utilized in this analysis follows a modified form of standardization that performed better than five other tested methods in a simulation study by Milligan and Cooper (1988). The superior form tested by Milligan and Cooper (1988) involved dividing one variable by the range of values of the other. The calculation of the 5.556 km ks−1 standardization factor was computed by dividing a typical temporal value by a typical spatial value. As a preliminary study, standardization was a secondary concern to the primary research question, and this portion of the methodology is open to refinement and improvement in future research.
To optimize computing efficiency (as the iterative nature of the clustering process makes it intensive), each day’s flashes were broken into packages of approximately 4000 flashes, with each package overlapping with the one preceding it by 5 min. The single-linkage clustering technique was then performed on each package until only 150 distinct clusters remained. This threshold was selected because it represented a cautiously large upper limit on the number of distinct storms that might have occurred within a 4000-flash package. Since single-linkage clustering allows a single flash to constitute a cluster, in a worst-case scenario when only one distinct storm truly occurred, it would still be apportioned 96.3% of its deserved flashes [i.e., (4000 − 149)/4000 = 96.3%]. Days with less than 300 flashes were investigated manually.
Once again, a much deeper level of complexity belies the above methodology. The selection of the number of final clusters, 150 in this case, is a difficult and complex task that ultimately requires some element of subjectivity (Everitt et al. 2001). Since the true number of thunderstorms on a given day is unknown, it is impossible to designate the number of final clusters beforehand. Consequently, statisticians have proposed numerous methods for estimating the most appropriate number of clusters prior to analysis (Everitt et al. 2001). Thirty such procedures were evaluated by Milligan and Cooper (1985), who documented a wide-ranging ability of the procedures to accurately ascertain the true number of clusters prior to any clustering. It should be noted then that this selection criterion holds some degree of influence over the subsequent results of the analysis. However, as a preliminary study, the termination point of the clustering technique was a secondary concern within the context of this study’s objective of assessing storm mode classification. This portion of the methodology is open to refinement in future research.
After all of the packages for a given day were clustered, the overlapping flashes were analyzed to identify clusters that were interrupted in the t dimension when the package was created. For instance, if during the 5-min period of overlap, 200 flashes in package A were grouped into cluster 10A, and the same 200 flashes were grouped into cluster 48B in package B, then all flashes in cluster 48B would be added to 10A, and all the redundant flashes would be removed. At the end of the process, all unique flashes on each day were seamlessly clustered across the entire time domain. In total, 56 697 clusters were generated for the period of consideration. This 5-min overlap period yielded acceptable results on a series of test packages, but future studies might consider alternative overlap values.
b. Characterizing the temporal and spatial properties of the clusters
Once lightning clusters were generated for 1 May–31 August 2012, each cluster was characterized in terms of beginning and end time, duration, IC flash total, CG flash total, IC + CG flash total, the maximum observed total flash rate, the time of the maximum observed total flash rate, and the average velocity of the cluster in two-dimensional (x–y) space. The average velocity was calculated by dividing the distance between the mean center of all flashes occurring within the first minute of the cluster and the mean center of all flashes occurring within the last minute of the cluster by the duration of the cluster. While 1-min mean centers appeared to satisfactorily describe storm motion, longer periods of averaging might also be considered in the future. The maximum flash rate was determined following the methodology outlined by Schultz et al. (2011). This method calculates the flash rate (flashes per minute) using the most recent 2 min of flash data. To reduce variability in the 2-min flash rate, the most recent 2 min are split into two 1-min intervals and averaged together.
All clusters with beginning times prior to 1805 UTC (1405 LT) and end times after 0100 UTC (2100 LT) were excluded from further consideration. While it was previously stated that all flashes between 1805 and 0300 UTC were included in the clustering process, the temporal window was narrowed so that clusters ongoing before 1800 or extending past 0100 UTC would be excluded. Artificially beginning all clusters at 1800 or ending all clusters at 0100 UTC might cause truly long-lived clusters to appear misleadingly short lived. Therefore, only clusters with natural begin times after 1805 UTC and end times prior to 0100 UTC would be retained for further consideration. Figure 2 displays the results of the clustering process for a single day at the conclusion of this stage of the analysis.
To describe the spatial characteristics of each cluster, the coordinates of all flashes representing each cluster were imported into a geographical information system (GIS). The GIS framework presented an excellent platform for analyzing the spatial properties of the clusters and proved to be a useful analytical tool in a recent total lightning study (Rudlosky and Fuelberg 2013). The minimum bounding geometry (MBG) tool within the Environmental Systems Research Institute’s (ESRI) ArcMap, version 10.0, GIS software was used on each cluster to draw the convex hull bounding polygon of 1) all flashes that occurred during the cluster’s 5 min of peak lighting activity (Fig. 3), from here on referred to as the maximum flash rate polygon, and 2) all flashes occurring across the entire lifetime of the cluster, from here on referred to as the total lifetime polygon. The MBG tool was also used to draw the circumscribing circle of each maximum flash rate polygon (Fig. 3), and its diameter was recorded for later use. The area of the maximum flash rate polygon was determined using ArcMap’s calculate geometry feature and recorded alongside the attributes described above. Since the calculation of a cluster’s areal extent was foundational to this stage of the analysis, all clusters with fewer than three flashes were discarded. By definition, these clusters did not possess a calculable area and, therefore, could not be analyzed spatially.
The total lifetime polygon was also used to eliminate any clusters that occurred within 2 km of the x–y domain boundary. If only partial lightning data for a cluster existed due to the restriction of the study domain, a large, long-lived cluster might appear small and brief. The 2-km buffer ensured that such clusters were eliminated from further consideration. Unfortunately, in a few instances, this buffering process removed very long–lived clusters that stretched across much of the study area. While eliminating these clusters was undesirable, it was judged to be better to exclude them than allow incomplete clusters to remain in the analysis. After removing all clusters beginning before 1805 UTC (1405 LT), ending after 0100 UTC (2100 LT), possessing fewer than three flashes, or existing within 2 km of the domain boundary, 8445 clusters remained for further analysis.
The final descriptive quantity generated for each cluster was a measure of “compactness” derived from the diameter of the maximum flash rate polygon’s circumscribing circle and the area of the maximum flash rate polygon. The calculation of compactness follows the formula summarized by MacEachren (1985), whereby the diameter of the circle with an area equivalent to the maximum flash rate polygon is divided by the diameter of its circumscribing circle. For a very circular polygon, the diameter of the equivalent circle and the circumscribing circle will be very similar, yielding a compactness value near one. Conversely, very linear shapes will have much larger circumscribing diameters and, thus, compactness scores closer to zero.
The compactness measure is an advantageous descriptor of the clusters because more circular clusters are distinguished from more linear clusters. Since a linear storm structure is often indicative of increased storm organization, and a non-single-cell mode, it is important to describe this spatial attribute of each cluster. It should be noted that storm motion can exercise some influence on compactness; however, the contribution of storm motion over the 5-min period used to calculate the maximum flash rate polygon would be minimal.
c. Quantifying thunderstorm organization
With a full suite of spatial and temporal attributes calculated for lightning clusters that occurred entirely within the bounding domain, the objective then shifted toward assigning each cluster a score representing its duration, areal extent, average velocity, and shape. To aid in this assessment, a storm index (SI) was developed to measure how well the duration, areal extent, shape, and motion of a cluster matched the scientific and operational knowledge of single-cell storm behavior. The goal of the SI was to assign large scores to small, nonlinear, short-lived, and slow-moving clusters while assigning smaller scores to large, asymmetric, long-lived, and fast-moving clusters. Storms adhering to the former expectations would be considered single cells while storms possessing the latter attributes would be considered multi-/supercells.
It has previously been hypothesized that convective mode is a continuous spectrum with no clear boundary between modes (Weisman and Klemp 1982, 1984). Therefore, it is necessary that the SI be a flexible measure that compares a storm’s overall behavior to single-cell expectations, but without using a hard cutoff for any single parameter. Since single-cell thunderstorms are the least organized form of convection, these storms served as the baseline by which more organized modes would be assessed. Therefore, for each parameter used to characterize a cluster/storm, two values were selected: one typical of a single-cell thunderstorm (where any storm meeting this mark could be safely considered single cell in regard to the selected parameter) and one typical of a non-single-cell thunderstorm (where any storm violating this mark could be safely considered non-single-cell in regard to the selected parameter). Any storm meeting the single-cell storm threshold would receive a one for the selected parameter. The score would then linearly decay to zero, with the zero value represented by the non-single-cell storm threshold. However, scores generated by this method were not restricted to between zero and one (Fig. 4). Instead, they were allowed to decay past zero so that storms with attributes grossly violating a parameter could be appropriately qualified. As a preliminary study, this analysis assumed a linear decay function between the two thresholds; however, it is possible that other decay models might be more appropriate or yield stronger results.
Unfortunately, there is very little literature to guide the threshold definitions for the different structural characteristics of ordinary and nonordinary thunderstorms. The primary source utilized is the 1949 Thunderstorm Project sponsored by the U.S. Weather Bureau (BB49), and it was the range of observations in this published study that provided the thresholds for duration, area, and motion used in this work. Though BB49 observed thunderstorms in both Ohio and Florida, this analysis elected to rely upon the Florida observations when available since it was more likely that these storms were embedded within a weakly sheared air mass. The study’s focus on comprehensively documenting nearly every spatial and temporal attribute of short-lived, single-cell thunderstorms qualified it as an excellent resource despite being conducted over 60 years ago.
1) Duration thresholds
According to the observations of BB49, isolated cell thunderstorms produced very short-lived radar echoes, lasting only 20 min on average. To place a conservative lower limit on the total duration of a single-cell storm, 20 min was selected as the lower threshold for the duration score. Any cluster lasting less than 20 min was safely within the limits of the single-cell storm regime. Additionally, the observations of the report placed the maximum combined duration of all three ordinary thunderstorm stages (cumulus, mature, and dissipating) at 75 min (BB49). This duration was selected as the upper threshold for the single-cell regime, meaning any storm lasting longer than 75 min would receive a duration score less than zero. The function relating storm duration to the numerical score assigned to storm duration is depicted in Fig. 4a.
2) Areal extent thresholds
Within the Thunderstorm Project’s discussion of the horizontal extent of thunderstorms, it is shown that an areal extent of 51.8 km2 was the average for storms with cloud tops greater than 25 000 ft (7620 m) (BB49). Since both older and more recent studies (Shackford 1960; Zhang et al. 2013) have found the 20 000–30 000-ft (6096–9144 m) cloud-top range to be the minimum for steady lightning production, 51.8 km2 was selected to represent the lower threshold. While BB49 did not include the individual observations, presenting group averages instead, the figures included in the report show that a cross-sectional area of 64.7 km2 represents roughly the 85th percentile for observed thunderstorms. This value was selected as the upper areal extent threshold, rather than the largest observed storm, since it is more likely storms nearer the maximum areal extent might be aided by atmospheric organization unbefitting a single-cell thunderstorm environment. The function relating the areal extent of a storm to the assigned score based on areal extent is depicted in Fig. 4b.
3) Motion thresholds
The lower threshold on storm motion was set equal to the smallest mean storm motion observed within the Thunderstorm Project (BB49) through the study’s data collection effort in Florida (5.9 kt, 3.0 m s−1). Similarly, the upper threshold (10.3 kt, 5.3 m s−1) was set equal to the largest mean storm motion observed in the data collection in Florida. Mean storm motions less than 5.9 kt were considered safely within the single-cell regime, while storm motions greater than 10.3 kt were deemed less likely to be characteristic of the single-cell mode. The function relating storm motion to the score assigned to storm motion is depicted in Fig. 4c.
It should be noted that supercell thunderstorms have been observed to remain nearly stationary while still maintaining supercell characteristics (Bunkers et al. 2000). In cases such as these, a highly organized storm mode would be mistakenly scored as a single cell with respect to the storm motion threshold.
4) Shape thresholds
While BB49 documented the shape of storms in a handful of cases, the observations were not great enough in quantity or detail to aid in the determination of the shape threshold. Therefore, the thresholds for uniformity in the shape of the lightning patterns were extracted from the total distribution of the compactness measure, the calculation of which was outlined earlier. However, the histogram of compactness values for the 8445 remaining clusters is bimodal with a clear maximum at compactness measures below 0.10 (Fig. 5a). As illustrated in Figs. 5b and 5c, this is a result of the large number of few-flash clusters possessing very low compactness values. Clusters with very few flashes will likely appear irregular and result in a large number of very low compactness scores.
Consequently, the upper and lower quartiles of the distribution are heavily influenced by these numerous, few-flash clusters that likely do not represent informative convective events. To prevent shape scores from being heavily influenced by what are quite possibly inaccuracies in the data, all clusters with fewer than 10 flashes were ignored (Fig. 5d). These clusters were not discarded; they were simply excluded from the distribution while establishing the shape thresholds. The lower threshold for shape was selected as the 75th percentile representing the top quartile of uniform clusters. Alternately, the 25th percentile was selected as the upper limit representing the lowest quartile of circularity measures. The function relating storm compactness to the score used to represent storm shape is depicted in Fig. 4d.
The thresholds selected for each cluster attribute are, of course, subjective, and arguments could be made for other threshold values. As such, several different threshold values were utilized prior to finalizing those defined above, and the results remained rather consistent. Under any selection of thresholds, small, short-lived, slow-moving, circular clusters of lightning flashes receive larger scores and large, long-lived, fast-moving, asymmetrical clusters of flashes receive smaller scores.
5) Relative parameter weights
With scores for each parameter established, the individual scores needed to be combined into a single SI score. Since the SI is intended to be a simple linear combination of the scaled parameters, only the relative weights of the parameters were needed to calculate the SI. With no published work available to guide this task, operational forecasters from the Blacksburg NWSFO were engaged to gain forecaster perspectives. While the forecasters had little experience with interpreting total lightning data, they were asked to compare the relative importance of each storm parameter as they would interpret them using radar imagery. Since the index is only concerned with the relative importance of each parameter, it was not necessary for the forecasters to be familiar with the lightning data source.
The forecaster perspectives were quantified using the analytic hierarchy process (AHP) as described by Wind and Saaty (1980). Each forecaster was asked to make pairwise comparisons between each unique combination of the four chosen parameters, and for each pairing, choose which parameter was more important in discerning a single-cell thunderstorm mode. A rating scale from one to nine was then used to describe the degree to which the more important factor was considered to be superior in the forecaster’s view. A rating of one meant that the “two activities contribute equally to the objective” while a rating of nine meant that “the evidence favoring one activity over the other is the highest possible order of affirmation” (Wind and Saaty 1980). Individual hierarchies were established for each of the five forecasters from the Blacksburg NWSFO before calculating the geometric means of the individual results to produce a group decision (Saaty 1994). Freely available software from Business Performance Management (http://bpmsg.com/) was used to perform the AHP. The geometric mean of the five results yielded the following group result: duration, 20.3%; motion, 18.0%; area, 17.4%; and shape, 44.3%. With a consistency ratio of 0.4%, the forecaster responses were excellent candidates for the AHP.
With relative weights established, the SI was calculated for each storm using the formula shown below. The terms within the expression are written in conditional notation in order to communicate the piecewise nature of the functions shown in Fig. 4. If the given parameter satisfies the condition listed first within the brackets, then it receives a score of one for that parameter score. Otherwise, the final expression within the brackets applies:
6) Filtering potentially misleading clusters
Although clusters of one and two flashes were discarded during the spatial analysis stage, further refinement was required to determine the minimum flash threshold at which a cluster’s SI score represented a meaningful commentary on its parent environment. For the parameters utilized within the SI to have a practical meaning, a minimum cluster flash threshold is needed. For example, it is difficult for a few-flash cluster to be anything but small, short lived, and display little mean movement. Therefore, few-flash clusters will inevitably possess large SI scores regardless of the convective environment in which they formed. These clusters are labeled “potentially misleading” because it is unclear whether they possess large SI scores as a result of storm parameters characteristic of a weakly sheared environment or large scores simply by virtue of being a few-flash cluster. As a preliminary study, this analysis chose to rely upon the spatiotemporal properties of healthy, robust convection (i.e., many-flash clusters) since these can either be long or short lived, large or small, stationary or mobile, and circular or linear. If spatiotemporal patterns of total lightning flashes are truly useful in identifying storm mode, their basic utility will be more readily apparent among many-flash clusters than few-flash clusters.
To establish the minimum flash threshold, the SI was compared to the conditions of the atmosphere in which each cluster formed. Traditional definitions assert that single-cell thunderstorms form in environments with weak vertical wind shear (Markowski and Richardson 2010). Thus, the median SI on a given day should be inversely proportional to the vertical wind shear. Median daily SI versus wind shear was chosen over individual cluster SI versus wind shear so that the relationship would not be skewed by single days with many clusters. The minimum cluster flash total required for a cluster to be included in the analysis was varied between 3 and 300 flashes at one-flash intervals. In each case, the Pearson correlation coefficient R was calculated between median daily SI and 0–6-km wind shear as measured from the 0000 UTC KRNK sounding. The results indicate a very clear step in correlation coefficients near a per cluster minimum flash threshold of 50 (Fig. 6). Thus, a minimum threshold of 50 flashes was selected to differentiate between potentially misleading clusters of lightning flashes and clusters that reasonably reflected the organization of their host environment.
d. Individual thunderstorm and thunderstorm day identification
As stated in the previous section, the intent of this preliminary analysis is to assess the storm mode discernment capability of total lightning flashes using the SI scores clearly established, healthy convection. These clusters were deemed to be the most dependable indicators of the parent convective environment. While a 50-flash minimum seemed to eliminate a potentially misleading subset of few-flash clusters with large SI scores, it could not yet be said whether all of the remaining clusters represented healthy, robust convective events. Thus, it was still necessary to identify the minimum cluster size that would best characterize a healthy, well-established thunderstorm. An additional threshold was then selected to ascertain the minimum cluster flash total that would best characterize the state of the atmosphere in which the clusters formed. Future studies might investigate an enhanced methodology that avoids the institution of a second minimum cluster flash threshold.
Since storm mode is largely dependent on the local wind field (Markowski and Richardson 2010), the state of the atmosphere in which the storms formed was used to calibrate the results of the SI. Additionally, yet to be determined was the number of individual clusters per day necessary to yield a robust description of the convective environment. For instance, it would likely be inaccurate to use the data for the regional atmospheric sounding (KRNK; Fig. 1) to characterize the environment as conducive to single- or multi-/supercell thunderstorms on a day marked by a single cluster of 50 flashes located far from the sounding site. To establish a number of daily discrete clusters necessary to provide robustness in atmospheric analyses, daily cluster minimums were instituted and the correlation between the median daily SI and 0000 UTC wind shear was recalculated while excluding days failing to meet each cluster minimum that was tested. However, instituting the daily minimums also decreased the sample sizes used to generate the correlation coefficients, and, thus, the value of the coefficient became increasingly sensitive to individual data points. To conserve sample size and the reliability of the correlation calculation, five series of correlations were calculated (Fig. 7)—one each for minimum daily number of clusters of 1, 2, 4, 6, and 8.
Between 50 and 100 flashes, two of the five series (six and eight minimum daily clusters), demonstrated a strengthening correlation between median daily SI and 0000 UTC 0–6-km shear. These two series both reached a maximum inverse correlation near 90 flashes with the correlation coefficient stabilizing near 100 flashes, indicating that clusters with at least this number of flashes most accurately and consistently characterized the environmental wind shear. Since only two data series possessed a stable minimum correlation near 100 flashes, it seems that the most accurate characterization of the convective atmosphere could be achieved using either a 100-flash, six-cluster minimum threshold or a 100-flash, eight-cluster minimum threshold. Since the SI–shear correlation when using the six-cluster minimum is more consistent across the entire minimum flash domain, the 100-flash and six-cluster daily minimums were instituted for the remainder of the analysis. This minimum cluster total agrees well with the observations of Goodman et al. (1988), who documented a robust, 116-total-lightning-flash thunderstorm within a weakly sheared environment.
From here on, clusters with at least 100 flashes will be referred to as events, and days with at least six events will be referred to as storm days. While it is possible that the 0000 UTC (2000 LT) soundings were influenced by convection occurring earlier in the day, wind shear measured at that time exhibited a stronger correlation to median daily SI than did wind shear measured at 1200 UTC (0800 LT; Fig. 8); thus, 0000 UTC 0–6-km wind shear was used as the standard to which the SI was compared, and 0000 UTC variables were chosen to characterize a day’s convective environment. Future generations of total-lightning-based classification efforts might seek to incorporate a smaller daily event threshold or eliminate the need for a threshold altogether.
4. Results and discussion
In this section, the events and storm days determined using the previously defined process are compared to the environments in which they formed. Days documenting predominantly lightning-defined (LD) single-cell storms are compared against the traditional expectations of a single-cell storm environment while days with predominantly LD multi-/supercell storms are compared against the traditional expectations of an organized storm environment. In some respects, this section serves as a validation of the above method.
a. Analyzing the SI distribution
Of the 123 days between 1 May and 31 August 2012, lightning was observed within the study area on 97 days. However, three lightning days recorded all flashes outside of the 1805–0100 UTC temporal bounds, and 23 days did not produce a cluster with at least 100 flashes. Of the 71 days with at least one 100-flash cluster, 35 days did not contain at least six such events. The result was a final pool of 36 thunderstorm days. The histogram (Fig. 9) of the median daily SIs was used to subjectively group the 36 days into two classes (Table 1) based on the distribution of the scores. The use of median daily SI tiers allows a comparison between days on which the simple majority of storms were quantified as single cells and days on which the majority of storms were quantified as multi-/supercells. Since storm mode is largely a function of the parent atmosphere, 0000 UTC 0–6-km wind shear should be smaller on predominantly LD single-cell days in order to remain consistent with traditional definitions. Thus, summary statistics for each tier were generated to characterize the atmosphere for the group of days in each tier (Table 2). It should be recalled when interpreting the results in Table 2 that 0–6-km wind shear was used to inform the minimum flash thresholds and minimum daily event thresholds that would constitute an LD storm and LD storm day.
A pooled Student’s t test for unequal variance was conducted (also referred to as a Welch’s analysis of variance procedure) to test the differences in the distributions of 0000 UTC 0–6-km wind shear (Fig. 10) for statistical significance. This test yielded a p value of 0.0102, suggesting that 0000 UTC 0–6-km wind shear on tier 1 days is statistically smaller than on tier 2 days with nearly 99% confidence. Based on these results, tier 1 is consistent with the expectations of a weakly sheared single-cell thunderstorm environment while tier 2 matches the characteristics for a more strongly sheared multi-/supercell environment. Though the mean shear on tier 2 days (26.5 kt, 13.6 m s−1; Table 2) is weaker than values typically associated with supercell thunderstorm environments (Thompson et al. 2003, 2007, 2012), the values satisfactorily correspond to wind shear supportive of multicell storms (Markowski and Richardson 2010). Since it is likely that organized storms occurring during the 1 May–31 August time period in this region were predominantly multi- as opposed to supercell, a mean shear of 26.5 kt is a reasonable characterization for the multi-/supercell category.
It is worth noting that one of the tier 2 days recorded a low 0–6-km wind shear value (<10 kt, 7.7 m s−1). While traditionally this day might have been characterized as a single-cell environment, comparison to radar revealed that these storms were closely tied to the heterogeneous topography of the study area. The highly variable terrain of the study domain frequently leads to preferential storm propagation along differential heating boundaries associated with ridgelines (J. Hudgins, NWSFO Blacksburg, 2013, personal communication), and the resulting lightning clusters consequently appear larger and more linear than anticipated. This topographic influence on storm propagation has been well documented in mountainous terrain (e.g., Hallenbeck 1922; Klitch et al. 1985).
Additionally, four tier 1 days recorded 0–6-km shear greater than 20 kt (10.3 m s−1) with the highest observation being 24.9 kt (12.8 m s−1). While still at the low end of the range of common shear values associated with organized, multi-/supercell storm modes (Markowski and Richardson 2010), these values might be somewhat larger than operationally expected for a single-cell thunderstorm day. However, as their high SI scores indicate, the median event on these days was nonetheless a relatively small, short-lived thunderstorm for which the spatial and temporal properties were indicative of a single-cell mode.
As previously stated, a potential advantage of lightning-based classification is the ability to automate a traditionally labor-intensive task that relies on radar analysis. As such, no effort was made in this preliminary study to individually classify all 470 LD single- and multi-/supercell storms as guided by existing radar-based classification techniques (e.g., Smith et al. 2012). However, each single- and multi-/supercell day was compared against radar to ensure they met the subjective operational expectations of their respective radar structures. As an example, Fig. 11 contains a comparison of two thunderstorms, one LD single cell and one LD multi-/supercell, against Weather Surveillance Radar-1988 Doppler (WSR-88D) imagery from the Blacksburg radar (KFCX) with the maximum flash rate polygons overlaid.
b. Synoptic atmospheric evaluation of SI categories
The stratification of thunderstorm days also allows for an evaluation of how the mean synoptic patterns for the LD storm modes compare to traditional expectations for single- and multi-/supercell convective environments. As with the immediately preceding section, the results of this evaluation seek to lend additional credibility to the lightning-based classification method presented herein. In a sense, it is a further validation of the methodology outlined in the previous sections.
For each category of storm days, daily mean composites of the synoptic-scale atmosphere were generated to illustrate the synoptic meteorological characteristics of days on which predominantly LD single- or predominantly LD multi-/supercell thunderstorms occurred. The publicly available Earth Systems Research Laboratory/Physical Science Division (ESRL/PSD) Daily Mean Composite tool (www.esrl.noaa.gov/psd/data/composites/day/) for depicting reanalysis data (Kalnay et al. 1996) was used to generate the maps. It should be noted that LD single- and multi-/supercell days were both evenly distributed throughout the 1 May–31 August time period, minimizing the influence of seasonal variation in the composite synoptic patterns.
Synoptic maps of 500- and 850-hPa geopotential heights (Fig. 12) illustrate greater values indicative of a warmer lower atmosphere across the eastern United States on single- compared to multi-/supercell days. A westward expansion of the eastern Atlantic Ocean Bermuda high is evident at 850 hPa on LD single- (Fig. 12b) relative to LD multi-/supercell days (Fig. 12c). This circulation feature has long been tied to eastern United States summer precipitation patterns (Stahle and Cleaveland 1992). During LD multi-/supercell days, a trough in the 850-hPa height pattern is evident in the Great Lakes region (Fig. 12c). This is in contrast to the pattern on LD single-cell days (Fig. 12d) and suggests a more dynamic atmosphere that could yield the forcing for sustaining larger and longer-lived storms on LD multi-/supercell days.
While it might be expected that moisture and instability would not differ greatly between single- and multi-/supercell days (Fig. 13), the link between storm mode and atmospheric wind shear leads to an expectation of differences in the strength of atmospheric flow in the mid- and low levels of the atmosphere on the two types of days. This appears to be the case when examining composites of zonal wind at various levels within the lower atmosphere (Fig. 14). Although zonal wind only captures the east–west component of the wind velocity, it was judged that this quantity would highlight any basic differences in the magnitude of the synoptic wind field. If differences in the combined u and υ components of the vector wind were shown, the magnitude of the difference vector would always be positive and it would be unclear which set of days experienced a stronger wind field. Vector wind maps were generated for each set of days to confirm the presence of primarily east–west geostrophic flow, ensuring meaningful representation of the difference in zonal wind between the two sets of days.
The large swath of negative values at the 500- (Fig. 14a), 700- (Fig. 14b), and 850-hPa (Fig. 14c) pressure levels indicates that winds were indeed stronger on LD multi-/supercell days than single-cell days. The stronger wind fields depicted in Fig. 14 would likely contribute to the greater 0–6-km wind shear values observed on LD multi-/supercell days. When combined with the results shown in Figs. 12 and 13, it is apparent that LD multi-/supercell days are potentially characterized by the approach of upper-level dynamical support that is not present during LD single-cell days. The relative similarity of precipitable water and lifted index values characterizing LD single- and multi-/supercell days suggests that LD multi-/supercell storm days occur in the warm, moist air mass ahead of a possible surface front. This would support the greater 500- and 850-hPa geopotential height differences located well northwest of the study area. Observations in these areas would reflect conditions west of the front, typically cooler (Figs. 12a,c), less humid (Fig. 13a), and more stable (Fig. 13b) than ahead of it.
Taken in their entirety, the synoptic composites indicate that humidity and instability are comparable between LD single- and multi-/supercell days across the study area. However, wind strength and 850-hPa geopotential height patterns (Fig. 12c) suggest that LD multi-/supercell days are more dynamically driven than LD single-cell days. These synoptic differences are supported by previous research findings. Fuelberg and Biggar (1994) observed very similar thermodynamic regimes between days with “strong convection” and “weak convection” while thermodynamic variables between days with convection versus without convection yielded statistically significant differences.
This study represents a preliminary attempt to classify thunderstorms into two basic categories of storm mode using the spatiotemporal distributions of total lightning clusters. The results of the work presented here suggest that spatiotemporal total lightning patterns can be used to differentiate between single-cell and more organized storm modes. Days exhibiting a majority of smaller, shorter-lived, slower-moving, and more circular storms were characterized by significantly weaker wind shear than days with a majority of larger, longer-lived, faster-moving, and more linear storms according to the statistical assessments. A series of synoptic atmospheric composites also suggested that LD single-cell storms formed in weaker-flow regimes than LD multi-/supercell storms. As total lightning data continue to gain prominence within the operational forecasting toolset, the integration of lightning data into traditional storm mode classification research methodologies could prove a fruitful addition.
Though radar-based classification possesses a much stronger pedigree than lightning-based classification, the potential benefit of this new classification tool is very attractive. LD storms can be easily analyzed for total lightning jumps (Schultz et al. 2009, 2011), a recent focal point of severe weather research. Additionally, a standard definition of LD storm modes could allow for the frequency and intensity of thunderstorms to be relatively easily tracked through time, helping clarify potential implications of climate change. The validations conducted as part of this analysis lend credibility to the potential utility of lightning-based classification.
Total lightning-based identification was limited in its ability to discern longer-lived terrain-driven events from convection supported by atmospheric dynamics. This appears to be a potential shortcoming of any purely lightning-based storm mode classification system. In addition, many aspects of the analysis are open to refinement and improvement by future researchers: most notably, the choice of clustering method, the time–space standardization factor, the single-linkage termination procedure, thresholds used to inform the SI, the relative weight of the SI components, the minimum cluster flash thresholds, and the minimum daily event thresholds. While the SI utilized in this study seems to have provided a satisfactory description of an event, it is very likely that the SI would need to be tailored for use outside of the central Appalachian region of the United States. The index was created by using responses from forecasters who monitor a relatively small region with a somewhat location-specific storm climatology. It is certainly possible that forecasters in different NWSFOs could very well diagnose storm mode using a different set of weights than those utilized in this study, or that superior results could be obtained using an equal weighting of the index’s inputs. The SI is also limited by the relatively short temporal span used to create the histogram that led to the delineation of tiers. A longer period of observation would likely provide a smoother, fuller histogram with different tier cutoffs. Further research could focus on expanding both the spatial and temporal areas of consideration.
This material was based upon work supported by the COMET Program of the University Corporation for Atmospheric Research (UCAR) and the National Oceanic and Atmospheric Administration’s (NOAA) National Weather Service (NWS) under Grant Z13-99434. Any opinions, findings, conclusions, or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the COMET Program, UCAR, NOAA, or the NWS. The authors thank M. Marston of the Department of Geography, Virginia Tech, for help with data analysis, and the staff at Earth Networks, Inc., for supplying the lightning data used in this project. The authors also express their gratitude to the anonymous reviewers, who helped substantially increase the strength and clarity of this text from its original form.