Climate Projections and Time Series Analysis over Roma Fiumicino Airport Using COSMO-CLM: Insights from Advanced Statistical Methods

1. Introduction

The evaluation of the effects of climate change on airports is important to ensure flexibility and safety of the infrastructure of air travel. Airports are unsafe for temperature oscillations and extreme weather events, which can disrupt operation and compromise passenger safety. For example, floods of runways and terminal areas may cause expensive damage and flight delay. High temperature can affect the performance of the aircraft, requiring a long runway for safe takeoff. Additionally, disruption of ground transport and power supply caused by extreme weather can obstruct the functionality of the airport. Understanding these risks enables better planning and adaptation strategies, such as re -designing infrastructure, increasing drainage systems and updating the operational protocol. The methodologies for quantifying the potential effects of climate change on aviation sector are still limited, especially in the Mediterranean area [1], although this region has been defined as a climate hot spot [2]. As explained in [3], a first step required to define a framework for risk assessment is the identification of relevant climate change projections over the XXI century and their potential impacts on airport infrastructures. Global climate models (GCMs) are powerful tools for the analysis of large-scale climate processes, but they lack the high resolution required to capture fine-scale features for regional decision-making. Regional climate models (RCMs) are able to remove this gap by providing details into localized climate dynamics, such as temperature variations, precipitation and extreme weather events. These simulations support a deep understanding about the climate trends in specific areas, considering factors like topography, land use and local climate conditions. This localized information is essential for sectors such as agriculture, urban planning and disaster management, to drive informed strategies for mitigation and adaptation to climate change impacts. In the last decade, several initiatives e.g., the WCRP Coordinated Regional Downscaling Experiment (CORDEX) project [4] were established to provide a global coordination of RCM downscaling for climate change adaptation and impact assessment. Over Italy, detailed analyses of the impacts of climate change in the aviation sector are still limited. In [5], the expected climate extreme variability over the Capodichino airport, located in Naples (Italy) was analyzed. De Vivo et al. [6] applied a climate risk assessment framework to evaluate the impacts of extreme temperatures on several Italian airports.

ALBATROS (Advanced systems and solutions for better practices against hazards in the aviation system) is a Horizon Europe-funded project whose ambition is to maintain a high level of safety and resilience in aviation in view of extreme weather conditions. In this frame, CIRA has conducted high-resolution climate simulations using the RCM COSMO-CLM over a geographical area centered over Roma Fiumicino Airport (Italy), in order to provide climate projections according with IPCC RCP8.5 scenario, which is often adopted as a benchmark for understanding the potential impacts of unmitigated climate change. The Italian peninsula is a geographical area that requires high-resolution climatology (less than 10 km) due to its highly complex topography, which encompasses high mountain ranges and a diverse array of coastlines.

The main aim of this work was to provide climate projections over Roma Fiumicino airport, widening the range of projections already available in the literature for this area, but with higher resolution, employing the “state-of-art” latest version of COSMO-CLM. Due to the complex nature of climate phenomena, some studies [7] suggest performing specific analyses of climate time series by using advanced statistical methods, because of the lack of regularity of occurrences. An exhaustive review of the application of fractal analysis on wind speed time series is presented in [8]. The analysis of climate time series with tools like FFT or the Higuchi method [9] is critical for understanding the intrinsic complexity and dynamics of climate systems. These methods allow to quantify patterns and structures, revealing properties like persistence or fractality. For instance, the dimension obtained through Higuchi’s approach can distinguish between chaotic, random, or persistent trends, providing information on the underlying processes driving climate variability. This is essential for improving predictive models and assessing long-term trends, such as those related to climate change. The visualization and quantification of these dynamics can support policy decisions aimed at mitigating climate risks.

This paper is organized as follows: Sec 2 contains a description of the methodology employed, including RCM and simulations set-up; Sec 3 reports the comparison of the model results with respect to the observations whereas climate projections are analyzed in Sec 4. Discussion and conclusions are reported in Sec 5.

2. The Methodology Employed

2.1. The RCM and Simulation Set-Up

The three-dimensional, non-hydrostatic regional model COSMO-CLM [10] is the climate-adapted version of the operational weather forecasting model COSMO-LM [11], originally developed by the German Weather Service (DWD). It is based on the Navier–Stokes equations for a compressible flow and assuming a hydrostatic base state at rest. The atmosphere is modeled as a multicomponent fluid consisting of dry air, water vapor, liquid and solid water, all of which obey the ideal gas law. This fluid is influenced by gravitational and Coriolis forces. Small-scale phenomena that cannot be explicitly resolved are incorporated using statistical approaches implemented through various parameterization schemes. The COSMO-CLM version used in this work is 6.0_clm2. The development of the COSMO model has ended in 2021, in fact COSMO 6.0 was released on 15 December 2021 and it was the last release of the COSMO model. In this version, several diagnostic output variables were added, for example for wind sector classes and for sunshine duration. Also new tuning variables have been introduced for several parameterization schemes. A new additional hydrology scheme (groundwater and runoff) has been developed.

In this study, a preliminary experiment was carried out by using ERA5 as forcing [12], which is a high-resolution global reanalysis dataset by the European Centre for Medium-Range Weather Forecasts (ECMWF). It is characterized by spatial resolution of 0.25° (about 31 km) and temporal coverage from 1940 to the present. The quality of this dataset makes it ideal for providing perfect initial and boundary conditions to RCMs, enabling the study of localized climate phenomena with great accuracy. The simulation covers the period 1979-2010 and uses the computational domain (2°-22°E; 35-49°N) shown in Figure 1. The grid is characterized by a spatial resolution of 0.0715° (about 8 km) and has 185 x 195 grid points. The time integration is based on a third-order Runge–Kutta scheme, convection is parameterized using the Tiedtke scheme [13], while the number of vertical levels in the atmosphere and the number of soil levels are set at 40 and 9, respectively. The sponge zone is made up of 15 grid points.

The main experiment was driven by the GCM EC-Earth at 1.125° (about 128km) resolution. The GCM EC-Earth (Earth System Model) [14] is a state-of-the-art global climate model developed by a consortium of European institutions. It simulates the global climate system, integrating atmospheric, oceanic, land surface and sea ice processes. It was built upon the ECMWF Integrated Forecasting System (IFS) and was designed for weather prediction and climate research. Specifically, it is used for the study of past periods and future climate scenarios, contributing to several projects (i.e., IPCC assessments) and advancing the understanding of climate dynamics and variability. This experiment used a double nesting approach: the first simulation was performed over the period 1979-2050 using a coarser grid (resolution 0.44°, about 50 km) on the domain (30°W-54°E; 22-66°N) that includes 95 x 95 grid points; the number of vertical and soil levels are set at 40 and 9, respectively. The second simulation was performed using the same computational grid already adopted in the preliminary experiment (resolution 0.0715°), with initial and boundary conditions provided by the previous 0.44° simulation, over the period 1980-2050 (the year 1979 was removed to reduce the spin-up period influenced by initial conditions).

2.3. Observational Datasets

Model evaluation for the entire Italian domain was conducted using the E-OBS (version 29) gridded dataset [15] at resolution of 0.10°. E-OBS provides daily mean values for several variables, such as precipitation and temperature, over the period 1950–2023. It was developed by interpolating one of the most comprehensive collections of station data across Europe. The dataset advantage lies in its extensive spatial and temporal coverage, making it a standard reference for climate research and a widely used tool for evaluating RCMs across Europe. In this study, seasonal bias maps of mean, 2-m above ground air temperature and total precipitation have been analyzed, obtained by comparing the values of simulations output with the E-OBS considering the time period 1981–2010. To this aim, the model data were interpolated onto the dataset grid using the bilinear interpolation. Furthermore, in order to have a more detailed performance evaluation over the Roma Fiumicino airport, the annual cycle of the investigated variables have been computed for all the simulations performed considering the nearest grid point to the airport location. In this case, in addition to the E-OBS dataset, the station data from the SCIA system (National System for the collection, elaboration and diffusion of climate data) developed by ISPRA (Istituto Superiore Protezione e Ricerca Ambientale) [16] has also been used for the comparison.

2.2. Statistical Methods for the Analysis of the Results

The Fourier transform is a powerful tool for signal analysis, but it is subject to certain limitations. For example, a chaotic regime is characterized by a continuous power spectrum, however a Fourier spectrum that looks continuous cannot be automatically attributed to a chaotic signal, since it can be also representative of a quasi-periodic signal with a high number of frequencies (white noise) [17]. For this reason, a Fourier analysis (e.g., FFT) could not be sufficient to characterize a chaotic signal. The measurement of the fractal dimension of an attractor provides important information about the kind of temporal regime that characterizes a dissipative system [18], since the fractal dimension measures the number of relevant degrees of freedom. The Higuchi method [9] is a powerful tool that can be used to capture the intrinsic complexity of time signals. It is based on the construction of several time sub-series, starting from the original one (N elements). For a given parameter k, the algorithm splits the time series in sub-series whose distance is k, evaluates the length L of them and calculates the average value:

$$L_k\left(m\right)={\displaystyle \frac{1}{k}}\sum _{i=1}^{{\displaystyle \frac{N-m}{k}}}\left|x\left(m+i·k\right)-x\left(m+\left(i-1\right)·k\right)\right|·{\displaystyle \frac{N-1}{\frac{N-m}{k}·k}}$$

$$L\left(k\right)={\displaystyle \frac{1}{k}}\sum _{m=1}^k{L}_k\left(m\right)$$

With m =1,2… k. The process is repeated for several values of k and then the fractal dimension d is evaluated as the slope of the log-log relationship between the average length of the curve and k. Values of d close to 1 are representative of relatively simple and regular time series, while values of d close to 2 are related to complex and non-regular signals. For atmospheric variables, a fractal dimension close to 1 indicates a strong positive persistency, which is representative of regular time series that reflects foreseeable and stable climate cycles (e.g., annual seasonality). An intermediate value of d (about 1.3) indicates a limited persistency or a non-linear dynamic, in other words the series features a deterministic and a random component. Finally, higher values imply that the series shows non-regular fluctuations, associated with climate changes or extreme weather phenomena. The analysis of changes in the value of d over future periods with respect to past ones could help to record if the system is becoming more complex or less foreseeable.

3. Model Evaluation

This section is devoted to the quality assessment of climate model data. The mean 2-m temperature (T_2m) seasonal bias distributions over the period 1981-2010 are shown in Figure 2, for winter (DJF) and summer (JJA), respectively for ERA5 (a-c) and GCM EC-EARTH (b-d) driven simulations. The color scale indicates temperature bias in °C, with blue showing areas where the model underestimates temperature and red for overestimations. For ERA5 driven, a pronounced underestimation is recorded in winter over Alps, potentially due to difficulties in simulating complex topography or snow-albedo feedbacks. Smaller biases are observed in southern Italy and islands. In summer, a slight temperature overestimation is visible in southern Italy and Sardinia, which might reflect challenges in simulating heatwaves [19], soil moisture dynamics and urban heat effects [20]. An underestimation persists on northern and central Italy. For GCM EC-EARTH driven, panels b-d suggest that COSMO-CLM generally is affected by a cold bias across Italy, with intensity depending on the season. In winter, the cold bias is particularly pronounced in northern Italy whereas minor overestimations (red areas) are visible in limited regions, such as north-eastern Italy. The cold bias is stronger during JJA, suggesting that the model has difficulties to simulate warmer temperatures; coastal regions and islands show smaller biases compared to inland and mountainous areas. This general cold bias is inherited by the driving GCM, in fact the temperature underestimation by the EC-Earth model over Italy (and other specific regions) has been reported in some analyses and scientific publications. As regards bias in summer temperatures, recent studies (e.g., [21]) have indicated that EC-Earth tends to underestimate temperatures in Italy and other areas of the Mediterranean, and that this could be due to an incomplete representation of the feedback between dry soil and high temperatures. Regarding the topographic influence, the complex orography of Italy may be too simplified in the model, leading to an underestimation of temperatures in certain areas (e.g., plains or coastal regions). The latest version of EC-Earth has improved spatial resolution and the representation of certain key processes, reducing some biases. However, the extent of these improvements may vary depending on the model configuration (e.g., atmosphere-ocean coupling, resolution used). A recent study [22] analyzed the very high-resolution configuration of the EC-Earth model for “The High-Resolution Model Intercomparison Project (HighResMIP)”, comparing it with lower-resolution versions. The results indicate that, despite improvements in the representation of certain climatic phenomena, biases persist in the simulations of summer temperatures in Europe, highlighting the need for further refinements in the model.

Figure 3 shows the seasonal precipitation bias (mm/month) of the ERA5 driven simulation (a-c) compared to the E-OBS dataset for the period 1981–2010, with blue indicating overestimation and red indicating underestimation. In winter, a positive bias is recorded in the northern part of the domain, whereas negative biases affect Apennines and southern Italy. In summer, precipitation is slightly overestimated in the north, which might be related to excessive convective activity simulated by the model. For GCM driven simulations (b-d) results highlight strong spatial and seasonal variations, with notable differences in northern Italy compared to central and southern regions. In DJF, an overestimation is observed in the northern regions, particularly along the Alps whereas some underestimations occur in southern and central Italy. In JJA, northern Italy shows overestimation, although the magnitude is reduced compared to spring. Underestimations are more widespread across central Italy, likely due to the model underestimating convective precipitation during the summer.

The annual cycles of 2m daily mean (T_2m), maximum (Tmax_2m), minimum temperature (Tmin_2m) and total precipitation (Tot_Prec) for Roma Fiumicino airport location are shown in Figure 4. Data are provided by ERA5 and EC-Earth driven simulations at both the resolutions considered (0.44° and 0.0715°), averaged over the period 1981–2010, obtained considering the nearest grid point to the location considered. They have been compared to annual cycles provided by the E-OBS dataset (nearest grid point) and ground station data from the SCIA system. The ERA5 driven simulation is characterized by great accuracy, while EC-Earth driven exhibits a systematic cold bias for mean, maximum and minimum temperatures, with the bias increasing during warmer months. Generally, the coarse resolution underperforms compared to the finer one, with the exception of the maximum temperature in summer. Regarding precipitation, during spring the EC-Earth driven simulations are characterized by a general overestimation, more pronounced for the coarser one, whereas the ERA5 driven simulation shows a slight underestimation. From September to December, instead, all the simulations underestimate the observed values, especially the run at 0.0715° forced by EC-Earth. The annual cycles have been computed also considering the average value over a 3x3 grid box centered on the nearest point to the airport location. The results (not shown) are generally similar to the ones obtained considering the single closest point, with the exception of the simulation at 0.44° resolution, since the analyzed area is much larger due to the greater distance between the points included in the box.

4. Analysis of Climate Projections

Climate projections have been performed under the IPCC RCP8.5 (Representative Concentration Pathway 8.5) scenario [23]. It is a high greenhouse gas concentration scenario, often referred to as a “business-as-usual”, which assumes a continued growth in fossil fuel use and limited mitigation efforts, resulting in a radiative forcing of 8.5 W/m² by 2100. This scenario projects significant increases in global temperatures, rising sea levels and increasing extreme weather events, serving as a benchmark for understanding the potential impacts of unmitigated climate change.

Figure 5 shows the yearly time series (5-year running mean) covering the period 1981-2050 of daily (mean, maximum and minimum) 2m temperature and total precipitation at Rome Fiumicino (nearest grid point). For temperature, a positive trend is projected (0.028, 0.027 and 0.029°C/year for mean, maximum and minimum respectively), which could exacerbate heatwaves, drought conditions and energy demands for cooling, as well as poses health risks. For precipitation, a modest positive trend is also detected (0.109 mm/year), with annual precipitation fluctuating between 490 mm and 673 mm; this large variability could result in more intense rainfall events, with possible implications for flooding, droughts and water resource management.

4.1. Analysis of Daily Time Series and FFT

Daily time series of 2m maximum and minimum temperature, total precipitation and 10m wind speed, covering the periods 1981-2010 and 2021-2050, have been processed by using an FFT to produce a power spectrum of the corresponding variables at regularly spaced bins, on order to measure the amounts of variability occurring in different frequency bands. The frequency resolution is 1/T₀ = 0.912 × 10−⁴ (given that T₀ is 10957 d). The values of the main periods recorded (in days) for the variables considered are shown in Table 1. Figure 6 (a-c) shows the time series of daily 2m maximum and minimum temperature over the periods considered, while Figure 6 (b-d) shows the corresponding mean power spectra. For Tmax_2m, a regular pattern in the time series can be observed, while three main peaks are visible in the power spectrum (f = 2.74 × 10−³, 5.48 × 10−³, 8.21 × 10−³), corresponding to periods of 365, 182 and 121 d respectively. The first one corresponds to the annual cycle of temperature, dominated by seasonal variations. The second one refers to variations that are repeated every 6 months, reflecting temperature changes related to seasonal transitions (e.g., winter-spring and summer-autumn). Moreover, higher harmonic frequencies may arise if heating and cooling are not uniform over the year. Both past and future periods show a similar behavior at low frequencies, while at high frequencies (>0.005) a slight attenuation of the amplitude is observed in the future period. This could be related to a reduction in rapid (day-to-day) fluctuations, linked to a greater persistence of warm events or a more uniform climate on short time scales. For Tmin_2m, in both periods a similar trend is recorded for low frequencies (T= 365 d) with comparable amplitudes. At higher frequencies, amplitudes are slightly different for the two periods, suggesting a potential reduction of the fast variations in Tmin_2m in the future, with a more stable climate.

Figure 7 (a-c) shows the time series of daily total precipitation and 10m wind speed over the periods considered, while Figure 7 (b-d) shows the corresponding mean power spectra. For precipitation, significant peaks in the amplitude spectrum are observed around low specific frequencies. These peaks indicate the presence of dominant periodic cycles in precipitation, probably related to seasonal or interannual climate phenomena. The main peaks are present in both past and future periods, suggesting that some periodic features persist over time. The future series generally presents a slightly larger amplitude than the historical series in some frequencies, especially around the main peaks, indicating an intensification of some cyclical components of precipitation in the future. The spectral noise appears slightly more pronounced in the future, reflecting greater climate variability in future precipitation. These results project a change in the intensity or periodicity of future precipitation events. In particular, an increase in spectral amplitude at certain frequencies may reflect more intense weather events, considering that any changes in the identified cycles could be associated with global climate processes, such as changes in atmospheric circulation. For 10m wind speed, the same peaks in the amplitude spectrum are observed for both time series around low specific frequencies. Moreover, the spectra (log-log representation) at high frequencies generally show a declining trend with noise superimposed. At higher frequencies, the observed distribution has been compared with the power law described by Kolmogorov for homogeneous isotropic turbulence (depicted as blue line in the figure). According to Kolmogorov’s theory [24], the velocity spectrum of turbulence can be divided into distinct frequency ranges, often classified as the turbulence source region (large scale), the inertial subrange (medium scale), and the dissipation region (small scale). This law [25] asserts that, assuming constant energy dissipation rate, the energy of turbulent flow depends solely on the scale of turbulence raised to the -5/3 power. From the figure, it is clear that at higher frequencies, the wind power spectrum aligns closely with the Kolmogorov distribution. However, it is crucial to emphasize that the Fourier power spectrum is a second-order statistical measure, providing insights into medium-level fluctuations but not fully characterizing the scaling process.

4.2. Evaluation of Fractal Dimensions

The values of fractal dimensions (d) represent a measurement of the complexity and variability of time series. Values of d associated with the daily time series already described were calculated by means of the Higuchi method, which allows the assessment of irregularity and seasonality of the variables considered, influenced by several atmospheric systems (such as topography) that behave at various scales and periodicity. Figure 8 shows the log-log plots of the relationship between L(k) and k for the four variables considered, for both past and future periods, which is the basis of the evaluation of the fractal dimension. Table 2 contains the corresponding values of the fractal dimensions. For Tmax_2m, the complexity will slightly increase in the future, which could be related to a growth of the daily thermal variability or to a larger number of extreme events. A larger increase is recorded for Tmin_2m, which could be associated with phenomena such as warmer nights or a larger nocturnal variability. For total precipitation, the value of d remains almost unchanged, implying that climate changes do not have a significant impact on the temporal variability of total precipitation. For 10m wind speed, the expected variation of d is small too, revealing that changes in the temporal structure of wind are limited, with a stable complexity over the time. Overall, the expected increases of d for temperature in the future are related to the effects of climate change, which introduce new patterns in the temporal distributions of extreme events, such as heat waves.

5. Discussion and Conclusions

In the present work, climate simulations performed with the RCM COSMO-CLM model over a geographical area centered over Roma Fiumicino airport have been performed and analyzed. The study emphasized the critical role of high-resolution regional climate modeling in understanding localized climate impacts. As regards the ERA5 driven simulation, a general slight underestimation of temperature is evident in all the seasons, stronger in winter. Negative precipitation biases in central Italy suggest difficulties in capturing Mediterranean storm dynamics and convective rainfall. As regards EC-Earth driven simulation, COSMO-CLM consistently underestimates temperatures in all seasons, which is mostly inherited by the global model, but could be also related to errors in simulating topography or snow cover. These biases reveal challenges in accurately capturing temperature dynamics, possibly due to parameterizations of surface energy fluxes, vegetation or atmospheric processes. Overall, the present analysis highlights that in some areas model improvements are needed, particularly in handling winter cold extremes, summer heat events and the representation of topographic influences, implying that bias correction techniques might be required when using RCMs for impact studies or decision-making processes.

For the specific area of Roma Fiumicino airport, the analysis under RCP8.5 scenario projects a consistent increase in daily mean, maximum, and minimum temperatures, with trends of 0.028°C/year, 0.027°C/year and 0.029°C/year, respectively, which could exacerbate heatwaves, increase energy demands for cooling and heighten health risks. Warmer nights could lead to reduced nocturnal cooling, affecting energy demands and increasing the heat stress on ecosystems and urban populations, while higher maximum temperatures could intensify heatwaves and increase risks for health, agriculture and energy systems [26,27]. Similarly, a modest positive trend in annual precipitation (0.109 mm/year) is detected, highlighting increased variability and a potential rise in extreme rainfall events. The area may experience a combination of warmer conditions and precipitation variability, leading to urban planning challenges such as heat stress, flood control and infrastructure resilience [28]. These changes underline the importance of developing adaptation strategies for aviation (e.g., runway reinforcement, heat-tolerant materials, fuel efficiency innovations) and other sectors, such as agriculture, water management, urban planning and public health.

These outcomes are consistent with previous studies [29], confirming the broader trends of climate change impacts over the Mediterranean regions. Moreover, the results provide novel contributions by offering high-resolution projections specific to the Fiumicino airport area, enabling a more detailed assessment of potential risks and adaptation needs. Furthermore, the use of advanced statistical methods, such as the Higuchi method and FFT analysis, revealed insights into the fractal dimensions and frequency components of climate variables, showing an increasing complexity and variability of future climatic patterns. Then, the methodological advancements demonstrated in this study provide a framework for assessing climate impacts in other regions with similar topographical and climatic features. Nevertheless, some limitations should be acknowledged, such as the inherent uncertainties associated with climate modeling and the reliance on a single IPCC scenario. To overcome these limitations in the future, other methodologies will need to be considered. For example, the integration of fractal analysis with AI techniques, such as artificial neural networks and neuro-fuzzy systems, has proven effective in improving wind speed prediction accuracy. A study by Petković et al. [30] highlighted the synergy between fractal extrapolation and AI methods in addressing the complexities of wind speed time series. Fractal metrics, when integrated with machine learning, could enhance model performance by quantifying complexity and refining data decomposition [8]. These metrics aid in the tuning of hyperparameter and the trend identification, improving predictive accuracy and reducing the computational overhead. In applications like wind speed forecasting, fractal analysis efficiently allocates resources by identifying components that significantly impact forecasting errors. Furthermore, the integration of physics-informed models, combining domain expertise with fractal techniques and machine learning, could offer a more accurate approach to climate change projections.