1 SIMULATION OF SEA STORMS INCLUDING MULTIVARIATE STORM EVOLUTION

A new method for the simulation of storms is proposed which takes into account the multivariate evolution of the storms, allowing to innovate in the form of each simulated storm, for all the variables involved. The method is based on two novel aspects: (a) measured storms are grouped using clusters techniques and a set of average evolution forms is defined for each cluster, one for each of the variables involved, and (b) a Vector Autoregressive model is fitted to the differences between the average evolution of each variable and the actual measured evolutions. The ability of the methodology to properly reproduce the joint probability distribution of all the variables involved is demonstrated for a case study at the mid Rio de la Plata northern coast.


INTRODUCTION
The simulation of storms of met-ocean variables is a tool commonly used in probabilis tic design and probabilistic coastal risk assessments (see e.g.Martín-Hidalgo et al. 2014, Li et al. 2014, Davies et al. 2017, ROM 1.0-09 and ROM 1.1 to be published).Traditionally, storms simulation methods have been based on the simulation of just a few values per storm, typically the peak significant wave height reached during a storm (i.e.Hs), other simultaneous or concomitant variables (e.g.Tp concomitant with Hs) and the duration of the storm (see e.g Solari et al. 2014).In those cases were it is also required to have information on the storm evolution, as is the cas e for analyzing damage progression or beach erosion, most common approach has been to impose some kind of equivalent form (see e.g.Borgman 1973, Castillo et al. 1977, Boccotti 2000, De Michele et al. 2007).However, with this approach it is not straightforward to reproduce the multivariate joint distribution of the values that are simultaneous to Hs, and at the same time to also reproduce the multivariate joint distribution of all the values that are not simultaneous to Hs (think, for example, on a triangular form for Hs and a rectangular form for the rest of the variables).
Fig. 1 presents a typical example of the temporal simulation process.Once simulated the number and point in time of the storms throughout the year (not shown in the figure), which is usu ally done using a Poisson model, the next step is to simulate the values taken by each one of the variables for each storm.Starting from the distribution of the significant wave height at the storm peak, the duration of the storm and the characteristic values of the other variables (e.g.mean direction, mean period, sea level, wind speed, etc., all simultaneous with Hs,peak) are simulated from joint distributions that relate these variables with Hs,peak.Then, a standard form for the evolution of each of the variables is assumed, typically a triangular or rectangular shape for Hs and a rectangular shape for the other variables, obtaining the time evolution of all the variables along the storm (see e.g.Payo et al 2004, 2008, Baquerizo and Losada 2008).This kind methodologies ensures the correct reproduction of the distribution of Hs,peak, as well as the reproduction of the joint distributions of Hs,peak and the concomitant values of the other variables; however, this methodology does not ensure that the distribution of Hs (all values above the threshold that defines the storm) or the joint distribution of Hs and the other variables (i.e. the joint distribution of all values, not only those that coincide with the peak of Hs).At the same time, it is evident that these methodologies do not innovate in the form (i.e.time evolution) of the storm, which limits the combinations of values of the variables that the method can generate.

OBJECTIVES
The objective of the present work is to devise and test in a case study a methodology for the simulation of multivariate sea storms that fulfills the following requirements.As is the case in many of the existing methodologies, the new one should properly reproduces the multivariate distribution of the variables at storm peak.However, it should also be capable of reproducing the multivariate distribution of the values of all sea states included in the storm (i.e.not only sea states corresponding to the peak of the storm) and capable of innovate in the storm evolution for all variables involved (i.e.storm forms or time evolutions should differ for every simulated storm).

METHODLOGY
The proposed simulation methodology considers that for each storm and for each variable, the temporal evolution of the variable is of the form given by Eq. ( 1) where  is any one of the variables being simulated (e.g.Hs, Tp, etc.),  (  ) is the time evolution of the variable during a given storm,  ̂( ) is the expected (or average) time evolution of the normalized variable,  ̅ ans   are mean value and standard deviation of the values taken by the variable within a given storm and () is a time varying, auto-correlated, innovation process .Fig. 2 shows an outline of the model given by Eq.(1).Thus, the following models need to be defined for performing simulations based on Eq. ( 1): a model for the expected or average time evolution of each variable  ̂( ) , models for mean and standard deviation of each variable,  ̅ and   , and a model for the innovation process ().These models are introduced and discussed below.

Expected evolution and storm duration model
First, a bivariate distribution is fitted for the variables Hs,peak and Duration, assuming that Hs is the main variable and the one defining the occurrence of a storm.For applications where this is not the case, this step could be redefined, fitting a bivariate model for Duration and any other variable or variables of interest.
Secondly, a vector is defined for each storm containing the time evolution of every variable, previously normalized using the mean and the standard deviation of the within storm values, where time was also non-dimensionalized using the storm duration.Then, these vectors are used to define average time evolutions; to this end k-means clusters technique is used.In this way, storm whose normalized and non-dimensional time evolution is similar are grouped into the same cluster, and the average or expected evolution is defined in the normalized and non-dimensional space as the cluster centroids.
Lastly, the occurrence probability of each cluster (i.e. of each average evolution) in the Hs,peak-Duration space is estimated.In this way, given a storm defined by the pair of values (Hs,peak,Duration), it is possible to estimate the probability of that storm having one of the previously defined average evolutions.
Fig. 3 shows, for a case with two cluster, an outline of the probability density function of each cluster in the Hs,peak-Duration space: in that outline the red bivariate distribution corresponds to storms with average time evolution given by  ̂1(  ) and the green bivariate distribution corresponds to storms with average evolution given by  ̂2(  ) .For those storms with duration equal or less than 2 sea states, no average time evolution is defined, resorting to simplified forms.

Innovation model
The difference between the expected evolution of the storms and the actual evolution is modeled using an auto-regressive vector model (VAR) of order one and zero mean, such as that presented in Eq. (2): where   is a vector with values taken at time t by random variables that follow a multivariate normal distribution, that depends on p previous values of  − (with i=1,…,p; in this case p=1 is used) trough linear coefficients contained in matrix   , and on the values taken by Ut , that are i.i.d realizations following a normal multivariate distribution with zero mean.The first step for fitting a VAR model is the estimation of the marginal distributions of the variables (), that are the difference between the average and the actual evolution of the storms (i.e.() is a multivariate variable).Then, following Solari and Van Gelder (2011), these distributions are used to transform () to the normal multivariate variable   .Lastly, a order 1 VAR model is fitted to   .Fig. 4 shows an outline of this process.

Mean and standard deviation models
It is necessary to define models for the mean  ̅ and the standard deviation   for every variable.To this end there are several possible approaches and the selection of the most suitable one will depends on the case study characteristics.In this case it was opted to define the mean  ̅ conditional to Hs,peak for every variable, and to define the standard deviation   conditional to the corresponding  ̅ .The only exception was sea level, for which the mean  ̅ was defined conditional to the mean value of the wave direction.
For the definition of all the conditional distributions a combination of lineal and no -linear regressions was use, along with a Gaussian for the case of the average sea level conditional to the average mean wave direction.

Steps for simulation
The steps followed for the simulation are: 1. Simulate the number of events per year and they occurrence in time (not discussed in this work).then, for every event: 2. Simulate Hs,peak from its marginal distribution (extreme value distribution obtained from a Peak-Over-Threshold analysis).3. Simulate Duration of the storm conditional to Hs,peak. 4. Given Hs,peak and Duration, determine the probability of the storm coming from each cluster and randomly choose one accordingly.5. Simulate mean and standard deviation of each variable conditional to Hs,peak.6. Use Duration, mean and standard deviation to obtain expected storm evolution for each variable 7. Simulate innovation process with the VAR model.8. Add innovations to the expected storm evolution of every variable to obtain the actual values of the variables.

CASE STUDY
For testing the performance of the proposed a 22 years hindcast of three-hourly wave and sea level parameters is used.Data corresponds to a location in the northern coast of the Rio de la Plata estuary (see Fig. 4).Mean water depth in the area is about 7 m; wave conditions are mainly sea, particularly during storms, and total sea level is mainly constituted by storm surges.The variables analyzed are significant wave height (Hs), mean wave period (Tm), mean wave direction (Dm) and total sea level (SL).

RESULTS
Fig. 6 shows the three multivariate, normalized, storm evolution forms identified by means of the kmeans algorithm.The distribution of the events corresponding to each one in the Hs,peak-Duration space is shown in Fig. 7.
In total 100,000 storms were simulated.From them the following results were obtained (all results are discussed in the next section).Fig. 8 shows the comparison of the marginal distributions obtained from the simulated storms with the corresponding empirical ones; in particular it shows the marginal distribution of Hs,peak, and that of Hs and Dm (all values, not only those concomitant with the peak of the storm).Fig. 9 shows the bivariate distributions (Hs,Tm) and (Dm,SL) for both: data corresponding to the peak of the storm (i.e.concomitant with Hs,peak) and all data.Finally, an analysis of the time evolution of the storms is presented in Fig. 10 and Fig. 11.Fig. 10 includes the time evolution of observed severe storms and simulated severe storms with the same range of Hs,peak values.Fig. 11 shows the multivariat e time evolution for 100-year return period storms (for which there are no hindcast data available for comparison).

DISCUSSION
The three average storm evolution forms identified by the methodology (Fig. 6) are clearly different.While forms  1 ̂ and  2 ̂ differ mainly in the evolution of the wave direction: in  1 ̂ they rotates clockwise and in  2 ̂ they rotate counterclockwise,  3 ̂ clearly differentiate from the others in the evolution of the sea level.It is noted that none of the average form would be properly reproduced by a triangular, nor by a rectangular form.From Fig. 7 it is seen that identified clusters not only differ in terms of its sto rm evolution forms but also in terms of its distribution in the Hs,peak-Duration space, with average form  2 ̂ tending to occur with lager Hs,peak and lower durations that  1 ̂ and  3 ̂.The simulations properly reproducing the distribution of Hs,peak (Fig. 8, a) should not be a surprise as this variable is directly simulated from its marginal distributions.Moreover, it is shown that the simulations properly reproduce marginal (Fig. 8, b and c) and bivariate (Fig. 9) probability distributions that are not used in the simulations.In all these cases the variables are simulated indirectly, through the simulation of the storm evolutions, their mean and standard deviatio ns and the innovations.The fact that present methodology is able to properly reproduce these distributions without simulating directly from them is an indicator that the general model given by Eq. ( 1) is a reasonable representation of the actual behavior of the storms.
In some cases the simulation missed some particular features of the bivariate distribution, as is the case in Fig. 9 (lower right panel).This is believed to be due to simplifications made when modelling the relation between average direction (  ̅̅̅̅ ) and average sea level ( ̅̅̅ ).As mentioned, quite simple models were used for modelling  ̅ and   , leaving some room for future improvements of the results.From Fig 10 it is seen that average storm evolution obtained from the simulations properly re produce those observed, but with an increased variability, as evidenced by the confidence intervals.Moreover, from Fig. 11 it is seen that for the central sea states the average evolution of Hs is more or less triangular, but that the evolution of the other variables are not necessary rectangular, not even in average terms; also, when confidence intervals are taken into account, it is concluded that innovation capability of the methodology is able to produce storm evolutions that clearly depart from triang ular or rectangular standard forms.

Figure 1 .
Figure 1.Exam ple of typical storm simulation method (see main text for explanation).

Figure 2 .
Figure 2. Outline of the model proposed in Eq. (1): (a) expected evolution of the normalized variable during the storm (t* is a non-dimensional time, t*=t/D, being D the duration of the storm); (b) and (c) transformation from norm alized variables to actual variables (i.e.Hs, Tp, etc.); (d) adding of the innovation process to obtain the actual evolution of the storm (in red) from the expected or average evolution (in dashed blue).

Figure 3 .
Figure 3. Outline of the m odel proposed for duration and expected time evolution of the multivariate storms.

Figure 3 .
Figure 3. Outline of the determination and normalization of the innovation process.

Figure 5 .
Figure 5. Approximate location of the data used for the case study.

Figure 6 .
Figure 6.Average multivariate storm evolution; variables are normalized w ith mean and standard deviation of each storm; time is non-dimensionalized w ith storm duration.

Figure 7 .
Figure 7. Distribution of the three storm evolution clusters in the Hs,peak-Duration space.

Figure 11 .
Figure 11.Left panels: identification of 100-years range of values for Hs,peak (upper) and simulated storms in that rage in the Hs,peak-Duration space (lower).Right panels: Average and 80% confidence intervals for the evolution of Hs, Dm and SL, and frequency of reaching time lags relative to the storm peak, from simulated storms in the 100-years return period range of values.