On the role of wave climate temporal variability in bias correction of GCM-RCM wave simulations

This work presents the performance analysis of a multi-model ensemble of wave climate projections in the Mediterranean Sea against hindcast data. The wave projections were developed with the numerical model Wavewatch III forced by surface wind fields of 17 EURO-CORDEX GCM-RCMs providing time series of the main wave parameters on a 3-h and 10-km resolution. The performance of the wave GCM-RCM simulations during the baseline period (1979–2005) was assessed by means of the deterministic metrics RMSE and Bias. Different bias correction methodologies were analyzed by means of the application of the widespread Empirical Quantile Mapping method considering different time periods of significant wave height in order to analyze the ability of the bias-correcting methods to capture the different wave climate temporal scales ranging from storm events, monthly, seasonal and interannual variability. The results show that the use of the EQM method for the full-time series without taking into account other timescales, can lead to increased biases in some regions and seasons and that the use of time-dependent bias-correction techniques leads to an improved accurate characterization of biases considering the interannual temporal variability of significant wave height. More specifically the use of the EQM method for monthly data provides a good performance in capturing the correlation and interannual temporal variability of wave climate.


Introduction
Coastal regions are amongst the most vulnerable areas subjected to climate changes impacts due to extreme coastal water level (ECWL) events, which include the combination of different physical processes acting at varying spatial and temporal scales, such as relative sea-level change, tides, waves, storm surge, and swash (Woodworth et al. 2019;Gregory et al. 2019). Swell and wind-sea waves are important contributors to ECWLs and crucial for coastal impact assessments and damage to coastal infrastructure. Despite being the predominant factor in many nearshore processes and coastal impacts, the contribution of waves has often been neglected given the difficulty to predict the large longterm variability in surge and wave conditions and the lack of local coastal topography and bathymetry data (Bricheno and Wolf 2018;Lira-Loarca et al. 2021a;Almar et al. 2021;Nicholls et al. 2021).
Projected changes in wave climate arising from changes in atmospheric circulation vary on a regional and local scale. Global and regional climate models are the primary tools to investigate the response of the climate system to different emission scenarios and are used as forcing conditions for modeling wave climate projections. Regional Climate Models (RCMs) are used to downscale coarse Global Climate Model (GCM) data to a finer resolution for its use in regional 1 3 and local vulnerability, impact, and adaptation assessments (Giorgi 2019;Lemos et al. 2020a).
Studies that present future global or regional wave climate changes usually do so in terms of relative changes in statistics between the future and baseline simulations of a given GCM-RCM or a multi-model ensemble of wave projections modeled using forcing atmospheric conditions from the corresponding GCM-RCM models (Bricheno and Wolf 2018;Morim et al. 2019Morim et al. , 2020. However, GCM and RCMs atmospheric simulations present systematic biases due to discretization and spatial resolution, simplified physics or parameterizations, internal variability and downscaling processes (Christensen et al. 2008;Teutschbein and Seibert 2012). Wave climate projections forced by GCM-RCMs wind field data, therefore inherit their systematic error and bias adjustment methods should be employed to address them and allow quantitative coastal impact projections (Lemos et al. 2020a, b). Additionally, wave climate presents a high temporal variability with timescales that range from decadal, interannual, seasonal, and monthly variability to hours, seconds, and minutes via storm events, swell and wind-waves and swash processes linked to individual waves (Melet et al. 2018;Lira-Loarca et al. 2021a).
Bias-correction techniques are extended in studies dealing with climatic and hydrological variables such as precipitation and temperature but their applications to wave climate still remains an open and challenging issue due to the multivariate behavior, varied temporal and spatial variability of waves as well as the influence of complex topography and thermal gradients on surface winds fields used as forcing conditions for wave models (Lemos et al. 2020a;Costoya et al. 2020). Several bias correction methods have been developed ranging from linear scaling of the mean through additive or multiplicative correction factors, parametric (distribution mapping) or non-parametric (empirical) quantile mapping that allows to correct the distribution function of simulated data to agree with the distribution function of the observed data to multivariate techniques (Teutschbein and Seibert 2012;Hempel et al. 2013;Parker and Hill 2017;Holthuijzen et al. 2021). The use of these methods leads to an improvement in the statistical distribution of the simulations to match that of observed data but can lead to misrepresentation of the temporal, spatial, and multivariate characteristics of the physical processes. Therefore, different techniques have been developed to account for the multivariate, temporal, and/or spatial variability of a given atmospheric simulation (Haerter et al. 2011;Maraun 2013;Vrac and Friederichs 2015).
The empirical quantile mapping (EQM) method, also known as distribution mapping, probability mapping, or quantile-quantile mapping, is one of the most widely used bias adjustment techniques in atmospheric variables due to its flexibility and ability to account for the extreme values of the distribution (Teutschbein and Seibert 2012;Déqué 2007). Although the application of different bias correction techniques to wave climate projections is still limited to a few studies, the EQM method has proven to be an effective tool to account for biases in wave projections data (Parker and Hill 2017;Lemos et al. 2020a, b;Lira-Loarca et al. 2021a;Lobeto et al. 2021). Nonetheless, to the authors' knowledge, wave bias correction studies have been done under a stationary approach without considering the temporal variability in wave climate. Therefore, there is still a need to address different bias correction techniques that take into account the varied temporal periods present in waves and estimate if stationary bias-adjustment methods could lead to increase bias in different timescales.
This work presents a multi-model ensemble of wave climate projections in the Mediterranean Sea modeled with the numerical wave model Wavewatch III forced by wind surface fields from 17 EURO-CORDEX GCM-RCMs. Their performance was evaluated against a validated hindcast during the baseline period . Then, different biascorrections methods were applied and the performance of the adjusted data was analyzed focusing on the ability of the GCM-RCMs to capture the multi-temporal variability present in wave climate. This paper is organized as follows. Section 2 provides a description of the wave hindcast and GCM-RCM simulations, the metrics used to evaluate their performance, and the bias adjustment methods applied to GCM-RCM simulations in the Mediterranean Sea. Then, Sect. 3 present the results of the performance of raw GCM-RCM data against hindcast and the application of different bias correction methodologies to improve the temporal representation of wave data. Finally, Sect. 4 presents the discussion and the main conclusions from this study.

Methods and data
This section provides a description of the wave hindcast in the Mediterranean Sea developed with the numerical model Wavewatch III and the multi-model ensemble wave climate projections under climate change scenario RCP8.5 with the same setup. Additionally, the different bias-adjustment techniques analyzed in this work are presented as well as the metrics used to evaluate their performance.

Wave hindcast and projections in the Mediterranean Sea
This work uses the wave hindcast developed by the Meteocean research group 1 of the University of Genoa (Italy) providing high-resolution wave climate data from 1979 to 2020 in the Mediterranean Sea using a regular grid of longitude 0.127 • and latitude 0.09 • , corresponding to ≈ 10 km (Mentaschi et al. 2013(Mentaschi et al. , 2015Cassola et al. 2016;Besio et al. 2016). It uses the third-generation wave model Wavewatch III (version 5.16; The WAVEWATCH III®Development Group 2019, hereinafter WW3) with the source terms of growth/dissipation ST4 (Ardhuin et al. 2010;Rascle and Ardhuin 2013). This database provides hourly time series of the main wave climate integrated parameters. In this study, we have focused on the significant wave height, H s . Using the same WW3 configuration, a multi-model ensemble of wave climate projections was obtained, forced by surface wind fields of seventeen different Euro-CORDEX (Jacob et al. 2014(Jacob et al. , 2020) models (GCM-RCM combinations, Table 1) with 6-h temporal resolution and 0.11 • ( ≈12.5 km) spatial resolution. Wave climate simulations for each GCM-RCM were obtained with a 3-h resolution for the base period from 1970 until 2005 and for the RCP8.5 high-emission scenario extending from 2006 until 2100 (De Leo et al. 2021;Lira-Loarca et al. 2021b). Additional details of the definition and performance of the different RCMs used in this work can be found on Strandberg et al. (2014) for the Rossby Centre regional climate model RCA4, Will et al. (2017) for the CLM-Community CCLM4-8-17 model, Christensen et al. (2007) for the Danish Climate Centre regional climate model HIRHAM5 and Leutwyler et al. (2017) for the COSMO-CLM accelerated version COSMO-crCLIM-v1-1.

Performance of GCM-RCM wave data
The performance of a GCM-RCM was evaluated by means of the Bias (Bias, Eq. 1) and root-mean-square error (RMSE, Eq. 2) metrics for the entire Mediterranean basin.
where Y is the analyzed variable, N is the length of the dataset and the subscripts hind and RCM hist correspond to the hindcast and GCM-RCM baseline simulations, respectively, from 1979 to 2005 (27 years). In this work we have analyzed the Bias and RMSE of the variables (Y): (1) significant wave height monthly means and (2) significant wave height monthly maxima. In this study we have analyzed 27 years of data and therefore the metrics are calculated for N = 27 years ⋅ 12 months = 324 or N = 324∕4 = 81 values, when performing a seasonal analysis.

Empirical quantile mapping method
The empirical quantile mapping (EQM) method (also referred to as probability mapping or quantile-mapping method) corrects the distribution function of the GCM-RCM projections to agree with the hindcast distribution as Déqué (2007), Teutschbein and Seibert (2012), where X * is the bias-adjusted value of the raw GCM-RCM variable, X and F hind and F RCM hist are the empirical distribution functions of the hindcast data and GCM-RCM, respectively, during the baseline period .
In this work, the EQM method is applied to significant wave height ( H s ) data and the empirical distribution functions are determined by a series of regularly spaced quantiles from the 1st to the 99th and additional quantiles from the 99th to 99.99th to account for the upper tail of the distribution. Then, a correction is applied for each quantile and linear interpolation between them is done to obtain the bias-corrected time series of the variable (Lemos et al. 2020a).

Temporal empirical quantile mapping methods
Wave climate presents a large range of temporal variability varying from sea-state, storm events, monthly, seasonal to interannual and mutidecadal timescales. In order to analyze the importance of these temporal dependencies and propose temporal bias-adjustment methods that account for them, we propose the same approach as the EQM method except that the correction is done by grouping the hindcast and GCM-RCM historical data into different timescales. Then, for each group, the bias-corrected quantiles are obtained following Eq. (3) and the bias is considered stationary within each selected temporal scale. In this work, we have analyzed the following temporal-bias corrections methods: • EQM-seasonal: The data is grouped by seasons considering December-February (Winter), March-May (Spring), June-August (Summer), and September-November (Fall). • EQM-monthly: The data is grouped by months, therefore, obtaining twelve different groups to be corrected. • EQM-dayofyear: The data is grouped by day-of-year considering a 365-days year.

Performance of bias-corrected wave data
The performance of the bias-adjusted data against hindcast was evaluated by means of the RMSE and Bias metrics presented in Sect. 2.2. In order to assess the skill of the proposed bias correction methods, the Skill Score index (SS) is used to quantify the improvement of the bias-corrected GCM-RCM with respect to the raw GCM-RCM data.
where E * is the value of the error metric for the bias-adjusted data, E raw corresponds to the error for the raw (not biascorrected) GCM-RCM data, and E opt is the optimal value for the corresponding metric. In the case of the root-meansquare-error (RMSE), a perfect fit to observations is given by E opt = 0 . Therefore, the higher the skill score, the better the performance of the bias-adjusted data.

Fig. 4 Bias between the raw GCM-RCMs and hindcast of the monthly maxima of significant wave height, H s [m]
Regarding the performance of the Bias metric, a simple delta method was used, where Bias * is the Bias value for the bias-adjusted GCM-RCM versus hindcast data and Bias raw corresponds to the error for the raw GCM-RCM versus hindcast data. Then, negative ΔBias values indicate an improvement in the performance of the bias-corrected GCM-RCM with respect to the raw GCM-RCM.

Performance of GCM-RCMs
As a first step of the study, the performance of the different GCM-RCMs was evaluated against hindcast data by means of the seasonal RMSE and Bias for the monthly means and maxima values of significant wave height H s . Figures 1 and 2 present the RMSE of the monthly means and maxima, respectively, for all analyzed GCM-RCMs. Figures 3 and 4  It can be observed in Fig. 1 that, in general, all GCM-RCMs, except CCLM4-CanESM2 present a similar RMSE spatial distribution for the monthly means ranging from 0.2 to 0.7 m, with the maximum values depicted in the Western Mediterranean. CCLM4-CanESM2 presents higher values up to 0.9 m in the Tyrrhenian Sea. Regarding the monthly maxima ( Fig. 2) all GCM-RCMs present the same spatial distribution between them and with respect to historical values of significant wave height and mean wave energy flux (Besio et al. 2016;Sartini et al. 2017). Higher values are depicted for the Western and Central Mediterranean and the lowest values are observed in the Adriatic and Aegean Seas as well as near-coast regions. Then, as expected, higher RMSE for the monthly maxima are obtained in places with higher waves where larger errors are expected to have a larger effect on the RMSE values and bias-correction techniques can be used to improve the models' performance in the upper quantiles, crucial for coastal impact assessments. In general, for the monthly maxima, RMSE values of [1.15, 1.6] m are obtained for most of the Central and Eastern Mediterranean basins as well as the Tyrrhenian and Alboran Seas for all GCM-RCMs, with the exception of the models CCLM4-MIROC5, HIRHAM5-NorESM1-M, and HIRHAM5-CNRM-CM5 for which the Central Mediterranean presents higher RMSE = [1.6, 1.75] m. In all models, the Western Mediterranean and Balearic Sea present the highest RMSE = [1.75, 2.05] m for the monthly maxima. The    the distribution. The highest underestimations are obtained for the HIRHAM5 simulations in the Gulf of Sidra and the Levantine Sea. This indicates that the performance of the GCM-RCMs simulations is highly dependent on the RCM regardless of the GCM forcing.
Given the dependence of the results to the RCM simulation and results obtained, the performance of the models CCLM4-MIROC5, RCA4-HadGEM2-ES, and HIRHAM5-CNRM-CM5, for the rest of the seasons, is further analyzed. Figures 5 and 6 present the metrics for the monthly means of significant wave height for models CCLM4-MIROC5 and RCA4-HadGEM2-ES, respectively, whereas Figs. 7 and 8 present the metrics for the monthly maxima for models CCLM4-MIROC5 and HIR-HAM5-CNRM-CM5. The highest RMSE are obtained during winter (Figs. 5A, 6A, 7A, 8A) in the Western Mediterranean, due to the higher occurrence of extreme waves and therefore higher errors, with values of up to 0.7-0.8 m for the monthly means and 2.5-2.75 m for the monthly maxima. This region, where higher waves are expected (Sartini et al. 2017), presents the largest RMSE for the remaining seasons. The RMSE of the monthly means presents the lowest values during summer for both analyzed GCM-RCMs (Figs. 5E, 6E), for which the lowest hindcast wave energy is reported in previous studies (Besio et al. 2016;Lira-Loarca et al. 2021b). The CCLM4-MIROC5 simulations present, for the monthly maxima, maximum RMSE values of 2.5 m for the Central Mediterranean during summer (Fig. 7E) and for the Western Mediterranean during winter (Fig. 7A) maxima, the remaining seasons present values up to 2.1 m and a higher spatial variability (Fig. 7C, G). Regarding the Bias behavior, a higher spatial variability is noted for all seasons (Figs 5, 6, 7 and 8; panels B, D, F and H) and different behaviors are obtained for the different analyzed GCM-RCMs. CCLM4-MIROC5 presents overestimations in the central parts of the Mediterranean for both monthly means ( Bias = 0.3 m) and maxima ( Bias = 1.1 m) and underestimations in the near-shore regions. The root-mean-square-error is more sensitive to outliers due to the second-moment sensitivity to large errors, therefore it can be explained that the high RMSE for the monthly maxima could be due to a poorer ability of the models to adequately capture the behavior of extreme waves.

Performance of bias correction methods
In order to assess the importance of multi-temporal variability of wave climate, the following four bias correction  methods are applied to the GCM-RCMs raw significant wave height, H s , data: • EQM: Empirical Quantile Mapping (Eq. 3) using a range of quantiles from the 1st to the 99.99th for the baseline time series of wave data   HIRHAM5-CNRM-CM5 presented large metrics and varied behaviors, all the results presented henceforth will be regarding the performance of the bias-adjusted data for these models allowing to assess the robustness of the biasadjusting methods for different GCM-RCMs simulations. Figures 9 and 10 present the the Winter (December-February) RMSE of the H s monthly means with respect to the hindcast of the raw and bias-adjusted data for the CCLM4-MIROC5 and RCA4-HadGEM2-ES simulations, respectively, for the different bias-correction techniques. Figure 11 presents the RMSE results for Fall (September-November) for RCA4-HadGEM2-ES and Fig. 12   Detailed caption as in Fig. 9 rows present, on the left (panels B, D, F, and H), the RMSE of the bias-adjusted data with respect to hindcast and, on the right (panels C, E, G, and I), the skill score that assesses the performance changes. Positive SS values indicate an improvement in the performance of the bias-adjusted data with respect to raw data. It can be observed that in Winter for CCLM4-MIROC5 (Fig. 9), for the Western Mediterranean and Gulf of Lion, which presented the highest RMSE, all four methods present an improvement in biases with maximum RMSE values of 0.6 m, leading to a skill score SS ≈ 10% . The largest improvement with respect to raw data is obtained, for all methods, for the Central and Eastern Mediterranean and Alboran and Ligurian seas with SS ≈ 20-25%. In line with this, for RCA4-HadGEM2-ES winter conditions (Fig. 10) the largest improvements are obtained for the Alboran and Adriatic Seas with SS ≈ 15-20% whereas the Western Mediterranean presents very small improvements ( SS ≤ 5% ) for time-dependent methods (Fig. 10E, G, I). In general, although all bias correction methods present overall similar behaviors, for CCLM4-MIROC5 the EQM (full-time series) and EQM-seasonal present slightly larger skill scores SS = 12-16% (Fig. 9C, E) with respect to SS = 8 − 12% for EQM-monthly (Fig. 9G) since the use of the complete time series provides a large  Fig. 12 Summer (June-August) Bias of the H s monthly means, between the raw CCLM4-MIROC5 and hindcast (A) and different bias-adjustment simulations (B, D, F, H) and delta, ΔBias , using the raw model as reference (C, E, G, I). Detailed caption as in Fig. 9 sample of data for the less probable higher quantiles, which leads to a smaller variance of the errors. Nonetheless, for RCA4-HadGEM2-ES, which presented lower RMSE values, the EQ-full did not provide significant improvements whereas the largest improvements were obtained in the Eastern and Levantine Seas where RMSE ≈ 0.4 were observed. This behavior is also observed for the fall months (September-November) for RCA4-HadGEM2-ES (Fig. 11) where the temporal-dependent bias-adjusted data (Fig. 11, panels D, F and H) perform better than the EQM method with the full-time series (Fig. 11B).
Furthermore, the advantages of the time-dependent methods are noticeable when examining the summer months (June-August) results. Figure 12 presents in the first row (panel A), the Bias of the monthly means between the raw CCLM4-MIROC5 data against hindcast, which is used as a reference to assess the performance of the different bias correction methods. Then, the following rows present, on the left (panels B, D, F, and H), the Bias of the biasadjusted data with respect to hindcast and, on the right (panels C, E, G, and I), the ΔBias that assesses the performance changes. Negative ΔBias indicates an improvement  in the performance of the bias-adjusted data with respect to raw data. The best performance is given by the EQMseasonal, EQM-monthly and EQM-dayofyear methods with |Bias| ≤ 0.1 (Fig. 12D, F, H) m for the entire Mediterranean basin, therefore reducing the spatial variability of bias given by the raw data and providing improvements in the performance of the bias-adjusted data against raw data. Indeed, the larger improvements in the performance of the temporal-dependent bias correction are obtained for the Balearic and Levantine Seas with ΔBias ≈ −0.4 m and slightly lower improvements ΔBias ≈ −0.1 in the Ionian and Adriatic Seas.
The EQM method using the complete time series leads to lower improvements in the bias in the Balearic and Levantine Seas with ΔBias ≥ −0.1 and even a loss in performance in some areas of the Alboran and Ligurian Sea and the Strait of Sicily with ΔBias ≥ 0.05 (Fig. 12B).
Regarding the results of the monthly maxima, Fig. 13 depicts, for CCLM4-MIROC5, the Winter RMSE where a higher spatial variability in the skill score results is obtained. The largest improvements are obtained for the EQM-seasonal and EQM-dayofyear methods (Fig. 13D-H)   performance is obtained in the southern part of the Strait of Sicily with SS ≈ −20% for the EQM-dayofyear method. The EQM-full (Fig. 13B) which considers a larger sample of extreme values presents overall increases although lower in magnitude that the EQM-seasonal method. Figure 14 presents the RMSE for CCLM4-MIROC5for Summer which presents SS ≈ 70% in the Central Mediterranean region for the time-dependent methods leading to RMSE ≈ 0.5 m with respect to the raw RMSE ≈ 2.3 m values. It is noted that the EQM-monthly and EQM-dayofyear methods also provide improvements in the Western Mediterranean with SS ≈ 20% where the bias-adjustment provided maximum values of RMSE ≈ 1.6 m. The EQM-full presents low improvements in this region and some loss of performance in the Alboran, Adriatic, and Aegean Seas with SS ≈ −10% . Indeed, as the EQM with the full-time series performs the bias correction using a larger sample, the adjustment done for some quantiles includes data from different seasons that could lead to an increased bias when considering seasonal-dependent res ults.
In order to better understand the behavior of the different bias-adjustment methods and the ability to capture wave temporal variability, we analyze the results in selected locations for the different representative regions of the Mediterranean Sea depicted in Table 2 and Fig. 15. Figures 16, 17, 18 and 19 present the Q-Q plots and ECDF of the full-time series of hindcast versus raw and bias-adjusted CCLM4-MIROC5 and HIRHAM5-CNRM-CM5 significant wave height data for selected locations of the Mediterranean Sea. The results for RCA4-HadGEM2-ES are presented in the Supporting Information (Figs. SI-27 and SI-28). For CCLM4-MIROC5, it can be observed that for West-1 (Alboran Sea), Centre-1 (Tyrrhenian Sea), North-3 (Adriatic Sea) and South-1 (Gulf of Sidra) the raw CCLM4-MIROC5 data presents a good performance against hindcast, and the bias-adjusted simulations present a similar performance. On the other hand, for West-2 (Western Mediterranean), North-1 (Ligurian Sea), Centre-2 (Central Mediterranean), Centre-3 (Aegean Sea), and East-1 (Levantine Sea) the raw data presents a deviation from the hindcast, more noticeable in the higher percentiles that is then corrected by all EQM methods, with a slightly worse performance given by the EQM-dayofyear method. Then, the EQM method is an adequate method that allows the correct characterization of the extreme wave values in agreement with previous studies (Lemos et al. 2020b, a). These results are also largely reproduced for HIRHAM5-CNRM-CM5 (Figs. 18,  19), where the EQM-method is able to bias-correct the H s data in different locations. It is noted that for this model, the worse performance is given by the EQM-seasonal method for almost all locations. In order to understand the ability of the bias-adjusted methods to accurately capture the temporal variability of the different percentiles and for different GCM-RCM simulations, Figs. 20 and 21 present the seasonal Q-Q plots for CCLM4-MIROC5 and RCA4-HadGEM2-ES, respectively, for the locations West-2 (Western Mediterranean), Centre-2 (Central Mediterranean), Centre-3 (Aegean Sea) and East-1 (Levantine Sea).   The raw CCLM4-MIROC5 data presented, for the location West-2 (Western Mediterranean), a slight overestimation in the quantiles of significant wave height when considering the complete time series, and the bias-adjusted data showed a good performance for all methods. When comparing the seasonal quantiles (Fig. 20, 1st row) it can be observed that the raw data slightly overestimates the upper percentiles (99 and 99.9th) for Winter and Spring and gives a large overestimation of all quantiles during summer. Although all bias-adjustment methods are able to correct the Winter and Fall bias, for Spring and Summer, the EQM with the full-time series, is not able to accurately adjust the bias in the percentiles 99 and 99.9th for both seasons and also the lower percentiles during summer whereas the EQM-monthly and EQM-seasonal are able to accurately correct the biases. This behavior is consistent when considering a different GCM-RCM as observed in Fig. 21. When analyzing the day-of-year quantiles (Figs. 22, 23), it can be clearly observed that using an EQM that does not account for shorter temporal scales in waves (2nd row), could lead to incorrect characterizations of the upper percentiles during the periods with lower significant wave heights as can be observed for the two GCM-RCMs and for all points, where the summer quantiles are not correctly adjusted or present even higher biases than the raw data for the EQM-full method. The use of time-dependent bias-correction methods allows a better assessment of wave climate temporal variability for shorter timescales (month or day-of-year). It can be highlighted that, for the shortest timescale used in this work, the EQM-dayofyear method, although it provides the best non-stationary quantile characterization, it can lead to overfitting of the data. Therefore, the choice in bias-adjustment methods needs to take into account not only the different timescales present in wave climate but also the physical dependencies and correlations in waves.
For CCLM4-MIROC5, the locations Centre-2 (Central Mediterranean) and East-1 (Levantine Sea) provide a similar behavior with overestimations of the raw data for the upper percentiles when considering the full-time series and seasonal Q-Q plots. In these cases, the EQM-monthly presents the best performance, able to correct seasonal and shorterterm biases in all the considered percentiles. The need to account for wave temporal variability is clearly depicted when analyzing the results of the Centre-3 (Aegean Sea) location where the raw CCLM4-MIROC5 data presented an underestimation of all the percentile values when considering the full-time series Q-Q plot (Fig. 16). When considering the seasonal distribution, the raw data depicts an underestimation with respect to hindcast for Winter, Spring, and Fall and an overestimation during summer. Therefore, the bias-adjusted EQM data done with the full-time series, provides good performance for Winter and Spring, whereas for Summer and Fall, it provides overestimation with even higher biases than the raw data during Summer. Indeed, it can be observed for the non-stationary quantiles (Fig. 22, 3rd column) that the bias-adjusted EQM data presents higher errors than raw data for the higher quantiles during the months with lower significant wave heights due to the fact the bias-adjustment for the higher quantiles is dominated by the Winter values where an underestimation was given by the raw data and the overestimation during Summer is not considered in the adjustment. This behavior is less noticeable for the bias-adjusted EQM-seasonal data whereas the bias-adjusted EQM-monthly accurately captures the temporal variability of wave climate during the year.

Discussion and conclusions
This work presents a multi-model ensemble of wave climate projections in the Mediterranean Sea modeled with the wave propagation numerical model Wavewatch III forced by surface wind fields from 17 EURO-CORDEX GCM-RCMs. For the performance and bias-correction analysis of the wave projections, the simulations for the baseline period  were used providing data of the main wave North-2 (Adriatic Sea)  Fig. 15. The colors represent the hindcast and the different CCLM4-MIROC5 datasets: raw, bias-adjusted EQM, bias-adjusted EQM-seasonal, bias-adjusted EQM-monthly and bias-adjusted EQM-dayofyear parameters on a 3-h and 10 km, temporal and spatial resolution, respectively. This work focuses on the performance of GCM-RCMs against hindcast for significan t wave height and the use of different bias correction methods taking into account wave temporal variability. Under the assumption that the statistical properties of the biases of the different GCM-RCMs are maintained in the future, these biases can be corrected by applying the obtained bias correction model to future projections, which will be the aim of future work. In this study, the performance of the raw and bias-adjusted GCM-RCM data was analyzed by means of the RMSE and Bias metrics for the significant wave height monthly mean and maximum values and the gain/loss of performance of the bias-adjusted against raw data was assessed by means of the Skill Score and ΔBias . The results of the performance of the different GCM-RCMs during the baseline period show that all models present a similar spatial distribution of RMSE according to the distribution of significant wave heights, H s . Therefore, a larger RMSE are obtained in the Western and Central Mediterranean where larger H s are expected and lower biases are obtained in the Adriatic and Aegean Sea as well as near-coast regions. On the other hand, the Bias presents a large spatial variability for all GCM-RCMs although with overall lower biases ranging from − 0.2 to 0.2 m.
For the assessment of the role of temporal variability, the widespread Empirical Quantile Mapping bias-adjustment method was used under different temporal periods for the correction of significant wave height data against a validated hindcast for the baseline period 1979-2005. Although the EQM method is widely used and has proven to outperform other methods (Switanek et al. 2017), there are some limitations that should be considered. More specifically, when considering a large number of quantiles in the calibration process, there could be an over-fitting of the transfer function in the non-robust upper-tail of the distribution which, when applied to a different period with a different tail behavior could lead to undesirable results of the bias-adjusted tail distribution (Berg et al. 2022). To further explore this, we perform a split-sample or cross-validation test for the EQM-full method whose results are presented   Figures SI-31 and SI-32. Indeed it can be observed that due to the inherent noise in the upper tail of the distribution, the performance of the upper-tail reduces in the validation period when the calibration is performed in a different period. Therefore, for the calibration of the tail, it is advisable to decrease the number of quantiles when calibrating the distribution or fit a theoretical distribution valid for the upper-tail (Berg et al. 2022). In this work, we have not fitted a theoretical distribution to the tail and applied EQM to accurately compare the performance of the method when used for different temporal windows. Then, the EQM method is applied (1) using the complete time series of significant wave data, (2) seasonal data, (3) monthly data and (4) day-of-year data in order to analyze the ability of the bias-correcting methods to capture the different temporal scales present in wave climate ranging from storm events to monthly, seasonal and interannual variability.
The results show that the use of the EQM method for the full-time series without taking into account other timescales can lead to increased biases in some regions and seasons. The EQM bias correction provides generally better performance in the RMSE in the Mediterranean Sea than the timedependent methods during Winter where the use of a larger sample allows a reduction in the variance of the errors. On the other hand, the time-dependent bias-adjusted data provides better performance for the RMSE during the remaining seasons and the Bias for all seasons. This is clearly observed in the seasonal and non-stationary Q-Q plots where the traditional use of the EQM method for the complete time series leads to increased biases in comparison with the raw data in some locations.  Fig. 15. The colors represent the hindcast and the different HIR-HAM5-CNRM-CM5 datasets: raw, bias-adjusted EQM, bias-adjusted EQM-seasonal, bias-adjusted EQM-monthly and bias-adjusted EQMdayofyear It is concluded that the use of time-dependent bias-correction techniques leads to an improved accurate characterization of biases considering the interannual temporal variability of significant wave height. Nonetheless, the analysis should consider not only different temporal scales but the physical constraints within waves and the sample size for statistical significance. Then, the EQM-dayofyear, which is the shorter time scale considered in this work, although it provided an excellent characterization of non-stationary quantiles, led to slight deviations when analyzing the seasonal statistics due to the lack of correlation in day-of-year annual wave data. Then, the EQM-monthly was able to capture the temporal variability of waves and provide statistics compared to Hindcast data depicting |Bias| ≤ 0.1 m and RMSE ≤ 1.2 m.
This work addresses the need to account for temporal variability within bias-correction methods by applying the recognized and widely used EQM approach to different temporal scales crucial for wave analysis. Nonetheless, wave climate also presents large spatial correlations due to the link to synoptic and mesoscale atmospheric dynamics and a multivariate structure linked to the correlations between significant wave height, wave period, and wave direction. Then, further research should be done to include the intervariable and spatial correlations in waves.