Evaluation of the Spatio-temporal Variation of Extreme Cold Events in Southeastern Europe Using an Intensity-duration Model and Excess Cold Factor Severity Index

1. Introduction

Over millions of years, the climate system behavior has been dominated by natural fluctuations at different time scales, but during the last 150 years, anthropogenic forcing substantially influenced natural variability [1,2]. In the context of climate change scenarios, an increase in global temperature of 2 °C compared to the pre-industrial period is considered an upper limit above which the risk of dangerous and irreversible climate shifts increases sharply [3,4,5]. Given that the rise in global temperature has already exceeded half of this limit, it can be crossed in the next 20-30 years [1,6,7]. Heat and cold waves are among the most impactful large-scale natural phenomena strongly related to climate change. Prolonged periods of extreme temperatures cause thermal stress and may have adverse or fatal effects on human health [8,9,10]. They directly and indirectly influence labor productivity and well-being [11,12]. Extreme temperatures considerably increase the risk of energy supply outages and damage to the infrastructure [13,14]. They can even represent a limiting factor for agriculture in the changing climate due to the amplified likelihood of frosts and prolonged hot spells during the vegetation period [15,16,17].

Against the background of the expected rise in global temperatures, a decrease in cold extremes seems intuitively sound [1,18,19]. The impact of heat waves will significantly exceed those of cold waves till the end of the century. Across all studies, the authors point out the multifold increase in the duration, frequency and intensity of future extreme heat events as highly likely [1,20,21]. In Europe, since the beginning of the century, heat waves have caused more deaths than all other natural disasters [22]. Schär et al. [23] revealed that the record summer of 2003 cannot be explained by solely shifting the mean of the temperature distribution to higher temperatures but with an additional increase in variability. The changed shape of the temperature distribution contributes to the occurrence of more hot extremes but also suggests that future cold extremes cannot be inferred from changes in mean and variance alone [24,25]. In this sense, recent severe cold events in Europe, such as in 2010, 2012 and 2017, have received particular attention [26,27,28,29]. Cold extremes would still be possible in the current climate [30,31] and will continue to impact various socioeconomic sectors, particularly healthcare, due to the complex relationship between human thermal adaptation capability and warming rates. Recent studies indicate increased tolerance to hot weather and significant differences among age groups in response to cold events [10,32,33,34,35].

After the deathly European heatwave of 2003, heat extremes have become a priority topic in climate change research [36]. Despite remarkable progress in recent years and the apparent upward trend in heat wave intensity, duration, frequency and spatial extent, there is still no universally agreed definition for this phenomenon [37,38,39]. Therefore, criteria and procedures for identifying extreme heat at a national level are very diverse [40,41]. This lack of consensus also extends to cold wave definitions and warnings.

An extreme event occurs when a climate variable is above or below a defined threshold at the upper or lower end of the range of its observed values. This threshold can be calculated as absolute (fixed) or relative (percentile-based) value [1]. Since around 60% of the national weather services that have criteria for heat/cold waves use absolute thresholds when considering heat and cold events, the World Meteorological Organization (WMO) Task Team for the Definition of Extreme Weather and Climate Events (TT-DEWCE) specified the term heatwave as unusually hot weather relative to the local climate lasting at least two consecutive days in the warmest months of the year. Similarly, cold waves denote unusually cold weather during the cold season, lasting at least two days and often preceded by a rapid and marked drop in air temperature caused by the advection of cold air [37]. These definitions allow flexibility in reporting heat/cold events when different criteria and approaches are used (from a simple descriptive approach to classification schemes based on local climatology, human perception of heat/cold stress, or increased mortality risk [41]).

Extreme heat events can be quantitatively assessed by indices describing their particular characteristics [42,43,44], and similar indices can be defined analogously for extreme cold events (ECEs). Some basic requirements should be satisfied when developing indices of heat and cold extremes [45]. Generally, they should be suitable for various economic and social sectors or climates, but they should also allow comparative analysis. Indices should be easy to understand, calculate, represent and forecast. They should give the public and decision-makers direct, meaningful information on the expected severity of the heat and cold waves. The advantage of percentile-based indices is their global applicability, but fixed-threshold indices are indispensable when monitoring high or low temperatures, which negatively impact health or the economy, is important [39]. Thresholds can be combined to define cascading schemes accounting for the spatiotemporal evolution of extreme events [46]. The rarity of heat and cold extremes should be determined relative to a standard climate period, so a balance must also be found between the threshold levels and the number of events identified as extreme [37,43].

The most popular and widely used set of climate indices consists of 27 core climate change indices (17 temperature-based) developed by the WMO Expert Team on Climate Change Detection and Indices (ETCCDI), initially designed for scientific purposes (http://etccdi.pacificclimate.org/list_27_indices.shtml) (accessed on 28 January 2025). Subsequently, the WMO Expert Team on Sector-specific Climate Indices (ET-SCI) expanded this list to include custom-defined threshold indices and new heat/cold wave indices based on definitions accounting for excessive heat/cold accumulation and short-term thermal acclimatization (https://climpact-sci.org/indices) (accessed on 28 January 2025). In recent years, there has been growing attention to indices based on health-related definitions, particularly for early warning purposes in large cities [41,47,48].

The region around the Mediterranean basin, including Southeastern Europe, is among the most vulnerable to climate change, experiencing more frequent and intense severe weather events and heat extremes [1,49,50,51]. Besides heat waves, cold events have also been the subject of various regional studies in recent decades. Tringa et al. [52] conducted an extensive climatological and synoptic analysis of winter cold spells in the Balkan Peninsula from 1961 to 2019. They found a decreasing trend in the frequency of cold spells toward the end of that period. Planchon et al. [27] performed a multi-scale agroclimatic analysis for Romanian vineyards in the meteorological context of the early 2012 European cold wave. Demirtaş [28] analyzed the anomalously cold January 2017 in Southeastern Europe through the connection between cold spells and atmospheric blocking using an objective coldwave indicator and a two-dimensional blocking indicator. Kostopoulou [53] investigated the potential effects of the January 2017 cold spell on mortality attributed to cardiovascular disease. Faranda et al. [54] conducted an attribution study on the winter storm Filomena in 2021, marked by historic snowfalls in the inland regions of the Iberian Peninsula and a prolonged 14-day cold spell. Kysely et al. [55] demonstrated that cold stress significantly impacts mortality in Central Europe, representing a public health threat comparable to that posed by heat waves. Díaz-Poso et al. [56] analyzed climate change scenarios using simulations from the EURO-CORDEX project, applying the Excess Cold Factor (ECF) index for the Iberian Peninsula and the Balearic Islands. Their results indicated that despite the expected decrease in cold waves, they would continue to pose a serious local threat due to the population’s acclimatization to higher temperatures. Piticar et al. [57] found substantial changes in the indices of excess heat/cold factor and concluded that heat waves have become more frequent, longer, and more intense, while cold waves have become less frequent but more intense in Romania.

The present study can be viewed as continuing previous research on spatiotemporal variations of hot spells in Southeastern Europe [58]. The indicators proposed here for quantitative assessment of cold-weather phenomena could be helpful in sector-specific climate assessments. The study aims to test the suitability of developed indicators for an overall analysis of ECEs in Southeastern Europe. The ability of the cold spell duration indicator to detect extreme cold events has also been examined and compared with the ECF-based indices.

2. Data and Methods

We consider the same domain over Southeastern Europe (SEE), with coordinates 35°N–50°N and 15°E–30°E (the red frame in Figure 1), as in the previous study [58]. This area almost entirely falls within the Mediterranean region, as defined by the Mediterranean Experts on Climate and Environmental Change (MedECC) [59]. According to the updated Köppen-Geiger classification [60], the European part of the Mediterranean region mostly belongs to the midlatitude temperate climate zone, comprising Mediterranean climate types with a dry hot/warm summer (Csa/Csb) and large areas with no dry hot/warm summer (Cfa/Cfb), but also areas with continental hot/warm/cold climates (Dfa/Dfb/Dfc) in the north [61].

The same dataset, developed for evaluating the spatiotemporal variations of extreme heat events in the period 1961-2020 in SEE [58], is used here to investigate extreme cold events. This dataset consists of daily maximum and minimum air temperature data from 70 selected meteorological stations covering the study area relatively evenly (Figure 1). Data from 12 stations of the national meteorological network of Bulgaria are provided by the National Institute of Meteorology and Hydrology (NIMH). Data from stations outside Bulgaria were downloaded from three freely available sources: the European Climate Assessment and Dataset (ECA&D) project [62], the U.S. National Climatic Data Center (NCDC) Global Historical Climatology Network daily dataset (GHCNd) [63], and the U.S. National Centers for Environmental Information (NCEI) Global Surface Summary of the Day (GSOD) dataset [64]. Most stations have altitudes below 500 m and are distributed across four Köppen-Geiger climate zones [60] – Cfa (25%), Cfb (34%), Csa (21%) and Dfb (20%). The stations’ metadata are presented in Appendix A, Table A1.

We have demonstrated in the previous study the suitability of climatologically justified for Bulgaria thresholds of daily maximum temperature (32, 34, 36, 38 and 40 °C) and corresponding duration thresholds (6, 5, 4, 3 and 2 consecutive days) to classify extreme heat events in SEE. Given that the prevailing in Bulgaria Köppen’s climate subtypes (Cfa and Cfb) is almost 60% of the whole domain, we expect to obtain a sufficiently accurate picture of changes in cold extremes in the period 1961-2020 by combining thresholds for daily minimum temperature and cold spell duration, using the same dataset and methodology as in the estimation of extreme hot spells.

The threshold indicators defined in [58] were developed using statistical modeling of daily maximum air temperature data from 36 stations, including the 12 stations shown in Figure 1, representative of the non-mountainous regions in Bulgaria with at least 80 years of observations in the period 1931-2020, referred to as BG36 hereafter (stations’ metadata are presented in Appendix A, Table A2). In the present study, we strictly followed this approach and used the same meteorological stations to ensure the methodologically consistent analysis of minimum temperatures and extreme cold events. Data on daily minimum air temperature from the BG36 dataset were processed to obtain the distribution of minimum temperatures and corresponding quantiles characteristic of the climate in the country.

Data processing involved several steps: 1) construction of empirical distributions of annual minimum temperatures for 1931–1980 (the period is considered as less affected by global warming [2,65,66]); 2) calculation of quantiles corresponding to return periods 2, 5, 10, 20 and 50 years by fitting the Generalized Extreme Value (GEV) distribution to data, accounting for nonstationarity of time series [67]; 3) determination of appropriate thresholds for daily minimum temperature (tn); 4) determination of cold spells duration thresholds for the period 1961–2020 at the different tn-thresholds.

Until 1980, daily minimum temperatures in the coldest months in Bulgaria (January and February) fall mainly into the range from -27.4 °C to 11 °C across the selected stations, and the absolute minima are typically between -30.8 °C and -22.5 °C. From the estimated 2-, 5-, 10-, 20- and 50-year return levels by the GEV model (Figure 3, middle panel), the set of thresholds (-10, -12, -14, -16, and -18 °C) have been chosen in such a way to be achievable in more than 80% of stations. During the second period (1981–2020), a substantial shift of tn-distributions toward higher values is observed (Figure 3, left panel). Therefore, the cold spell duration thresholds were determined for the period 1961–2020 to provide a reliable analysis of ECEs in the contemporary climate. Cold spells last 2-3 days most frequently, and the 90th percentile varies between 3 and 5 days.

Assuming that the research focuses on applying threshold indicators to analyze extreme cold in a climatically diverse region, such as SEE, we explored the minimum temperature features by Köppen-Geiger climate zones. First, we calculated empirical quantiles, corresponding to the return periods used here, from the time series of annual minimum temperatures (1961–2020) for all 70 stations. As seen in the left panel of Figure 2, the median values remain well below the chosen tn-thresholds (blue dots in the chart). This suggests that more frequent and longer cold spells will be detected in at least 50% of the stations, determining the proposed cold spell duration indicator as an indicator of moderate extremes [42]. The chosen thresholds are well-suited for more than 70% of stations in the domain but not for stations in the Csa climate zone. The thorough analysis showed that for 2-year return periods, -5 °C is achievable in 35% of the Csa stations and 85% of all stations, while -20 °C is achievable in around 10% of stations. The right panel of Figure 2 shows the calculated percentiles of the lower tail of distributions of daily minimum temperatures by stations in the four climate zones as box-plot diagrams. One can see that only percentiles below the 10th, usually used in ECE definitions, fall within the range from -20 °C to -5 °C.

2.1. Climatologically Justified Threshold Indicators

Given the obtained estimates about low temperatures in the studied domain, the cold spell duration indicator was developed based on samples from the BG36 dataset and the extended set of tn-thresholds (-5, -10, -12, -14, -16, -18 and -20 °C). The right panel in Figure 3 shows the distribution of cold spell durations at different tn-thresholds for the period 1961-2020. Duration thresholds were selected as the 90th percentile of the respective distributions: 7 days for tn ≤ -5 °C, 5 days for tn ≤ -10 °C, 4 days for thresholds of -12 °C and -14 °C, and 3 days for thresholds of -16 °C and -18 °C. For cold spells with tn ≤ -20 °C, we selected a minimum duration of 2 days [37]. A recent study on cold waves in the non-mountainous regions of Bulgaria in the period 1952-2011 found that cold events, assessed using the ETCCDI Cold Spell Duration Index (CSDI), mostly last between 6 and 11 days, reaching minimum temperatures between -5 °C and -20 °C [68]. These results confirm that the proposed range of tn-thresholds and the 7-day duration threshold for lowest-intensity cold spells are relevant to capturing the evolution of cold spells in the current climate.

Analogously to the hot weather indicators in [58], three climate indicators of cold weather have been defined:

• Annual number of cold days (ncd-5) – i.e., the annual count of days when tn < -5 °C;
• Maximum number of consecutive cold days (ccd-5) – i.e., the longest continuous calendar period when tn < -5 °C;
• Cold spells duration at different tn-thresholds (csd-5/-10/-12/-14/-16/-18/-20) – i.e., the annual count of days when tn ≤ -5, -10, -12, -14, -16,- 18 and -20 °C for at least 7, 5, 4, 4, 3, 3 and 2 consecutive days, respectively.

The first indicator represents a common measure of unfavorable cold-related conditions during the year. The second indicator delineates the upper bound of the cold weather persistence. The third indicator allows combined (intensity-duration) evaluation of ECEs in seven sevirity categories (very low/low/moderate/elevated/high/very high/extreme). The technology of calculation of the csd indicator allows the separation of events and the analysis of their temporal evolution.

The results are summarized for the SEE domain by the Köppen-Geiger climate zones. The significance and magnitude of trends in proposed indicators were assessed through the non-parametric Mann-Kendall’s test [69] and Sen’s method [70] using the R-package ‘trend’ [71].

Figure 3. Left panel: Medians of calculated percentiles of tn-distributions by stations for 1931–1980 and 1981–2020. Middle panel: Median of calculated lower percentiles of tn-distributions and 2-, 5-, 10-, 20- and 50-year return levels (estimated by GEV model) by stations for yearly minima (1931–1980). Right panel: Box plot of distributions of cold spells duration at the different tn-thresholds for 1961-2020. The BG36 dataset was used for all calculations.

Figure 3. Left panel: Medians of calculated percentiles of tn-distributions by stations for 1931–1980 and 1981–2020. Middle panel: Median of calculated lower percentiles of tn-distributions and 2-, 5-, 10-, 20- and 50-year return levels (estimated by GEV model) by stations for yearly minima (1931–1980). Right panel: Box plot of distributions of cold spells duration at the different tn-thresholds for 1961-2020. The BG36 dataset was used for all calculations.

2.2. Excess Cold Factor (ECF)

Several studies have recently analyzed extreme heat events in Europe and worldwide using the excess heat factor (EHF) that accounts for short-term acclimatization (e.g., [72,73,74,75]). This approach can be analogously applied to extreme cold events by using the excess cold factor [56,76,77,78,79]. ECF measures cold wave intensity at each location with an additional component to account for adaptation [45]. The units of excess cold are °C² as the ECF is a factorization of two cold indices representing the degree of acclimatization to cold and the climatological significance:

ECI_accl = (T_i + T_i-1 + T_i-2)/3 – (T_i-3 + … + T_i-32)/30

(1)

ECI_sig = (T_i + T_i-1 + T_i-2)/3 – T₀₅

(2)

T_i = (Tmax_i + Tmin_i)/2

(3)

ECF = -ECI_sig × min(-1,ECI_accl)

(4)

where T_i represents the average daily temperature for day i and T₀₅ is the 5^th percentile, calculated within the given base period (1961-1990 in this study); Tmax_i and Tmin_i are the maximum and minimum temperatures for day i.

ECF is negative (i.e., coldwave is in progress) when ECI_sig is negative (because the three-day period is cold in an absolute sense, being below the 5^th percentile for daily mean temperature), but additionally, the three-day period is substantially cooler than the preceding 30 days (ECI_accl < -1 °C). Positive values of ECF mean a coldwave is not in progress [45].

The ECF severity (ECF_sev) metric permits the comparison of cold wave events and their impacts across the world and can be readily implemented within early warning systems [80]. For given day i, ECF severity is defined as follows:

ECF_sevi = ECF_i /ECF_p15, if ECF_i < 0

(5)

ECF_sevi = 0, if ECF_i ≥ 0

where ECF_p15 denotes the 15^th percentile of all negative ECF values calculated for the base period.

Nairn et al. [80] defined three levels of ECE severity:

- L1 (low-intensity), when 0 < ECF_sev < 1;

- L2 (severe), when 1 ≤ ECF_sev < 3;

- L3 (extreme), when 3 ≤ ECF_sev.

The annual summaries of ECF_sev characteristics, such as mean and maximum magnitude and maximum number of consecutive days of excess cold, were calculated by stations and presented by the Köppen-Geiger climate zones.

2.3. Software Products Used in the Research

All calculations on observational data were made in the free software environment R, version 4.4.2 [81] and RStudio 2024.09.1 [82]. The map shown in Figure 1 is prepared using QGIS Desktop 3.28.5 [83].

3. Results and Discussion

In recent years, many studies have focused on the rising temperatures in Europe since the second half of the 20th century. Twardosz et al. [84] found that the rate of increase is almost twice as rapid from 1985 to 2020 as in the whole period after 1951 (0.027 and 0.051 °C/year, respectively). This aligns with our results concerning trends in annual mean temperatures in the examined domain for 1961-2020 and 1985-2020 (0.032 and 0.055 °C/year, respectively) [58]. However, the authors conclude that the rise in temperature during the winter should be considered much less steady than in the summer. Sippel et al. [65] found that the observed rapid increase in cold-season temperatures in the late 1980s was followed by relatively modest temperature trends thereafter, which can be interpreted as a manifestation of internal variability superimposed upon the long-term warming component.

Table 1 shows magnitudes of trends in annual mean (Ta), annual mean minimum (Tmn) and annual minimum temperatures (Tn) in the period 1961-2020 by climate zones. The largest magnitudes are obtained for Tn in the Cfb and Dfb climate zones. Most stations with larger trend magnitudes are located in the central and northern parts of the SEE domain. Unlike Ta, with statistically significant trends for all stations, the percentage of stations with significant trends in Tmn and Tn varies between 82% (Cfa) and 100%(Cfb) for Tmn and between 6-7% (Cfa and Csa) and 50%(Cfb) for Tn.

Figure 4 illustrates the long-term variation of the averaged by climate zone indicators ncd-5 and ccd-5 (calculated separately for each station), which reveals a downward tendency in all climate zones except for Csa. Both indicators show a statistically significant negative trend for 18% of stations in Cfa and up to 88% in the Cfb climate zone (Table 2). The average magnitude of the ncd-5 trends ranged between -3.1 (Dfb and Cfb) and -1.5 (Cfa) days/10 years, and the maximum magnitude of -4.5 days/10 years was reached in the Cfb climate zone. The ncd-5 maxima vary by stations between 0 and 118 days, with lower values occurring in coastal Csa stations. Around 50% of maxima were recorded in 1963, 10% in 1985, and 6% in 1987 and 1993.

The average magnitudes of the ccd-5 trend ranged between -1.0 (Cfb) and -0.4 (Cfa) days/10 years, and the maximum magnitude of -1.3 days/10 years was reached in the Cfb climate zone. The ccd-5 maxima vary by stations between 0 and 50 days, with the lowest values occurring again in Csa stations. The distribution of maxima by years differs from this one of ncd-5 – most of them have occurred again in 1963 (41%) and 1985 (13%), but also in 2012 (10%), 1964 (7%) and 1996 (6%).

The comparison of multiyear means of ncd-5 and ccd-5 for 1961–1990 and 1991–2020 is shown in Figure 5. The number of cold days decreased during the second period by 14% and 21% in the Dfb and Cfb climate zones, respectively, but there was no substantial change in the Cfa and Csa zones. The most significant change in ccd-5, both in absolute terms (-2.4 days) and as a percentage change (-24%), was registered in the Cfb climate zone.

Figure 6 presents the evolution of the cold spell duration indicator (csd) by categories, averaged by climate zones. The highest values for csd-5 were reached in 1963 (Dfb and Cfb) and 1985 (Cfa and Csa). The csd indicator distinguished the severity of extreme cold weather well in most of the SEE domain and marked some extremely severe cold events in the Csa climate zone, such as in 1968, 1985 and 2001 [85,86]. The calculated average trend magnitudes decreased with the csd-category increasing, varying from -2.7 to -0.2 days/10 years in Dfb and from -1.7 to -0.3 days/10 years in Cfb climate zone (Table 2).

Our results regarding the relatively weak-expressed trends in cold weather indicators can be viewed in the context of the less steady rise in winter temperatures compared to the rapid increase in summer temperatures in Europe [84]. Moreover, a recent study shows that little has changed with cold extremes since 1990, despite the intuitive suggestion that Northern Hemisphere cold extremes would become less intense and frequent not only with global warming but as the Arctic warms at an accelerated pace [87]. Concerning SEE, these findings align with the research of Buric and Doderovic [88], which reported a decrease in the number of cold nights (defined as days when tn < 10th percentile) during winter in Montenegro. This reduction varies from 0 to -2.5 days/10 years for stations with altitudes below 800 m, but it is statistically insignificant at 9 out of 11 stations.

Temperature thresholds underlying the proposed method for detecting and classifying ECEs through intensity-duration model include some physiological and technological limits, which crossing could have significant adverse health and economic consequences, such as increasing morbidity and mortality rates in cardiovascular, respiratory and some other illnesses [10,89], agricultural damage [15,17,90] and disruptions in energy production and supply [13,91]. On the other hand, ECF severity is helpful as an exposure index that scales well against human health impact [80] and an extreme cold indicator in climate studies [45]. Therefore, we choose the ECF_sev index to assess the suitability of the proposed threshold indicators. First, we examine the consistency between the most extended periods of consecutive cold days on a yearly basis, calculated as ccd-5 and ccd-ECF (where the latter is defined as the longest period of excess cold, ECF_sev > 0) using the SEE dataset. Long-term variations of ccd-5 and ccd-ECF, averaged by climate zones, are shown in Figure 7. It is evident that the -5 °C threshold is rather strict for Csa but quite soft for Dfb climates when assessing the excess cold that populations may experience. However, the correlation between ccd-5 and ccd-ECF remains high enough, even for the Csa zone.

During the studied period, each station reached a daily value of ECF_sev ≥ 3 at least once. The absolute maxima of ECF_sev range from 3 to 11.3, with the highest value recorded in December 2001 in Larissa, Greece. The number of days with excess cold varies slightly by station (12 to 17 days per year). Approximately 81% of the calculated daily values of ECF_sev are below 1.0 (severity level L1). Excess-cold days mainly occurred in January (45%) and February (25%). The estimated statistically significant trend magnitudes for ccd-ECF vary between -0.4 and -1.0 days/10 years, affecting 6-7% of stations in Cfa and Csa climates, 21% in Dfb and 75% in Cfb climates.

The indices averaged over the domain (ccd-5_d and ccd-ECF_d) indicate a strong linear relationship (r = 0.958). According to the definition of ccd-5, the longest cold period of the year falls typically into the coldest winter months, January and February. Thus, an interesting classification of winter severity can be derived by the scatter plot of ccd-ECF_d vs. ccd-5_d, accounting for the averaged yearly maxima of ECF_sev (Figur 8). Most years with mild winters (severity level L1) in SEE are associated with ccd-5_d below 7 days. For years with severe winters (L2), ccd-5_d is usually between 7 and 10. A few winters have prolonged cold periods exceeding 13 days, with ECF_sev above 2.5, categorizing them as significantly more severe. The only winter with the highest severity level (L3) is that in 1963, also characterized by the most prolonged cold period (ccd-5_d = 18 days). Generally, this categorization of winters is confirmed by various regional analyses and studies [e.g., 85,86,92,93].

Figure 9 presents an estimate of extreme cold weather in SEE based on the average ECF_sev, calculated for 1961-1990 and 1991-2020 on a yearly basis, over all periods of at least 7 days with tn ≤ -5 ℃ (in such a way being synchronized with the lowest severity category of the cold spell duration indicator, csd-5).

The results show that despite the decreasing number of ECEs in recent decades, the severity of cold weather in SEE has changed little regarding the excess cold that populations may experience. Percentage changes indicate that the Cfb climate zone is the most affected by regional warming.

The detailed comparison between yearly maxima of csd-categories (denoted by C1 to C7) and ECF_sev levels demonstrated good spatiotemporal conformity, as seen in the left panel of Figure 10. Although the relationship between them is complex, the percentage distribution of different categories by severity levels shows that the categories from C4 to C6 are more strongly related to L2. In contrast, the other categories are fuzzified between neighbor ECF_sev levels: C1 category mainly corresponds to L1, and over 30% of cases in C7 category align with L3 (Figure 10, right panel).

5. Conclusions

The subject of our study is the cold-weather climate phenomena related to the prolonged retention of extreme cold in Southeastern Europe. Although the cold spells are typical for the cold season in SEE, the research on their spatio-temporal changes is somewhat neglected compared to the extreme heat. This study can be viewed as a continuation of our previous work concerning the development of criteria and indicators for quantitative assessment of hot spells [58]. We used the same methodologies and datasets to investigate the spatiotemporal variations of extreme cold in the period 1961-2020. Three threshold indicators based on daily minimum temperature, presenting the unfavorable cold-season thermal conditions (ncd-5), cold weather persistence (ccd-5) and combined intensity-duration characteristic of cold spells in seven categories (csd-5/-10/-12/-14/-16/-18/-20) were analyzed and summarized by the zones of the Köppen-Geiger climate classification. All three indicators show statistically significant decreasing trends (p < 0.05, Mann-Kendall method) during the period 1961-2020 in the areas with typically cold winters in the northern and central parts of the studied domain. Comparing the two 30-year periods (1961-1990 and 1991-2020), we found that cold days decreased by 14% and 21% in the Dfb and Cfb climate zones, respectively, but there was no substantial change in the Cfa and Csa zones. The results show good consistency between the indicator of consecutive cold days (ccd-5) and the yearly longest period of excess cold. We presented the classification of winter severity in the SEE domain based on the averaged yearly maxima of ccd-5 and ECFsev. Excess-cold days in SEE are typically 12-17 per year, occurring mainly in January (45%) and February (25%). Most years with mild winters (severity level L1) in SEE are associated with a short retention of cold weather (below 7 days). A few winters have more prolonged cold periods exceeding 13 days, with ECFsev above 2.5, categorizing them as significantly more severe. The only winter with the highest severity level (L3) is that in 1963, also characterized by the most prolonged cold period (18 days).

Given the increased mortality risk, agricultural losses and disruptions in the energy sector linked to prolonged cold weather and extremely low temperatures, we demonstrated the suitability of the proposed intensity-duration model to classify extreme cold events in the studied domain based on the set of minimum temperature thresholds (-5, -10, -12, -14, -16, -18 and -20 °C) and corresponding duration thresholds (7, 5, 4, 4, 3, 3 and 2 consecutive days). The detailed comparison between categories of the cold spell duration indicator and yearly maxima of the ECFsev levels demonstrated good spatiotemporal conformity.

As far as we know, no survey has been conducted on the various aspects of using threshold indicators for a quantitative assessment of extreme cold events in the region. This research also delineates the aim of our future work to combine independently de-fined metrics, such as the threshold indicators presented here and the ECF metrics, when developing early warning systems.