Data-driven quantification and visualization of resilience metrics of power distribution systems

Outage and weather data description

This study utilizes two datasets: (i) outage records provided by a Midwestern US distribution utility, and (ii) open weather data and severe weather events from NOAA, including wind and precipitation measurements from multiple weather stations and all severe weather-related extreme event history within the utility’s service area28,29. Both datasets span a 20-year period from 2004 to 2024 and are temporally aligned to support joint analysis.

The outage dataset contains a total of 168,462 individual outage records. Each record corresponds to a discrete outage that affects a component in the distribution system and includes the start and end times of the outage, accurate to seconds; the duration of restoration in minutes; the number of customers interrupted; and the associated cause code. The cause codes in the outage dataset are detailed and associated with each outage. Examples include tree-related causes (e.g., “Tree limb on line,” “Tree outside trim zone”), animal interference (e.g., “Squirrel contact,” “Bird contact”), equipment failures (e.g., “Fuse blown,” “Switch failure,” “Underground cable fault”), and weather-related causes (e.g., “High wind,” “Ice/snow accumulation,” “Lightning strike”). Although the proposed models do not rely directly on cause codes, their specificity helps validate the event categorization and interpret outage patterns.

Before proceeding with analysis, the dataset was carefully preprocessed to improve reliability and internal consistency. Due to equipment malfunctions or manual reporting errors, historical outage data often contains incomplete or inaccurate entries, which can lead to misleading conclusions in data-driven analysis. To address this, records with missing values in key fields, such as timestamps, restoration duration, or customer count, were discarded. In addition, entries with logically inconsistent values were filtered out. These included cases where the restoration time exceeded the overall outage duration or where the values fell outside reasonable operational limits.

The weather sensor data includes hourly wind and precipitation data from 8 stations across the study region. Wind-related features include hourly average wind speed and fastest 2-minute wind speed, all recorded in meters per second. Precipitation-related features include hourly totals of rainfall, snowfall, and ground-level snow depth, measured in inches. All sensor data were preprocessed to ensure temporal alignment with outage records. The severe weather event dataset29 is obtained from a separate NOAA database and consists of documented weather incidents classified by event type. Each record includes a start and end time, accurate to the minute, event type, starting location, and affected area of the event represented as longitude and latitude, and a brief description of the event. The types of events considered in this study include tornadoes, high-speed wind events (classified as wind events), snowstorms, and flooding (classified as precipitation events). However, unlike the weather sensor dataset, severe event records do not contain direct numerical measurements of wind speed or precipitation.

GIS mapping of component outages

To visualize the spatial distribution of the outages in the utility service area, each outage record from the OMS was mapped to its corresponding geographical location. The OMS database provides information on damaged components or protective devices, such as tripped breakers or blown fuses that triggered each outage. Using the unique identifiers for these components, their GIS coordinates were retrieved and plotted to generate a heatmap of the outage density. Figure 3 presents the resulting heatmap, where the regions with higher outage concentrations are shown in red and yellow. This geospatial representation demonstrates clusters of vulnerability within the service territory. However, it is important to note that OMS records do not contain any weather measurements or meteorological context, the weather dataset mentioned in the previous subsection needs to be aligned to each OMS record.

Outage-restoration events extraction

To systematically analyze outage accumulation and restoration behavior, this work defines outage-restoration events based on the temporal progression of component outage records. Rather than independently evaluating each outage, outages are grouped into events that represent continuous disturbance intervals followed by full system recovery. As discussed in the Introduction section, restoration performance is not solely determined by the characteristics of a single outage, but is heavily influenced by the number and timing of coincident outages within a given time frame. When multiple outages occur in rapid succession, often triggered by severe weather, the restoration process becomes increasingly complex due to limited repair resources, unavoidable delays, and difficulty in arrangement. By aggregating temporally overlapping outages into outage events, this process will analyze how restoration time scales with the number of concurrent outages and provide a more realistic representation of system behavior under high system stress conditions.

An event begins with the first outage that occurs when all components are functional and ends when all outages have been restored, bringing the cumulative number of active outages, denoted as C(t), back to zero. Specifically, the process involves sorting all outage start and restoration times in chronological order and tracking C(t) at each time point. The value of C(t) increases when a new outage appears and decreases when an outage is restored. Each time C(t) returns to zero, the current event concludes, and any subsequent outage marks the start of a new event. For example, if an event includes n outage onsets, it must also contain n corresponding restorations to ensure that C(t) returns to zero at the end of the event window. This approach captures the cumulative behavior of coinciding outages and allows each event to be characterized by two primary features: the number of outages and the total restoration time, defined as the event restoration time from the first outage to the final restoration within the event. Figure 4 illustrates representative examples of extracted outage-restoration events. These plots illustrate how outages accumulate and are subsequently resolved during each event. The curves typically exhibit a sharp increase in the outage count at the beginning of an event, followed by a gradual decline as restorations are completed. The variation in shapes and durations reflects the diversity of event profiles, including both short, concentrated disruptions and longer, more complex restoration processes.

Note that in this study, an outage-restoration event is defined to end when the cumulative number of active outages returns to zero. This represents the strictest possible boundary condition and ensures that the entire period of disturbance and recovery is captured. We acknowledge that this definition may extend event durations and occasionally include background outages not directly related to the severe weather event, since utilities often experience a small number of routine outages even under normal conditions. In practice, utilities may prefer a more flexible threshold, such as defining event completion when 98% of outages have been restored or when the number of remaining outages falls below a small tolerance level. The proposed framework is compatible with such customizations, and the threshold can be adjusted according to utility preferences and operational standards. The strict criterion adopted here is used to demonstrate the methodology under the most conservative setting.

Delineation of weather zone via voronoi polygon method

To capture the spatial variation in weather conditions across the service territory, the studied area is delineated into localized weather zones based on the distribution of NOAA weather stations. Specifically, two wind and six precipitation stations located within the region are used to distinguish wind and precipitation weather zones for resilience analyses. Each weather station serves as a reference point for one zone, ensuring that the weather data associated with each outage are drawn from the most geographically relevant source.

The division is accomplished using the Voronoi polygon method, a spatial partitioning technique that assigns every point in the region to its nearest weather station. Given a set of weather stations \({s_1, s_2, \dots , s_n}\) with known geographic coordinates, the Voronoi cell for station \(s_i\) consists of all points in the plane that are closer to \(s_i\) than to any other station \(s_j\), where \(j \ne i\). This results in a tessellation of the area into non-overlapping polygons, each representing the area of influence of a particular weather station. The method is implemented using the QGIS software platform, which supports geospatial processing and visualization of shape files based on station coordinates.

Based on the spatial distribution of NOAA weather sensors and the measurement types available at each station, two wind zones and six precipitation zones were defined for the utility service region using the Voronoi polygon method. Wind zones are determined using the locations of stations that record wind-related features, such as hourly wind averages and gust speeds. The precipitation zones are defined using the locations of sensors that provide hourly rainfall, snowfall, and snow depth. The defined wind and precipitation zones are shown in Fig. 5. The left part shows the division of the service territory into two wind zones, each surrounding one of the two wind stations. The right part shows the division into six precipitation zones, each linked to one of the six precipitation stations.

Data-driven resilience models

This section introduces two models to quantify resilience across the distribution utility’s service territory. The first is a fragility model that correlates the number of outages with weather variables. The second is a restoration time model that relates the restoration time to the number of outages. These resilience metrics are subsequently quantified using real-world data from a utility.

Data-driven fragility models using severe weather outage events

A structured, multi-step process was followed to analyze system resilience against severe weather events. Each step was carefully designed to ensure the integration of weather records, outage data, and spatially defined weather zones into meaningful fragility models.

The first step involves compiling a dataset of historical severe weather records affecting the service area between 2004 and 2024. This dataset includes records classified by location, duration, and weather types such as high winds, tornadoes, floods, heavy snow, and thunderstorms. In the second step, all severe weather records are classified into two categories: wind-related events (e.g., tornadoes, high winds) and precipitation-related events (e.g., floods, snowstorms). This step classifies and associates the relevant weather variables, such as wind speed or precipitation levels, with each type of event. After that, each weather event is assigned to a defined weather zone based on the spatial coverage of the weather stations and the Voronoi polygon boundaries established earlier. Each event is matched to the zone in which it occurred, which supports the development of fragility curves reflecting the unique conditions of each zone. Then, the number of outages during each weather event is computed. This is done by filtering the outage dataset by event time window and geographic location, capturing only the outages that occurred within the event duration and the corresponding weather zone. In addition, we link each event to its corresponding weather measurements. For wind-related events, wind speed observations using the fastest 2-minute averages are extracted. For precipitation-related events, the hourly precipitation, snowfall, and snow depth are obtained. These values are aligned with the event temporal window to represent the intensity of the event. Finally, the collected data are aggregated at the zone level, yielding paired observations of weather intensity and outage impacts. These data are then used to construct fragility models, which express the number of outages as a function of weather severity. Separate fragility models are developed for wind and precipitation zones using wind and precipitation events, as their purpose is to characterize the relationship between weather intensity (e.g., wind speed or precipitation depth) and the number of outages resulting from those weather conditions. Exponential functions are used for model fitting, as they capture the nonlinear increase in outage likelihood with rising weather intensity.

Data-driven restoration time models using outage events

To characterize restoration performance across varying weather conditions, restoration time models are developed using the full set of outage-restoration events. These models aim to capture the general relationship between the number of outages and the total restoration time, regardless of the underlying cause, and therefore reflect the operational response characteristics across all types of outage conditions.

For each extracted event, as discussed in section “Outage-restoration Events Extraction”, two quantities are summarized: the total number of outages and the total restoration time. These values are plotted in scatter diagrams to explore the empirical relationship between the number of outages and the duration of restoration. Following this, events are grouped by weather zone, and separate restoration models are fitted for each zone. A combination of constant and double exponential functions is selected for its flexibility in capturing both the initial increase and the saturation effects observed in large-scale events.

This method produces a model for estimating restoration times across wind and precipitation zones using the number of outages predicted by the fragility models. By linking outage volume to restoration time in an event-specific manner, this method supports data-driven assessment of restoration performance under different severe weather conditions.

Fragility and restoration time models for wind events

To account for the impacts of wind events on the distribution system, zone-specific fragility and restoration time models are developed. In this study, fragility functions are defined at the zone level rather than for individual events. Each severe weather event that occurs within a zone contributes one aggregated observation, which represents the event’s weather intensity and the total number of outages it caused, to the zone’s fragility dataset. The resulting fragility curve for each zone is therefore fitted using multiple severe weather events accumulated over the 20-year record, reflecting the collective relationship between hazard intensity and outage volume rather than the behavior of any isolated event. Figure 6 displays the results for Wind Zone 0 and Wind Zone 1. The central map shows the spatial division of the two wind zones, with background shading representing the outage density throughout the study area. Wind Zone 0 covers the northern part of the region, while Wind Zone 1 covers the southern part.

The fragility model for Wind Zone 0 is shown in the upper left. The blue dots represent observed outages for individual wind events, plotted against maximum wind velocity (in m/s). The red curve shows the best-fit exponential function of the form \(y = 0.0002e^{0.3675x}\), where x is the wind speed and y is the number of outages. The curve demonstrates a rapid increase in outage counts once wind speeds exceed approximately 28 m/s, which is a threshold effect where system vulnerability intensifies significantly. This pattern reflects the relationship between outage numbers and wind speed in this zone. The corresponding restoration model for Wind Zone 0 is shown in the bottom left. Similarly, each blue dot represents a historical outage event, with the x-axis indicating the number of outages and the y-axis indicating the total restoration time in hours. The red curve represents the fitted model, a two-term exponential function that captures both the initial rapid increase in restoration time and the eventual saturation effect. This relationship shows that, as the outage volume increases, restoration becomes delayed, likely due to constraints in available repair crews and resource allocation.

The upper right panel presents the fragility model for Wind Zone 1. Compared to Zone 0, this zone exhibits a slightly different fragility profile. The fitted exponential model can be expressed as \(y = 2.9214e^{0.1058x}\), and this curve shows a more gradual increase in the outage count with wind speed. This indicates that Zone 1 may have relatively lower resilience against wind-related failures, or possibly differences in terrain or vegetation contributing to this trend. The figure on the bottom right shows the restoration model for Wind Zone 1. As in Zone 0, the restoration time increases nonlinearly with the outage volume, but the fitted curve here suggests a sharper increase at lower outage counts. This pattern may indicate less operational redundancy or more constrained access conditions in Zone 1. The model is defined by the expression \(y = 232.60 – 217.17e^{-0.001x} – 15.22e^{-0.041x}\).

The number of data points in each fragility curve is limited because severe weather events that satisfy the temporal and spatial matching criteria are relatively infrequent. Each such event yields only one valid observation per zone, as the outage dataset provides event-level total outages rather than continuous, asset-level failure counts. Only events with a clear overlap with NOAA severe weather records and the corresponding weather zone are included. This inherent sparsity of severe-weather observations is a well-recognized characteristic of empirical outage datasets and naturally results in a small number of points for fragility fitting.

Note that the blue data points in the fragility and restoration models represent different underlying datasets. For the fragility models, each point corresponds to a single severe-weather event within a specific weather zone and captures the relationship between weather intensity and the total number of outages caused by that event. Because severe-weather events are relatively infrequent and each event contributes only one aggregated observation, the resulting fragility curves contain a limited number of data points. In contrast, the blue points in the restoration-time models represent the full set of outages-restoration events extracted from the 20-year outage history, regardless of cause. These events capture the operational restoration behavior of the utility, including both weather-related and routine outages. Since restoration performance depends primarily on outage volume, rather than being exclusive to severe-weather events, the larger dataset enables the restoration model to generalize across diverse operating conditions.

Fragility and restoration time models for precipitation events

The utility sevice area has six precipitation zones, and Figure 7 presents the fragility and restoration models for each zone. The map in the lower right illustrates the boundaries of the six precipitation zones, and the background is the density of the outage records. In the fragility plots, blue dots represent observed outage counts during individual precipitation events, plotted against precipitation depth (in inches). The red curves correspond to the fitted exponential fragility models. The corresponding restoration time models are shown below or next to the fragility plots. Each blue point represents an outage event, with the x-axis showing the number of outages and the y-axis indicating the total restoration time in hours. The red curves represent fitted double-exponential functions. Although exact coefficients vary by zone, the general trend in all cases shows that restoration time increases nonlinearly with outage volume.

Although the fitted fragility curves for Zones 0, 1, 2, 3, and 5 share the same exponential form \(y = 0.0654e^{2.7683x}\), this does not imply that the outage levels in these zones are identical. The similarity in the curve shapes arises because many precipitation events affecting the area are spatially large and influence multiple precipitation zones simultaneously, resulting in similar precipitation-depth measurements for each affected zone. When aggregated at the event level, the paired observations of precipitation intensity and total outages fall along comparable ranges across these zones, causing the fitted exponential functions to converge. Nevertheless, the actual number of outages observed in each zone varies from event to event, as indicated by the blue data points in the restoration plots and the outage-density heat map. These differences reflect historical variation in exposure, infrastructure layout, and local vulnerability, even when the precipitation intensity is similar across zones.

Zone 4 shows a distinct fragility curve, modeled as \(y = 1.179e^{1.6159x}\), which grows more slowly than the fragility functions in the other zones. This may indicate lower sensitivity to precipitation depth, possibly due to reduced exposure, fewer vulnerable components, or differences in system configuration. The restoration model for Zone 4 also shows a steeper increase in duration with increasing outage count, showing longer recovery periods under stress conditions. For Zone 5, while the fragility model is identical to that of Zones 0–3, the restoration curve deviates. The fitted model indicates a slower increase in restoration time at lower outage volumes, followed by a steeper climb once the number of outages exceeds approximately 100. This behavior may reflect the threshold of repair resources.

It is important to mention that precipitation-driven outages can arise from different mechanisms associated with rain, snow, or mixed precipitation. However, the number of snow-only or rain-only severe events within each zone is insufficient to support independent fragility models for each precipitation type. In this work, all precipitation-driven events are therefore aggregated to construct a single precipitation fragility curve per zone. While this may introduce additional variability in the fitted functions, the aggregation is necessary to maintain statistical robustness given the available data. Future work could incorporate separate rain and snow fragility models when higher-resolution or longer-term datasets are available.

Model performance

To evaluate the accuracy of the proposed framework, we conduct a validation using historical severe-weather events. For each event, the recorded weather intensity is used as input to the fragility model to estimate the number of outages, and the predicted outage count is then used in the restoration-time model to estimate the total event duration. These predicted values are compared with their historical counterparts, and accuracy is quantified using goodness-of-fit metrics, including the coefficient of determination (\(R^2\)) for both the fragility and restoration models. The reported \(R^2\) values in Table 1 show that the fitted models capture the major trends in the outage data, providing a consistent and data-driven basis for resilience analysis.

Visualization of resilience metrics

The developed fragility and restoration time models can be combined to generate zone-specific predictions of outage restoration time under a given weather scenario. For each target weather condition, such as wind velocity (m/s) for a wind event or precipitation depth (inches) for a precipitation event, the fragility model is used first to estimate the expected number of outages in each zone. Specifically, the wind or precipitation level is input into each zone’s corresponding exponential fragility model, yielding an estimated outage count for that scenario. The predicted number of outages is then input into the corresponding restoration time model. The restoration model maps the outage count to an expected total restoration time. This two-step process estimates the restoration time in each zone, given a particular level of wind or precipitation. Applying this procedure across all zones generates a complete set of restoration time estimates for the entire service area. The results are visualized on a spatial map, with each zone shaded according to its predicted restoration time under the specified weather conditions. These visualizations directly compare resilience between zones and help identify areas more vulnerable to specific weather conditions.

Figure 8a presents an example of a wind event with a velocity of 35 m/s. In this visualization, darker shading indicates shorter predicted restoration durations and thus higher resilience. As shown in the figure, Wind Zone 0 exhibits a shorter expected restoration time compared to Zone 1 under the same wind intensity, indicating a stronger resilience performance. Figure 8b shows the restoration performance for a precipitation event with 2.5 inches of rainfall. As the map illustrates, Zones 0 and 2 demonstrate relatively shorter restoration times, indicating higher resilience to precipitation-driven outages, while Zones 4 and 5 exhibit longer restoration durations and greater vulnerability. This spatial comparison can help utilities target specific zones for resilience improvement efforts.