Wind Waves in the Mediterranean Sea: An ERA5 Reanalysis Wind-Based Climatology

A climatology of the wind waves in the Mediterranean Sea is presented. The climate patterns, their spatio-temporal variability and change are based on a 40-year (1980–2019) wave hindcast, obtained by combining the ERA5 reanalysis wind forcing with the state-of-the-art WAVEWATCH III spectral wave model and verified against satellite altimetry. Results are presented for the typical (50th percentile) and extreme (99th percentile) significant wave height and, for the first time at the regional Mediterranean Sea scale, for the typical and extreme expected maximum individual wave height of sea states. The climate variability of wind waves is evaluated at seasonal scale by proposing and adopting a definition of seasons for the Mediterranean Sea states that is based on the satellite altimetry wave observations of stormy (winter) and calm (summer) months. The results, initially presented for the four seasons and then for winter and summer only, show the regions of the basin where largest waves occur and those with the largest temporal variability. A possible relationship with the atmospheric parameter anomalies and with teleconnection patterns (through climate indices) that motivates such variability is investigated, with results suggesting that the Scandinavian index variability is the most correlated to the Mediterranean Sea wind-wave variability, especially for typical winter sea states. Finally, a trend analysis shows that the Mediterranean Sea typical and extreme significant and maximum individual wave heights are decreasing during summer and increasing during winter.


INTRODUCTION
Wind waves are a key factor of the Earth global climate system contributing to the modulation of the exchanges at the atmosphere-ocean interface (Cavaleri et al., 2012b). At the same time, wind waves can significantly influence human activities at global-to-local scales, as they can grow under moderate-to-intense winds provided a sufficient fetch length is available, therefore impacting offshore structures and navigation as well as coastlines and coastal recreational or productive activities. The largest individual waves of sea states, particularly during storm events, are recognized as a potential risk of structural damage and of goods and even human-life losses at sea. Many accidents have been reported in the past and attributed to the occurrence of very large, sometimes abnormal, wave heights (Dysthe et al., 2008). Therefore, the assessment of the wind-wave climate (i.e., the wind-wave characteristics averaged on a long-term temporal scale) and of its spatio-temporal variability and change is of fundamental importance for coastal and offshore engineering purposes (to design, for instance, littoral protections, oil rigs and wind farms), navigation (ship routing) and for all other activities related to the marine environment (DNV GL-Det Norske Veritas Germanischer Lloyd, 2017). Beside the most widely used significant wave height, information on the maximum wave individual height expected at certain locations is more and more required by the naval and offshore industries for the definition of the environmental loads over the lifetime of a ship or a structure (DNV GL-Det Norske Veritas Germanischer Lloyd, 2017).
In this study, we characterize the long-term wind-wave climate of the Mediterranean Sea (hereinafter MS; Figure 1), located in southern Europe, which is a primary source of food, ecosystem services and economic activities for all the countries in the region. It is estimated that one-third of the Mediterranean countries' population lives along the MS coastlines (46,000km long), and that 30% of the global international tourism is headed to this region, particularly the coastal areas, where large and populated cities have developed (United Nations Environment Programme/Mediterranean Action Plan and Plan Bleu, 2020). The MS is also one of the world's busiest shipping lanes, carrying 20% of seaborne trade (including 17% of the world's oil tank capacity; United Nations Environment Programme/Mediterranean Action Plan and Plan Bleu, 2020), 10% of world container throughput and over 200 million passengers (Piante and Ody, 2015). During MS storms wind waves can grow until threatening heights that may cause naval accidents, as happened to the passenger ships "Voyager" in 2005 (Bertotti and Cavaleri, 2008) and "Louis Majesty" in 2010 (Cavaleri et al., 2012a), both cruising in the western MS. Besides this, the MS basin, with its processes and dynamics, is a pivotal environmental factor for the climate not only of the whole region (Lionello, 2012), but also of the neighboring areas and countries. Given its fragility and exposure to anthropogenic pressure, the MS is considered an early responder to climate change and a hot-spot of its effects and impacts (Giorgi, 2006). Hence, it is essential to assess both the past and present climatologies in order to create a solid reference for studies on the future climate.
The purpose of this paper is to characterize the typical and extreme climate of the MS wind waves, together with its spatio-temporal variability and changes at a long-term spanning four decades . To cope with the requirements of such a complex environment and characterize the wave climate therein, we rely on a long-term and high-resolution wave hindcast, obtained by combining numerical wave modeling and atmospheric reanalysis. We have thus run the WAVEWATCH III model [Tolman, 1991;The WAVEWATCH III Development Group [WW3DG], 2019] with the ERA5 reanalysis wind forcing (Hersbach et al., 2018(Hersbach et al., , 2020 and performed a wave climate analysis. In doing so, we do not expect to cover all aspects of the wave climatology of the MS, some of which have been addressed in various studies and publications (an non-exhaustive list of them can be found in Section "The study area"). We rather aim at bestowing new results based on state-of-the-art numerical modeling and reanalysis data, using one of the longest and most recent datasets to date. A novel season definition for the MS states based on observed data is proposed and adopted. Additionally, and most importantly, for the first time, we present and discuss the MS climate of expected maximum individual waves and of the most relevant factors for their development.
The paper is structured as follows. The study area, the wave hindcast and the methodology we have used for the wave climate analysis are detailed in Section "Methodology and Data". In Section "Results" we present the results of this wave climate analysis, which first focuses on the significant wave height of sea states, showing its typical and extreme seasonal climate in the MS, its spatio-temporal variability and change, and then on the expected maximum individual waves and relevant parameters for their estimation (wave steepness and narrow bandedness). Section "Discussion" is a discussion on the limits and strengths of this dataset and on the results of the wave climate analysis, also compared to previous studies in the same region. Conclusions in section "Conclusion" close the paper.

METHODOLOGY AND DATA The Study Area
The MS is a mostly deep (1,500 m on average) semi-enclosed regional sea, connected to the west to the Atlantic Ocean by the narrow Gibraltar Strait and to the east to the Black Sea by the Bosporus and Gallipoli (Dardanelles) Straits. Geographically, it is located between the European and the African continents and, from a climatic perspective, at the transition between a midlatitude regime zone, to the north, and a subtropical regime zone, to the south. Meteorologically, the MS region is interested by an intense synoptic-scale activity, with a separate branch of the northern hemisphere storm track passing over the region (Lionello and Sanna, 2005). The MS has a complex morphology, with many islands and peninsulas dividing the basin into subbasins, some of them connected by straits (Figure 1). Together with the presence of steep mountain ridges close to the coast, this morphological complexity provides to the MS environment with a spatial heterogeneity, which also reflects on the MS region meteorology and climate, making their successful modeling far from being trivial (Lionello, 2012). In addition, and most important for the topic of this study, the orographic complexity of the MS region plays a fundamental role in the genesis of cyclones. In this context, cyclones typically originate (mostly in the western MS) from the interaction of synoptic-scale flows with the mountains, constraining the flow in the lower troposphere and producing local and often strong winds (Lionello, 2012). FIGURE 1 | Map of the Mediterranean Sea basin and sub-basins (separated for the climate analysis purposes of this paper by dashed gray lines). The most relevant regional winds are indicated by black arrows denoting blow-to direction.
In general, the MS wind climate is characterized by the occurrence of winds with regional spatial structure (200-1,000 km) and seasonal variability (Zecchetto and De Biasio, 2007;Lionello, 2012). In Figure 1 we have indicated the most relevant winds for the different sub-basins. Levanter and Vendaval winds blow in the westernmost part of the MS (Western sub-basin), respectively, from the north-east and the south-west: the first during all seasons (but most intensely in winter), the second mainly during autumn and spring. Mistral is a cold dry wind jet, stronger during winter, blowing from the north/north-west over the western and Tyrrhenian sub-basins and occasionally reaching the African coasts. Libeccio and Sirocco are southerly winds blowing in the Tyrrhenian and Ionian-Meridional sub-basins: the first mainly in winter from the south-west, the second mainly in autumn and spring from the south-east. Sirocco is also channeled into the Adriatic sub-basin, and it contributes to the intense storm surges that take place in the northern Adriatic Sea and in the Venice lagoon . The other dominant wind in the Adriatic climate is Bora, a cold and dry north-easterly wind producing jets over the sea after interacting with the local orography. Northerly winds also dominate the climate in the Levantine-Aegean sub-basin: Etesian winds blow the from north/north-west and are particularly strong during summer.
The climate of wind waves in the MS, mostly driven by these regional winds, is subject to the morphological complexity constraints described above, resulting in limited fetches in large parts of the basin (with some exceptions) and in coastal shallow regions relatively proximal to the generation areas. Some regions present peculiar characteristics, for instance the northern Adriatic Sea or the Aegean Sea, where wave propagation experiences either the evolution over a relatively shallow flat bottom or the shadowing by numerous islands, respectively. Most of the studies on the MS wave climate (including a wind and wave atlas of the MS, Western European Union and Western European Armaments Organisation Research Cell, 2004) based on hindcast or observed data, focused on the extreme significant wave height (e.g., the monthly 99th percentile: Sartini et al., 2015a;Sartini et al., 2015bSartini et al., , 2017De Leo et al., 2020;Morales-Márquez et al., 2020) with only a few exceptions that focused on the average significant wave height (for instance, Lionello and Sanna, 2005;De Leo et al., 2020). Others were dedicated to the assessment of the wave energy climate at specific locations or over the whole basin, with energy harvesting purposes (Barbariol et al., 2013;Liberti et al., 2013;Arena et al., 2015;Besio et al., 2016). A similar approach to the one presented in this paper, although with a different purpose, can be found in von Schukmann et al. (2021), where Benetazzo et al. have used a shorter portion of our dataset for an assessment of the climate of the maximum individual waves in the MS in the last quarter of a century, focusing on the 2019 anomalies and their physical interpretation. The complexity of the MS climate has also been exploited to test novel approaches for wave climate analysis (Sartini et al., 2015b;Portilla-Yandún et al., 2019). Besides providing the statistics over decadal periods (generally one to three decades), some of the previous studies used statistical and data analysis approaches to obtain accurate return level estimates of the significant wave height (e.g., Sartini et al., 2015a), and to describe the spatio-temporal variability and the long-term change (trends) of the wave climate (e.g., Sartini et al., 2017;De Leo et al., 2020;Morales-Márquez et al., 2020). To interpret inter-annual and inter-decadal variabilities, some studies also tried to relate them to the modes of atmospheric variability, looking for teleconnection patterns and using climate indices (Lionello and Sanna, 2005;Morales-Márquez et al., 2020).
Finally, some studies were devoted to providing projections of the future wave climate under climate change scenarios over the whole basin or sub-basins (see for instance, Benetazzo et al., 2012;De Leo et al., 2021).

Wave Model Setup
The climate of MS wind waves, its spatio-temporal variability and change presented in this paper are grounded on a numerical model hindcast, produced by running the wave model WAVEWATCH III R (hereinafter WW3; version 6.07; 1 ) for a 40-year long period from 1980 to 2019 over the whole MS basin (Figure 1). WW3 is a third-generation spectral windwave model that solves the random-phase spectral action density balance equation for wavenumber-direction spectra [Tolman, 1991;The WAVEWATCH III Development Group [WW3DG], 2019]. Various source-term packages allow computing wave generation by wind, decay by dissipative processes, nonlinear wave-wave interactions, and wave transformation near the coast. For this study, WW3 has been run over a high-resolution structured curvilinear grid with a horizontal spacing of 0.05 • (∼5 km). WW3 has been forced with the 10-meter wind speed horizontal components (u 10 , v 10 ) provided hourly by the ERA5 atmospheric reanalysis, over a 0.25 • resolution horizontal structured grid (Hersbach et al., 2018). To represent wave growth and decay in WW3, we have used the state-of-theart ST4 source-term parameterization of Ardhuin et al. (2010), relying on the default values, but adjusting some coefficients in agreement with the results of TEST405 of the Ardhuin et al. study [The WAVEWATCH III Development Group [WW3DG], 2019] to values that are supposed to perform well for younger seas (β max = 1.55 and z 0,max = 0.002). Wave propagation has been computed using a third-order accurate scheme, together with the discrete interaction approximations for nonlinear wavewave interactions, and accounting for subgrid-scale obstructions (propagation-based approach; Tolman, 2003), very common in some parts of the MS (e.g., Adriatic and Aegean seas). In shallow waters (i.e., coastal areas and the northern Adriatic Sea) bottom friction has been modeled using the JONSWAP parameterization with default values (Hasselmann et al., 1973), while the depthinduced wave breaking has been parametrized following Battjes and Janssen (1978). To define the bottom topography and coastlines we have used the ETOPO1 relief model (doi: 10.7289/ V5C8276M).

Wave Model Outputs
The model outputs considered in this study, produced at hourly rate, are the significant wave height H s , the maximum expected crest height C max , and the maximum expected crestto-trough height H max . Results for H s and H max are presented in the paper, while results for C max are in the Supplementary Material. The maximum individual crest and crest-to-trough heights herein considered are computed using the space-time extreme WW3 model implementation distributed since version 5.16 [The WAVEWATCH III Development Group [WW3DG], 1 https://polar.ncep.noaa.gov/waves/ 2016] which assumes nonlinear (for C max ) and linear (for H max ) constructive interference of directionally spread waves as the leading mechanism for the generation of the largest waves in a sea state. In particular, the WW3 implementation is based on the Tayfun (1980) approximation of the second-order nonlinear model for C max and on the linear Quasi-Determinism model for H max , such as they can be defined over a given space-time domain as: (1) In Equations (1) and (2), ξ 0 is the mode of the probability density function of linear space-time extremes (Fedele, 2012), σ = H s /4 is the standard deviation of sea surface elevation, N 3 , N 2 and N 1 are the average numbers of waves in the spacetime domain and over its boundaries , γ = 0.5772 is the Euler-Mascheroni constant, µ is the wave steepness and ψ * is the narrow bandedness parameter. Steepness and bandwidth are sea state characteristics that contribute to the maximum individual wave height characterization. ψ * is defined as the absolute value of the first minimum of the autocovariance function of the sea state (Boccotti, 2000). As both µ and ψ * are not direct outputs of WW3, we have a posteriori derived them, as µ = H s /(gT 2 m /2π), with T m being the mean zero-crossing spectral wave period produced by WW3, and from H max = C max,1 √ 2(1 + ψ * ) (Boccotti, 2000) by deriving the linear estimate of the maximum expected crest height C max,1 from C max , using the Tayfun (1980) quadratic equation [Equation (36) of Benetazzo et al., 2017]. Further details on the spacetime extreme theoretical framework implemented in WW3 can be found in Barbariol et al. (2017); Benetazzo et al. (2021b) and references therein. A discussion on the role of wave steepness and bandwidth on the maximum individual wave heights during extreme events is given by Benetazzo et al. (2021a).
Our climatology of maximum individual waves deals with the maximum crest height C max and the maximum crest-totrough height H max that might be expected in 20-min sea states and over a 100 m × 100 m area. This duration and spatial extent point to a customary sea state duration [e.g., measured by in situ instruments as a buoy, World Meterological Organization [WMO], 2018] and to the horizontal sea surface area covered by an offshore structure (e.g., a large oil/drilling rig), respectively.

Assessment Methodology
The MS wind-wave hindcast has been assessed using satellite altimeter observations of H s , starting from 1991. Altimeter H s data have been retrieved from the IMOS platform 2 Young and Ribal, 2019), which provides calibrated and cross-validated data (we have only retained data flagged as "good data, " which according to the authors of those studies proved very reliable). Satellite missions considered in the August 1991-December 2019 period, for a total of 4066142 collocations, are listed in Supplementary Table 1. The longest dataset is 14-year long (ERS-2), the shortest is almost 4-year long (SENTINEL-3A). To inter-compare modeled and observed data for validation purposes, wave hindcast data have been bi-linearly interpolated in space and linearly interpolation in time on the satellite data points location and time. The statistical parameters and error metrics used to evaluate the hindcast performance are the model-observation bias BIAS = H s,WW3 −H s,altimeter (overbar denoting ensemble average), the relative bias RBIAS = (H s,WW3 /H s,altimeter −1) × 100, the Pearson cross-correlation coefficient CC, the Mean Absolute Error MAE = |H s,WW3 −H s,altimeter | and the slope of the best-fit line (linear regression).
To assess that the ERA5 climatological signal transferred from the surface wind to waves can represent the intra-and interannual climatological variability of MS waves, as well as the trends, we have estimated the correlation between the monthly 50th and 99th percentile H s from satellite altimeters and colocated WW3 hindcast and compared the linear trends of observed and modeled data. For these specific latter tasks, we have retained only satellite datasets longer than 10-year (see Supplementary Table 1).
Modeled maximum individual waves are generally verified at the short-term (i.e., sea state) time scale (see e.g., Barbariol et al., 2017Barbariol et al., , 2019Benetazzo et al., 2021a,b), as there are no observational platforms in the global oceans that are continuously collecting this type of data over a long-term period.
In addition, the model-to-observation comparison would require using observations with the same space and time domain sizes used for the hindcast (i.e., 20 min and 100 m × 100 m), while the sizes used for observations (application dependent and fixed) may be different from those required. Hence, the assessment of the WW3 hindcast in terms of C max and H max can be only done indirectly. We rely here on the assessment of H s , which is the pivotal parameter for determining the height of maximum expected individual crest and crest-to-trough heights (Benetazzo et al., 2021b).

Wave Climate Analysis
The climatology of wind-wave sea states in the MS is presented by using the 50th (i.e., median) and the 99th percentiles of the significant wave height H s , which may represent typical and extreme sea-state conditions, respectively. Similarly, the climatology of maximum individual waves is based on the 50th and 99th percentiles of the maximum crest-totrough height H max (and maximum crest height C max in the Supplementary Material).
In order to highlight the effects of intra-annual variability, the 50th and 99th percentiles are estimated at seasonal scale and averaged over the length of the hindcast to provide an empirical estimate of the expected typical and extreme conditions during the different seasons. In this study, we propose and adopt a novel definition of stormy seasons for the MS windwave climate based on the observed H s data in the MS, rather than using the canonical meteorological definition based on the air temperature (i.e., winter: December, January, February, DJF; spring: March, April, May, MAM; summer: June, July, August; JJA; autumn: September, October, November, SON). To this end, satellite altimeter observations of H s are used (see Section "Assessment Methodology"). Yearly-averaged monthly H s values are first computed for each satellite altimeter and then the eight satellites are combined to obtain the ensemble average and standard deviation.
To let the interannual climate variability and climate change due to external forcings emerge, we remove the source of intrinsic intra-annual variability due to seasonality by taking separately into account the different seasons, as done for instance by Morales-Márquez et al. (2020). Therefore, in the result presentation (climate variability and climate change), we shall focus on the two most representative seasons only, i.e., winter and summer. The interannual climate variability is shown in terms of climate anomalies, which are defined as deviations of the winter and summer 50th and 99th percentile values of a generic wave height H in a specific i th year H pp i (pp = 50, 99) with respect to the winter and summer yearly-averaged 50th and 99th percentiles H pp (overbar denoting average over the years), i.e., H pp i − H pp . Further, to assess the spatial distribution of the inter-annual climate variability and highlight the MS subbasins and regions that have experienced the largest variations in the study period, we have also normalized anomalies with respect to the local yearly averages, thus showing the relative anomalies, i.e., (H pp i − H pp )/H pp . In order to motivate the inter-annual climate variability of the MS waves, we look for possible teleconnections between the wave climate variability and the principal modes of atmospheric variability, taking into account four northern hemisphere indices (downloaded as monthly averaged values from 3 ): North Atlantic Oscillation (NAO), Scandinavian (SCAND), East Atlantic (EA) and East Atlantic-Western Russian (EAWR). Monthly mean data are then combined to obtain seasonal averages. Finally, in order to assess any potential long-term change in the typical and extreme sea state climate over 1981-2019, the linear trend of winter and summer H s , C max and H max are estimated over the whole MS and sub-basins through the Sen's slope of the 50th and 99th percentiles (Theil, 1950;Sen, 1968), testing the statistical significance at the 90% confidence interval with the Mann-Kendall test (Mann, 1945;Kendall, 1948).

RESULTS
The characteristics of the MS wind-wave climate obtained from the ERA5 wind forced WW3 hindcast are presented in this Section stemming from the yearly-averaged seasonal 50th and 99th percentiles of H s (the conventional variable used to describe the wave climate), together with their variability and change. Same analyses are also conducted to investigate the distribution of the crest-to-trough height of the maximum individual waves H max (results for the crest height of the maximum individual waves C max are in the Supplementary Material). Before presenting the MS wind-wave climatology a novel definition of seasons for the MS states and the assessment of the wave hindcast performance are introduced.

An Observed Data-Based Definition of Seasons for the Mediterranean Sea States
The monthly 50th and 99th percentiles of H s , as observed over the MS by the eight satellite altimeters used for the hindcast assessment (Supplementary Table 1) are shown in Figure 2A.
Observations are also normalized with respect to the maximum monthly H s percentile ( Figure 2B), which occurs in January for both the 50th and 99th percentiles. Both typical and extreme sea states clearly show two dominant regimes: a stormy period, from November to March (with values within 20% of the January H s ; Figure 2B), and a relatively calm period, from June to August (with values within few per cent of the minimum H s , occurring in June for the 50th percentile and in August for the 99th percentile). In between, there are two transitional seasons. Following this variability, we define the four sea-state "seasons" as winter (NDJFM), spring (AM), summer (JJA) and autumn (SO). Compared to the classic meteorological definition, we therefore propose a longer winter season, starting in November and ending in March (both showing significantly larger H s compared to the other spring and autumn months), and shorter transitional seasons (2-month long), with summer remaining unchanged. This definition applies to the whole basin but still holds if only northern, southern or western sub-basins are considered. If only the eastern part of MS is accounted for, the relative difference between summer and transitional seasons reduces (not shown here), most likely due to the effect of Etesian wind-waves. This is in agreement with findings of Queffeulou and Bentamy (2007) who pointed out differences in western to eastern MS monthly mean H s from satellite altimetry, especially during summer season. November and December H s values clearly belong to the winter season that accounts for the early months of the following year; therefore, to have an equal number of seasons, we start counting seasons from winter 1981 (beginning with November 1980) to autumn 2019 (ending in October 2019), thus analyzing only 39 years in the whole period (hereinafter 1981-2019).

Assessment of the Mediterranean Ses Wind-Wave Hindcast
The MS wind-wave hindcast performance in the 1991-2019 period against satellite altimeter observations is summarized in Figure 3, which collects information about the relative model-observation bias RBIAS in the MS and in the subbasins (Figure 3A), the H s scatter ( Figure 3B) and the relation between modeled and observed climate variability and change (Figures 3C,D). Modeled H s is simulated in the MS with 0.87 best-fit line slope, BIAS = −0.12 m and MAE = 0.26 m. The relative model-observation bias at the 50th percentile of modeled H s (0.93 m; blue dot in the Figure 3B) is −8%, while at the 99th percentile (3.78 m; red dot in the Figure 3B) it is −14%, denoting a decrease in model performance with H s , in agreement with the MS wave hindcast study of Mentaschi et al. (2015; verified against in situ wave buoys). The model-observed bias might be partially corrected by calibrating the input/dissipation source term parameterization in WW3. However, the spatial heterogeneity of the MS (with narrow, semi-enclosed sub-basins or/and surrounding orography) makes it difficult to achieve a set of parameters that is optimized to obtain good performance in FIGURE 2 | Yearly-averaged monthly 50th (black dots) and 99th (red diamonds) percentiles of H s over the MS, from eight satellite altimeters (see Supplementary  Table 1  the whole MS basin and in the whole H s range (Lionello and Sanna, 2005;Mentaschi et al., 2015;Sartini et al., 2017;Morales-Márquez et al., 2020). This can be seen by comparing the slope of best-fit lines and the relative bias at different percentiles for the MS with those for the sub-basins. Best performance is achieved in the widest sub-basins (Ionian, Levantine), with larger slope and smaller relative bias compared to those obtained in the whole MS. Worst performance are in the Adriatic Sea, where the ERA5 wind is less effective in correctly reproducing wind and then waves due to spatial resolution effects in a narrow basin. These effects are partially visible also in the Tyrrhenian and Western sub-basins (both with islands and narrow parts), although performance in these sub-basins is comparable to that for the MS.
The average (and standard deviation of) correlation coefficients between the monthly 50th and 99th percentile H s from satellite altimeters and co-located WW3 hindcast are: 0.98 ± 0.004 for the 50th percentile and 0.97 ± 0.009 for the 99th percentile. An example of the monthly percentiles time-series for JASON-1 is shown in Figures 3C,D, with correlation coefficients larger than or equal to 0.98 for both 50th and 99th percentiles. Also, the linear trends of H s estimated from the same satellite altimeters and co-located model data are in agreement, with consistent (i.e., both positive or both negative) Sen's slope signs and a general tendency for larger model trends for the 50th percentiles and smaller for the 99th percentiles (compared to JASON-1, the relative trend differences are 15 and −16% for the 50th and 99th percentile, respectively, while the largest relative difference is −43%, between Envisat and model 99th percentile H s trends).

Mediterranean Sea Wind-Wave Climate
The yearly-averaged seasonal 50th percentile of H s is shown in Figure 4. As expected, the largest values occur during the winter (NDJFM) season when wind speeds are typically stronger (Zecchetto and De Biasio, 2007). Spatial patterns of the 50th percentile H s are rather different in winter and summer (JJA), on the one hand, and very similar in spring (AM) and autumn (SO), on the other hand. The most energetic sea states occur over the western and southern sub-basins in winter (50th percentile up to 1.5 m) and spring (up to about 1 m), and over the eastern and western sub-basins in summer (up to 1.25 m) and autumn (up to 1.0 m). The western basin experiences the largest typical conditions over all seasons except summer, particularly the Sardinian Sea (amid the Gulf of Lyon, Sardinia, the North African coast and the Balearic Islands) where the north-westerly Mistral wind persistently blows from the French inland (Zecchetto and De Biasio, 2007). In summer, the largest wave heights are attained in the Aegean and Levantine basins (particularly southeast of Crete, with a typical jet-like pattern) with values comparable to those of winter and produced by the northwesterly Etesian winds, typical of late spring, summer and autumn seasons but reaching their maximum intensity during summer (Zecchetto and De Biasio, 2007).
The yearly-averaged seasonal 99th percentile of H s is shown in Figure 5. The spatial distribution of extreme sea states is different compared to that of typical conditions (Figure 4). Indeed, the 99th percentile of H s in all seasons have similar spatial patterns, and in particular the sea severities during spring and autumn seem comparable (despite some minor differences, for instance in the Alboran Sea and the Levantine basin). The largest values occur in the Gulf of Lyon and Sardinian Sea (up to 5.5 m in winter, up to 4.0 m in summer, up to 4.5 m in spring and autumn). Other regions with large, though smaller, values are located in the southern (Ionian Sea) and eastern (Levantine Sea) sub-basins, while the summer waves southeast of Crete island (the largest in the summer 50th percentile maps; jet-like pattern) are surpassed by the extreme waves generated by the Mistral wind in the Gulf of Lyon. In the Adriatic Sea, the spatial  pattern changes compared to the 50th percentile, with large values more homogeneously distributed over the whole basin (including in the northern part) and with the largest ones (in the southern part) off the eastern coast, as a consequence of the south-easterly Sirocco wind storms hitting the Croatia, Albany and Montenegro shorelines.

Variability of the Mediterranean Sea Wind-Wave Climate
The first assessment of the wave climate variability is made by analyzing the intra-seasonal extreme-to-typical condition variations. In Figure 6 we show the winter (NDJFM) and summer (JJA) ratios 99 p,H s /50 p,H s of the 99th to the 50th percentile of H s . This ratio provides a measure of the width of the empirical probability distribution of H s at seasonal scale and it is yearlyaveraged to obtain an empirical estimate of the expectation. We observe that sea states are much more variable during summer than winter, with the largest values of 99 p,H s /50 p,H s . There is also a larger spatial variability of the ratio during summer with respect to winter. Indeed, during winter extreme conditions are generally up to four times greater than the typical ones, over all the MS basin, with only spatially limited exceptions in the northern Adriatic, southern Tyrrhenian, western Aegean seas and eastern Levantine sub-basin, where they can be up to six times larger (locally, more than six times). During summer, extreme conditions can be more than six times larger than typical ones (up to eight times) over large parts of the central MS, including the Sardinian, Tyrrhenian and Adriatic Seas. The lowest summer variability occurs in the Levantine basin, where the sea states generated by the steady and low-variability Etesian winds are intense but constant over the season.
We now proceed by assessing the interannual variability of the wave climate in the MS region by showing the standard deviation and range of the anomalies of typical and extreme H s during winter and summer. We first focus on the standard deviation of the H s anomaly, which is used to point out the regions with the largest climate variability over the years. Then, we focus on the standard deviation and range (maximum-minusminimum) of the relative anomaly (i.e., with respect to the yearly-averaged percentiles of H s , see Section "Wave Climate Analysis" for the definition) to show the regions with the largest variability with respect to the local wind-wave climate (and quantify it), allowing to inter-compare regions with different H s percentiles. The interannual variability of the 50th percentile of H s is presented in Figure 7. During winter the largest variability occurs in the Western sub-basin (Mistral wind region, up to 0.25 m standard deviation; Figure 7A); in summer, in the Aegean and Levantine sub-basins (Etesian wind region, up to 0.2 m standard deviation; Figure 7D). These are also the regions with the largest 50th percentile of H s (Figure 4) although the spatial distribution of the standard deviation of the H s anomaly does not perfectly reflect that of the yearly-averaged 50th percentile of H s . Looking at the standard deviation of the relative anomalies (Figures 7B,E), it emerges that the regions with the largest variability with respect to the local wave climate are in the Tyrrhenian, Adriatic and Aegean sub-basins during summer, and in the Western, Tyrrhenian, Adriatic, and Levantine subbasins during winter. The maximum standard deviation of the relative H s anomaly, up to 20% of the 50th percentile H s in both summer and winter, occurs in the northern Adriatic Sea and the southern Tyrrhenian Sea, where during 1981-2019 winters the range of variability has reached 100% of the local 50th percentile H s .
The interannual variability of the 99th percentile of H s is presented in Figure 7, and show different characteristics with respect to that of the 50th percentile of H s . Indeed, in summer the maximum standard deviation of the H s anomaly is in the Western sub-basin (Sardinian Sea, 0.6 m; Figure 7J), where also the maximum 99th percentile occurs, whereas in winter the maximum standard deviation is in the southernmost (Ionian sub-basin, offshore Libya; Figure 7G) and westernmost (Western sub-basin, southeast and southwest of Balearic Islands) parts of the MS, which are not the locations with the maximum 99th percentile of H s (Sardinian Sea and Gulf of Lyon, Figure 5). Similarly, the largest standard deviation of the relative anomaly (20% of the 99th percentile H s ) is in the Western, Tyrrhenian, Adriatic, Aegean and Ionian sub-basins during summer (as for the 50th percentile H s ; Figure 7K), while in winter the maximum is offshore Libya and southwest of the Balearic Islands ( Figure 7H). In these regions, the range of the relative anomaly is up to the 100% of the local 99th percentile H s .
To provide a synthetic description of the interannual climate variability at the regional scale (i.e., MS basin and sub-basins), in Figure 8 we show the time series of the spatially averaged relative anomaly of the 50th and 99th percentile of H s over the MS (spatial averages over the sub-basins are shown in Supplementary  Figures 1-5). In general, during 1981-2019 there is an alternation of positive and negative anomalies that, more frequently during winter and especially in some sub-basins, tend to aggregate forming clusters of years with (almost) only positive or negative anomalies. An example of this aggregation is the 1989-1994 anomaly of the winter 50th percentile H s . During this period the largest negative anomaly in the MS have occurred (winter 1990, −18%), followed by four more seasons with negative, often considerable, anomalies. The largest negative anomaly at sub-basin scale have occurred in the Tyrrhenian Sea in the same period (winter 1989, −36% of the 50th percentile of H s ; Supplementary Figure 1), but similar pattern can be observed in the Ionian, Adriatic and Western sub-basins, even if with smaller anomalies. During 1989-1994 winters, the atmosphere in the MS region has been very stable, in particular over the Western, Tyrrhenian, Adriatic and Ionian sub-basins ( Supplementary  Figures 1-4; as shown in Supplementary Figure 6, displaying the 500 hPa geopotential height, sea level pressure and wind intensity anomalies over western Atlantic, Europe and western Asia from ERA5 atmospheric reanalysis, for some seasons of interest). The cause of this stability can be found in the localization of the pressure fields on a synoptic (European) scale (Supplementary   summer and winter (for the Tyrrhenian Sea, 11.3 and 13.6% standard deviation, respectively; see Supplementary Figure 1), in agreement with results about typical climate variability in Figure 7. Extreme H s appears to be most variable in the Tyrrhenian and Adriatic seas during summer, and in the Ionian and Levantine sub-basins during winter, (for the Adriatic Sea summer 11.2% and for the Ionian Sea winter 8.9% standard deviation, respectively; see Supplementary Figures 2, 4), in agreement with results for the climate variability shown in Figure 7.
The inter-annual climate variability of the MS waves observed in Figure 8 (and Supplementary Figures 1-5 for sub-basins) can be partially explained by possible teleconnections between the wave signals and the principal modes of atmospheric variability. To this end, we compare the time series of the relative anomaly with the time series of the NAO, SCAND, EA and EAWR indices. Thus, we have (i) averaged the monthly values of 1981-2019 indices over the winter and summer seasons, respectively, (ii) computed cross-correlation coefficients between the four indices and the MS (and sub-basins) H s relative anomalies and then (iii) plotted in Figure 8 (and Supplementary Figures 1-5 for sub-basins) the time series of the index showing the highest correlation (in search for teleconnection). Generally, with only a few exceptions, relative anomalies in the MS and sub-basins are mostly correlated to SCAND positive phases, particularly, but not only, for the 50th percentiles. In the MS, this is observed for winter and summer typical and extreme sea states, except for winter extreme H s , which is related to the EA negative phases. The highest correlations are found for the winter 50th percentiles, although the correlation coefficients are generally small, reaching 0.54 at most for the MS (SCAND and winter 50th percentile, Figure 8) and 0.65 maximum value for the sub-basins (SCAND and Tyrrhenian Sea winter 50th percentile, Supplementary Figure 1). The impact of the SCAND index on the Mediterranean atmosphere dynamics and hence on the MS wave climate can be observed, as an example, in the winter 1989 and 1990 atmospheric configuration (Supplementary Figure 6). The persistent low 500 hPa geopotential height and SLP fields over Scandinavia and Eastern Europe (with −1.2 and −0.8 seasonal SCAND index, respectively; Figure 8), induce negative anomalies of wind speed over the whole MS, which, in turn, generate negative anomalies of the 50th percentile H s (−15 and −18%, Figure 8).

Change of the Mediterranean Sea Wind-Wave Climate
We assess here what might have been the climate change of the winter and summer MS sea states over 1981-2019 showing    Averages over the points where the trends are significant (black dots in Figure 9; tested with Mann-Kendall test and 90% confidence interval) and over all the points. Value between parentheses is the percentage of points with statistically significant trends in a (sub-) basin with respect to all the points in the (sub-) basin.
the maps of the H s linear trends (using the Sen's slope, in cm decade −1 , evaluated on 39 yearly data) in Figure 9 and their spatial averages over the basin and sub-basins in Table 1. We show both the trends at the points where there is statistically significance (at 90% confidence interval according to the Mann-Kendall test) and at all the points in a (sub-) basin, including where there is no statistical significance. The winter 50th percentile of H s (Figure 9A and Table 1) shows positive trends over large part of the western and northern sub-basins of the MS, and negative or no trends in the eastern and southern sub-basins, resulting in a substantially net zero trend for the MS typical winter sea states (0.3 cm decade −1 ; Table 1). The largest increases are in the Sardinian Sea and offshore Algerian coasts (up to 5 cm decade −1 ) and largest decreases in the southern MS (down to −3 cm decade −1 ). However, only trends offshore the Algerian coasts, in the central Adriatic Sea (Sirocco and middle Bora jets regions) and in small portions of the Aegean and Balearic seas are statistically significant, bringing the spatial average to 1.1 cm decade −1 ( Table 1). The summer map ( Figure 9B) is far different, showing generally smaller Sen's slopes, except in some limited areas in the eastern part of the basin, with statistically significant increases (up to 4 cm decade −1 ) in the Aegean Sea and eastern and western Levantine basin and decreases (down to −4 cm decade −1 ) in the central Levantine basin where the summer Etesian winds blow. As a whole, the MS presents a very mild decreasing trend during typical summer sea states (−0.8 cm decade −1 ; Table 1).
The 99th percentile H s trends, shown in Figures 9C,D, are significantly larger than the 50th percentile trends shown above. During winter the largest statistically significant changes are in the southern Levantine sub-basin with increases up to 12 cm decade −1 offshore Egypt and Israel (11.0 cm decade −1 on average in the sub-basin), in the Gulf of Lyon with increases up to 12 cm decade −1 (Mistral jet region, despite a smaller average trend in the Western sub-basin), in the northern Adriatic Sea (northern Bora jet, increases up to 8 cm decade −1 ) and in the southern Tyrrhenian Sea, the only statistically significant area of negative trend for winter extreme sea states (down to −12 cm decade −1 , −10.0 cm decade −1 on average in the Tyrrhenian Sea). In the MS the winter 99th percentile H s has increased 6.5 cm decade −1 (Table 1) as a whole. In summer, the trend distribution and sign are quite different from the winter ones. Indeed, globally the 99th percentile H s has decreased (−1.6 cm decade −1 ; Table 1), and the spatial distribution is quite similar to the distribution of the summer 50th percentile H s trends. In particular, a statistically significant increase is observed in the Ionian (up to 12 cm decade −1 , 7.9 cm decade −1 on average in the sub-basin) and Aegean seas only. The Levantine basin shows a decrease down to −10 cm decade −1 and of −6.6 cm decade −1 on average, while the Balearic, Tyrrhenian (including the Gulf of Genoa) and Adriatic seas show moderate decreasing trends in the summer extreme sea states. Concluding, on average, in 1981-2019 the MS wind-wave climate has decreased during summer seasons and increased during winter seasons (both 50th and 99th percentiles) with the largest changes in the extreme sea states.

Mediterranean Sea Maximum Individual Wave Climate
In this section we change the focus of the analysis from the climate of the significant wave height of sea states to the MS climate of the maximum individual wave heights and the other (beside significant wave height) relevant parameter for their estimation, i.e., wave steepness and narrow bandedness.

Maximum Individual Waves
The yearly-averaged winter and summer 50th and 99th percentiles of H max are shown in Figure 10 (and in Supplementary Figure 7 for C max ). As expected, seasonal and spatial patterns closely mirror those of H s , which is the driving parameter determining the height of maximum individual waves (others, like the steepness for C max and the narrow bandedness parameter for H max , will be discussed in the following). The largest typical crest heights occur in the Sardinian Sea (Western sub-basin) and in the southern MS during winter, reaching up to 3.5 m height, whereas during summer the largest values are found southeast of Crete island in the Levantine basin with heights up to 3 m. The largest extreme wave height (99th percentile) can be found in the same parts of the MS with values up to 12 m in winter and 8 m in summer.
The variability of the maximum individual wave climate (C max and H max ) well matches the variability of the H s climate, shown in section "Variability of the Mediterranean Sea Wind-Wave Climate." The extreme-to-typical condition variations at seasonal scale for H max and C max are shown in Supplementary  Figures 8, 9, respectively. Compared to the H s variations in Figure 6, they only show minor differences, which do not substantially change the conclusions drawn about the seasonal and regional characteristics of the wave climate variability in the MS, which are also valid for the inter-annual variability.
What deserves more attention is the change in the maximum individual wave climate, expressed using the Sen's slope of the linear trend over 1981-2019. Trends are shown in detail for H max in Figure 11 (in Supplementary Figure 10 for C max ) and summarized as spatial averages over the MS and sub-basins in Table 2 (in Supplementary Table 2 for C max ). On average, maximum individual waves change more than H s ( Table 1) with slopes larger than H s slopes: e.g., on average for the MS in winter 2.5 and 13.9 cm decade −1 for H max , instead of 1.1 and 6.5 cm decade −1 for corresponding H s values. More than doubled slopes for H max are also observed in Figure 11, with winter H max increasing up to 10 cm decade −1 in the Western sub-basin (50th percentile) and up to 25 cm decade −1 in the Levantine sub-basin (99th percentile). In general, however, the spatial distribution of the trends for H max (and C max too) reflects the distribution of H s trends. The meaning of these trends and how they can be interpreted compared to the H s trends will be examined in the Discussion section.

Wave Steepness and Narrow Bandedness
We have shown that, while the maximum individual wave climate displays characteristic typical and extreme values, spatial patterns and variability of the maximum individual wave climate closely mirror those of the significant wave height. The reason is that  (1) and (2)]. Therefore, we expect that regions where H s is large are also the regions where C max and H max are the highest. This is confirmed by the comparison of Figure 10 and Supplementary Figure 7. At the same time, though with a smaller effect with respect to H s , C max and H max are also monotonic functions of the wave steepness µ and the narrow bandedness ψ * , respectively [see Equations (1) and (2)]. Therefore, maximum individual wave heights in Figure 10 and Supplementary Figure 7 are the result of the combined effect of H s (both C max and H max ), µ (C max ) and ψ * (H max ).
In Figures 12A,B we show the yearly-averaged winter and summer 50th percentile of µ. As younger waves are generally steeper, largest steepness (up to 0.07) can be found close to the shores where prevailing Mediterranean onshore winds blow from: e.g., in the Gulf of Lyon (Mistral), in the Aegean Sea (Etesian), along the eastern Adriatic coast (Bora), but also on the leeward side of straits (i.e., Sirocco to the north of Otranto Strait in the Adriatic Sea, Vendaval to the east of Gibraltar Strait, Mistral to the east of Bonifacio Strait in the Tyrrhenian Sea, Etesian to the west of the Bosporus Strait). Steepness decreases moving off-shore along the wind direction as waves develop. Although the steepest waves occur in winter, the summer Etesian wind waves are as steep as in winter. The 99th percentile of µ (Figures 12C,D), again larger in correspondence of younger waves, is up to 0.075 over the whole MS except in the Adriatic Sea (up to 0.08) and in very narrow coastal regions where can be up to 0.09. Wind-sea dominated seas typically have steepness values larger than 0.03-0.04, while swell dominated seas may have halved values compared to wind-seas .
The steepness values we find in the MS (generally larger than 0.03) denotes the predominance of wind-seas in the MS climate, especially in winter. During summer smaller values (0.02; 50th percentile) can be found in the Gulf of Genoa and in the southern Tyrrhenian Sea, where mature waves driven by south-westerly winds can develop after frequency-and direction-dispersion along the large distances from the Alboran Sea to these areas (among the longest in the MS).
The yearly-averaged 50th percentile of ψ * , shown in Figures 13A,B, generally ranges between 0.625 and 0.675 in winter, while can be up to 0.75 in summer. These are values representative of a unimodal wind-sea (0.65 ≤ ψ * ≤ 0.75; see Boccotti, 2000), denoting the prevailing nature of typical MS sea states. Smaller values can be found in the wide western and southern basins where a local wind-sea may combine with more mature seas coming from different directions. It is worth mentioning that a proper long-term statistical analysis of the sea state characteristics should account for different wave systems separately. That is, wind-sea and swell (if any) should be separated (e.g., using spectral partitioning techniques accounting also for wave direction) and treated independently as they have a different origin and belong statistically to different populations. However, the range of superimposed wind-sea and swell (ψ * ≤ 0.6; Boccotti, 2000) is rather rare in the MS, even if crossing sea states with waves with different wave direction may occur. The 99th percentile values of ψ * (Figures 13C,D), generally larger than 0.7 and up to 1, denotes the presence of narrow banded sea states, able to produce the largest crest-to-trough excursions according to the Quasi-Determinism theory (Boccotti, 2000;Benetazzo et al., 2017Benetazzo et al., , 2021b, especially in the Adriatic, Tyrrhenian and Aegean seas.

DISCUSSION
In this section, we first highlight and discuss the limits and strengths of the dataset and of the related climate analysis. Then, TABLE 2 | Spatial averages over the MS and sub-basins of the 50th and 99th percentile H max trends shown in Figure 11 (cm decade −1 ), during 1981-2019 summer and winter seasons. Averages over the points where the trends are significant (black dots in Figure 11; tested with Mann-Kendall test and 90% confidence interval) and over all the points. Value between parentheses is the percentage of points with statistically significant trends in a (sub-) basin with respect to all the points in the (sub-) basin.
we comment on the results we have obtained, by also comparing them to those obtained in previous studies for the same area. The climatology of MS waves presented in the previous section relies upon the climatology of the ERA5 reanalysis wind. Even if ERA5 horizontal resolution may be poor for representing the dynamical structure of the wind fields in some regional seas, as the narrow semi-enclosed Adriatic or Aegean seas, it currently represents the state-of-the-art of long-term atmospheric reanalysis and hence the best forcing option for the present work. In this context, using a dynamical downscaling of ERA5 winds, Vannucchi et al. (2021) showed an improvement of wind and wave hindcast performance against wind stations and buoys in the MS compared to the original ERA5 data mostly at coastal stations. We have decided to use ERA5 wind and not ERA5 waves, as the latter may be not suitable for long-term climate analyses, at least for studies extending before and after 1992, i.e., the beginning of satellite altimeter data assimilation in the wave models. The introduction of spurious trends is a well-known issue pointed out by Aarnes et al. (2015) for ERA-Interim reanalysis waves (Dee et al., 2011) that may apply also to ERA5. To cope with the resolution issue of ERA5 in the MS and with the potential 1992 singularity in wave products, we have therefore (i) relied on the ERA5 wind forcing (on which we expect a smaller or at least indirect data assimilation effect) and on an ad hoc high-resolution (about 5 km) wave hindcast and (ii) verified the wave hindcast results against satellite altimetry. While this has not increased the wind forcing resolution, it has allowed us to produce one of the longest consistent Mediterranean wave hindcasts up to now, providing a detailed representation of the wave climatology in a morphologically complex enclosed sea and assessing its uncertainty compared to observations in the region. Concerning the hindcast performance, the MS wave climate, its (intra-and inter-annual) variability and the trends have shown to be well reproduced by the ERA5 wind forced hindcast. However, the negative modeled-observed H s relative bias shown in Section "Assessment of the Mediterranean Ses Wind-Wave Hindcast" should be regarded when using the estimates of the 50th and 99th percentile of the MS wave heights.
Results of the wave hindcast and following wave climate analysis have shown marked regional characteristics of the wave climate, related to the peculiarities of the regional winds blowing over the Mediterranean sub-basins. Results for the typical and extreme wind-wave climate are in agreement with other studies both for spatial distribution and intensity, albeit with some differences [see e.g., Morales-Márquez et al. (2020) for extreme winter H s and Lionello and Sanna (2005) for mean seasonal H s , both presenting slightly smaller maxima]. In particular, results have highlighted the MS regions where the largest waves occur, which are in the Western, Ionian and Levantine sub-basins, depending on the season. Indeed, compared to previous studies, in this paper we have provided results at a seasonal scale. We have shown that the variety of the MS wave climate characteristics may be well represented by winter and summer seasons, in terms of intensity and patterns of wave heights, with spring and autumn in-between. This is in agreement with the results of Lionello and Sanna (2005), who argued that the MS climate, including wave climate, is characterized by two main seasons, i.e., winter and summer, with spring and autumn only being transitional seasons.
The definition of seasons we have adopted in this study is based on the wind-wave climate and on the altimeter H s observations in particular, rather than on monthly air temperatures, as in the meteorological definition of seasons. This choice is supported by other studies on the MS region climate. Indeed, Trigo et al. (1999) suggested that the traditional four meteorological seasons do not fit well the cyclone occurrence patterns in the MS region. In this context, attempts to objectively define season start/end and duration in the MS region, based on different indices than the largely spatially variable temperature, led to definitions of seasons that generally correspond to the traditional ones, but with longer winters and shorter springs and autumns (Alpert et al., 2004;Kotsias et al., 2021) in agreement with our definition. Also Lionello and Sanna (2005), based on principal component analysis of the modeled H s , defined four seasons, again with longer winters (and summers) and shorter springs and autumns.
Following the conclusions of previous studies that found that the climate (including wave climate) variability of part of the MS (in particular the northern part with a mid-latitude climate regime) can be partially related to midlatitude climate indices (e.g., Lionello and Sanna, 2005;Morales-Márquez et al., 2020)  we have searched for teleconnections between wave climate variability and the principal modes of atmospheric variability. We have found that MS typical (50th percentile) wave climate variability is generally positively correlated to the SCAND index (i.e., positive H s anomalies generally corresponding to positive SCAND phases, and viceversa). This is in accordance with the SCAND positive phases being generally associated with enhanced cyclogenesis in the MS region (Xoplaki, 2002). Extreme wave variability (99th percentile) is also positively correlated to the SCAND index in summer, while it can be negatively correlated to the EA index in winter. The negative phases of EA were also related to extreme waves in the MS by Izaguirre et al. (2010) using satellite data. Over the sub-basins there can be differences, for example the Adriatic Sea, whose wave climate variability is always positively correlated to the SCAND index during winter and to EAWR during summer. In general, our results are in agreement with the results from another study (Morales-Márquez et al., 2020) and suggest that the MS wave climate variability is mildly-to-slightly related to the variability of large-scale atmospheric structures, which have a more clear influence on the wave climate in other seas (e.g., the North-Atlantic Ocean; Morales-Márquez et al., 2020). We thus share the conclusions of Lionello and Galati (2008) that even if some patterns seem to exert a larger influence with respect to others, this is not sufficient to explain the dominant variability of the wave climate in the MS and its sub-basins, which is the result of a combined effect of large-scale atmosphere dynamics and orography (forcing regional winds).
The observed trend of decreasing H s during summer and increasing H s during winter (both for typical and extreme wave climate; Table 1) seems to disagree with some other past climate change studies in the region, although based on different wind forcing datasets, wave models, trend estimate, and season definition. Indeed, Lionello and Sanna (2005) found a statistically significant decreasing trend of the mean H s in the MS during winter months (December to March, −2.4 cm decade −1 ), while we have found a statistically significant increasing trend of the median H s in the MS during November to March (Table 1, 1.1 cm decade −1 ). Also, the spatial distribution and intensity of the 99th percentile winter (December to February) H s trend obtained by Morales-Márquez et al. (2020) only partially resembles that shown in Figure 9C. There is instead remarkable agreement, at least in some specific regions, with the trend estimates for the annual maximum H s by De Leo et al. (2020). Comparing those estimates with ours for the winter 99th percentile H s (the closest to the annual maxima), we see many areas of corresponding increase/decrease: for instance, positive trends in the Gulf of Lyon and other parts of the Western basin, in the northern Adriatic and Tyrrhenian seas and in the easternmost part of the Levantine basin; negative trends in the southern Tyrrhenian Sea, Sardinian Sea and central Adriatic Sea. However, the change in the wind-wave climate we have found (Figure 9 and Table 1) reflects the change in the ERA5 wind forcing intensity and spatial distribution while for instance Lionello and Sanna (2005) study is based on the ERA-40 reanalysis (Uppala et al., 2005) and Morales-Márquez et al. (2020) and De Leo et al. (2020) on the CFSR reanalysis (Chawla et al., 2013). Hence, the reasons for different trends on the intensity of storms should be searched for in the wind forcings that have generated the wave model results. Changes in ERA5 wind characteristics might be due, for instance, to changes in storm trajectories, as suggested by numerous climatological studies, which investigated not only the intensity and frequency of Mediterranean cyclones but also the variability of their trajectories in the last decades and in the future climate projections (Lionello et al., 2002;Lionello and Giorgi, 2007;Cavicchia et al., 2014;Messmer et al., 2020). We have verified monthly H s percentile trends from the ERA5 wind forced wave hindcast at the MS basin scale against trends from satellite altimeters and found a general agreement in the increasing/decreasing trends as well as in the magnitude of the trends. For H max (and C max ), we have found trends that are larger than the H s trends: for instance, for the winter 50th percentile on the whole MS, on average 2.44 times for H max . This can be explained with H max being on average 2.36 H s . Hence, the linear trend of H max , i.e., H max / T, should be on average H max / T ≈ 2.36 H s / T, the difference between 2.44 and 2.36 being possibly explained with a change in the other wave parameters that H max depends upon.

CONCLUSION
In this paper, we have taken advantage of the ERA5 atmospheric wind reanalysis and spectral wave modeling to characterize the wind-wave climate of the MS, its spatio-temporal variability and change. The wind-wave dataset we have produced has allowed us to obtain one of the longest wave climate assessments in the MS to date. Also, for the first time we have bestowed an assessment of the maximum individual wave climate in the MS. Our main results are here summarized: • We have verified the ERA5 wind-based wind-wave hindcast against satellite observations of significant wave height with respect to its performance in reproducing the typical (50th percentile) and extreme (99th percentile) H s , and the variability and trend of the MS wind-wave climate. Despite a general tendency to underestimate H s (in particular in narrow basins as the Adriatic Sea and for the extremes) the dataset has been shown to properly reproduce the temporal variability and the trends of H s . • We have presented the typical and extreme MS wind-wave climate patterns and characteristics at seasonal scale and, to this end, we have proposed and used a definition of seasons based on the satellite observation of significant wave height over the MS, which presents a stormy season (winter) lasting 5 months, a calm season (summer) lasting 3 months and two transitional seasons (spring and autumn) lasting 2 months each. This definition is in agreement with other objective season definitions in the MS. • The largest typical waves (both as H s and maximum individual waves) occur in the western and southern MS during winter, and in the eastern MS during summer, whereas the largest extreme waves occur in the western MS in all seasons, with maximum values during winter. • The intra-seasonal variability of MS wind waves (expressed as the ratio of the 99th to the 50th percentile H s ) has been shown to be largest during summer and in the Adriatic, Tyrrhenian and Sardinian Seas, indicating these are season and the MS sub-basins that are more prone to the development of large extremes waves compared to the typical ones. On the other hand, the inter-annual variability has proven to be largest where the largest winter and summer waves occur. However, when the inter-annual variability has been expressed with respect to the local climate (i.e., in terms of relative anomaly) it has emerged that the Adriatic and Tyrrhenian Seas are characterized by the largest variability of typical and extreme H s , with the Ionian and Levantine sub-basins also showing large variability of extreme H s . During the analyzed period (1981-2019) these regions of the MS have experienced variations up to 100% of both the local 50th or 99th percentile H s . • We have motivated the largest positive and negative relative anomalies of H s in the MS basin (and subbasins) thanks to the geopotential height at 500 hPa, mean sea level pressure and wind intensity anomalies and we have related the temporal variability of the relative anomalies in the MS to the principal modes of atmospheric variability. The seasonal Scandinavian index seems the most correlated to the seasonal wind-wave variability in the MS, especially during winter and for the typical H s , with positive Scandinavian index phases associated to the largest typical winter H s (extreme winter H s are best correlated to the Eastern Atlantic index). However, correlations found are not remarkable (0.54 at most for the MS, 0.65 for sub-basins) and this suggests, in agreement with previous studies in the MS, that the wind-wave variability in the MS can only be partially motivated by teleconnections, most likely due to the effects of local orography that, interacting with synoptic scale atmospheric structures, generates cyclones and consequent winds with regional characteristics that partially lose the link with their largescale source. • The long-term trends found in the ERA5-wind based hindcast of the MS waves are negative for the summer season and positive for the winter season, both for typical and extreme sea states (and maximum individual waves). These trends are statistically significant though generally modest in magnitude (however, locally up to 12 cm decade −1 for extreme winter H s and 25 cm decade −1 for extreme winter H max ) and are larger for extreme sea states compared to typical sea states. This suggests the MS wind-wave climate has changed during 1981-2019, with decreasing H s (and maximum individual waves) during summer and increasing during winter. • The climate characteristics (patterns, spatio-temporal variability and change) of maximum individual waves closely mirror those of H s , although there are some differences that are motivated by the dependence of H max (and C max ) by other sea state characteristics, i.e., the narrow bandedness parameter and the wave steepness. We have presented their typical and extreme seasonal patterns and intensities, showing that MS sea states are generally dominated by wind waves.
Concluding, this study have proven that the ERA5 wind can be successfully used to hindcast the wind waves in the MS, a crucial task to assess the past and present climate in a regional basin providing environmental and economical services to the whole MS region and where the impacts of climate change are expected to be significant in the near future. As a recommendation for future reanalysis-wind based assessments in a semi-enclosed basin with similar characteristics (i.e., narrow sub-basins and surrounding orography), in order to reduce the spatial resolution effects in the narrowest basins (as the Adriatic Sea) and close to the coasts, a dynamical downscaling of ERA5 winds would be advisable. This represents one of the further improvements on the hindcast production side. As regards the wave climate analysis, although swell seas are rather rare in the MS, a partitioning analysis of the directional wave spectra determining the different wave systems that may mix in a sea state would allow treating them separately (if necessary) and also characterizing the crossing seas that often pose serious problems to ships in navigation. Finally, further investigations, for instance looking for a change in the direction of cyclonic paths and their intensity, would help disentangling the dominant mechanisms that drive the changing seasonal wave climate we found.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
FB, SD, and AB contributed to the conception and design of the study and produced the dataset. FB and SD performed the dataset assessment. FB and AR performed the analyses. SD, FF, RF, MS, and AB contributed to the analyses. FB wrote the manuscript. All authors contributed to manuscript revision, read, and approved the submitted version.

FUNDING
This study was partially supported by the Copernicus Marine Environment Monitoring Service (CMEMS) LATEMAR project. CMEMS is implemented by Mercator Ocean in the framework of a delegation agreement with the European Union. This work has also been conducted as part of the bilateral project EOLO-1 ("Extreme oceanic waves during tropical, tropical-like, and bomb cyclones") between CNR-ISMAR and the University of Tokyo (Japan). Antonio Ricchi funding has been provided by PON Ricerca e Innovazione 2014-2020 "AIM"-Attraction and international mobility program. EU Social Fund and Regional Development Fund; Ministero dell'Istruzione e della Ricerca, grant number AIM1858058. Hersbach et al. (2018) was downloaded from the Copernicus Climate Change Service (C3S) Climate Data Store. The results contain modified Copernicus Climate Change Service information, 2020. Neither the European Commission nor ECMWF is responsible for any use that may be made of the Copernicus information or data it contains. We would like to thank Luigi Cavaleri for fruitful discussions on the topic.