- Email: [email protected]
Long-term shoreline response of a nontidal, barred coast
Coastal Engineering 52 (2005) 79 – 91 www.elsevier.com/locate/coastaleng
Long-term shoreline response of a nontidal, barred coast Grzegorz Ro´yyn´ski* Institute of Hydroengineering, Polish Academy of Sciences, 7 Kos´cierska, 80-953 Gdan´sk, Poland Received 16 February 2004; received in revised form 22 July 2004; accepted 22 September 2004 Available online 18 November 2004
Abstract Long-term variations of shoreline positions along the southern Baltic coast were investigated using multichannel singular spectrum analysis (MSSA) to determine the most dominant long-term response patterns. The investigated beach is located at Lubiatowo on the Polish Coast and is mildly sloping with multiple bars. Data on coastal morphology have been collected at Lubiatowo including (1) bathymetric surveys since 1987 twice a year, and (2) beach topography surveys since 1983 every 4 weeks on the average, extending from the shoreline to the dune foot. Furthermore, several dedicated field campaigns have been carried out at Lubiatowo, as well as measurements of deep-water wave properties since 1998. MSSA was employed to the whole data set of shoreline position from all survey lines. In summary, three patterns emerged reproducing alongshore standing waves with different periods 7 to 8, 20+ and several decades. They represent long-term shoreline response, such that at some locations the longest wave is most predominant, at other locations the medium cycle predominates, whereas the shortest is the most prominent at yet other locations. However, all three can be detected at every location monitored, eliminating the confusion resulting from ordinary singular spectrum analysis (SSA) analysis, done previously for the same data set. D 2004 Elsevier B.V. All rights reserved. Keywords: Data-driven modeling; Shoreline evolution; Multichannel singular spectrum analysis; Lagged covariances; Reconstructed components
1. Introduction Coastal processes include complex, nonlinear phenomena and involve morphological features that are characterized by a wide range of spatial and temporal scales (De Vriend, 1991, Larson and Kraus, 1995). Despite substantial progress in the research of * Tel.: +48 58 552 2011; fax: +48 58 552 4211. E-mail address: [email protected]. 0378-3839/$ - see front matter D 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.coastaleng.2004.09.007
physical mechanisms that shape these morphodynamic processes, present knowledge is still too limited to properly understand, describe, and forecast the behavior of many coastal systems (Stive et al., 2002). Fortunately, high-quality data sets on coastal morphology become more easily available now, providing an opportunity for extracting valuable information on morphological properties by means of computerintensive analysis and data-driven modeling (Larson et al., 2003).
80
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
In analysis of morphological data, it is common practice to split the state vector (i.e., the morphological data set) into two subvectors: the signal vector which is described by a few characteristic patterns believed to represent the dynamics of the considered process, and the noise vector which contains all processes which seem to be unrelated to the signal vector. For example, linear techniques like empirical orthogonal functions (EOF) are optimal in representing the variance. Therefore, the EOF analysis has been widely used by researchers in coastal morphology to identify patterns in the data (e.g. Winant et al., 1975, Aubrey, 1979, Wijnberg and Terwindt, 1995, Larson et al., 1999). Another linear method, namely canonical correlation analysis (CCA) maximizes the correlation between two simultaneously observed vector time series (Barnett and Preisendorfer, 1987, Ro´z˙yn´ski, 2003). Currently, a lot of efforts is devoted to introduction of nonlinear techniques in order to establish better representation of key morphological processes that are essentially nonlinear. Such methods usually involve the combination of at least two signal processing techniques in which ordinary EOF method usually serves as a data compression tool. For example, Ruessink et al. (2004) used neural network based circular nonlinear principal component analysis to study similar bar behavior at coastal segments in the Netherlands, Japan and USA. Other less known methods, like principal oscillation or interaction patterns, satisfy certain dynamical constraints (Hasselmann, 1988, Ro´z˙yn´ski and Jansen, 2002), which might be identified with some a priori knowledge of a studied system. Ordinary singular spectrum analysis (SSA) was employed to a limited extent for investigating coastal morphology and its evolution. Southgate et al. (2003) utilized SSA to extract long-term trends from measurements of the variation in shoreline position at several beaches around the world with the overall purpose of detecting oscillatory behavior in the filtered signal (also, compare Stive et al., 2002). Ro´ z˙yn´ski et al. (2001) applied this method to determine characteristic patterns in the shoreline response at a sandy beach in the southern Baltic Sea seeking patterns related to forced or self-organized behavior. That study revealed three significant patterns. The most important always represented longterm shoreline change, but its period varied from
location to location, so finally periods of some 9, 16 and around 30 years were identified. The two other patterns were identified to represent deterministic chaos and thus were found to reflect self-organized shoreline behavior. The current study is a direct extension of that SSA analysis and is aimed at a more thorough analysis of the long-term behavior. In particular, the study is focused on understanding why ordinary SSA demonstrated different predominant cycles of long-term shoreline response at different points of the monitored area. Analysis of the whole data set, including its cross-covariance structure was thus intended to come up with a more comprehensive idea of that part of shoreline change. Long-term stability of the studied area provides excellent conditions for application of the MSSA method. First, the field site at Lubiatowo Beach in Poland is described together with the morphological data collected. Subsequently, description of SSA and MSSA is given. Then, MSSA is applied to the whole data set (27 survey lines analyzed jointly) on shoreline position from Lubiatowo and patterns of variation in time and space are established and discussed.
2. Field site and data collected The Coastal Research Station (CRS) of the Institute of Hydroengineering, Polish Academy of Sciences, is located at Lubiatowo on the Polish Coast facing the southern part of the Baltic Sea (see Fig. 1). The beach usually features multiple longshore bars and has a mild mean slope of m=tana=1–1.5% with a median grain size of D 50=0.22 mm. Since the resultant longterm energy flux is oblique to the shore, a predominant west to east littoral drift is produced. Multiple breakers are frequently observed in the surf zone. Long-term records indicate that for average storms the significant wave height outside the surf zone (in water depth h of about 20 m) usually reaches H s=2–2.5 m with a mean period of T=5–7 s. The incoming wave energy is dissipated as the waves propagate onshore, and for h=2–3 m the mean height typically equals – H =0.5–1 m with T=4 –5 s (Pruszak et al., 1999). Closer to the shoreline (h=1 m), the height of the waves during storms reduces to 0.3–0.5 m. Since the Baltic Sea is practically isolated from the Atlantic Ocean, there is no influence from tides. The average
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
81
Fig. 1. Location of Coastal Research Station at Lubiatowo.
surf scaling parameter, defined as e=[(H bx 2)/(gm 2)], where x=2p/T, g the acceleration of gravity, m the average beach slope, and H b the wave height at incipient breaking, in general falls between 1000 and 3000, indicating that the coast at Lubiatowo is a highly dissipative system. Due to the mildly sloping bottom, the values of e correspond to low p values offfi ffiffiffiffiffiffiffiffiffiffiffiffi the usual surf similarity parameter n ¼ m= H0 =L0, which falls below 0.1 for the quantities H 0 and L 0 representing deepwater wave characteristics. Long-term bathymetric surveys have shown that the beach at CRS Lubiatowo normally includes a system of four stable bars (Ro´z˙yn´ski et al., 1999). The first, innermost bar is situated about 120–170 m from the baseline with depth over crest of 1 m, the second 220–300 and 2 m, third 400–500 and 3–3.5 m, and the fourth, outermost bar 600–800 and 4–5 m. Occasionally, an ephemeral bar develops in the vicinity of the shoreline, typically located some 20–80 m from the
baseline. The bars do not migrate, but only oscillate around their average locations (Pruszak et al., 1997, 1999). This observation was confirmed by measurements of nearshore bathymetry commenced in the early 1960s, which gained geodetic compatibility since 1983, when the current baseline was established. The baseline itself is situated 10–50 m from shoreline on a dune crest. Very rarely, it can be destroyed by a series of extreme events. In such instances, exact baseline positions need to be reestablished by measurements referenced to a local state geodetic system. The surveys of the beach topography, from the shoreline to the dune foot, started in 1983. Although the measurements indicate long-term beach equilibrium, the shoreline is not a perfectly stable feature. It was observed that some parts of the shoreline advanced, while others receded or displayed no movement. However, after some time the direction
82
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
to obtain equal sampling interval in the MSSA, data points were derived for these gaps using linear interpolation, which showed little impact on longterm patterns though, which can be explained by the fact that the duration of gaps is much shorter than the identified long-term cycles. As previously pointed out, no net long-term erosion or accretion has been observed at CRS Lubiatowo, although the shoreline regularly moves back and forth over a decade or so (Pruszak, 1993). Thus, these measurements are highly suited for extracting temporal and spatial patterns of shoreline change for a beach where oscillations occur but without any net movement. The available data set covers the period between Sept. 1983 and Sept. 1999, which makes it unique in duration in a sense that there are only few shoreline observations, whose time span covers more than a decade (Dutch Jarkus data set or FRF Duck, NC, USA). The spatial extension of more than 2.5 km is also an important factor. Together they allow for comprehensive analysis in both time and space resulting in extraction of large shoreline standing waves having periods of many years.
of movement for some parts may reverse, whereas others may be stable or commence oscillating back and forth. Over long periods all these combinations of movements are possible for a particular shoreline stretch and the mechanisms governing these phenomena appear to be quite complex. The surveys are attached to a local geodetic system, so all measurements taken at different times and along different lines are fully consistent. In addition, shoreline positions at the profile lines are all related to the mean sea level. In total, the surveyed area covered 2600 m of alongshore distance with 27 cross-shore profiles uniformly spaced at 100 m. Enumeration of these profiles originates from the history of the CRS Lubiatowo geodetic base (cf. Fig. 2). Initially, eight baseline points east of point 3 were fixed (points 3, 4, . . ., 10), after which the baseline was extended westwards to fix points 11, 12,. . . , 29. Thus, the baseline consists of two segments that meet at fix point 3 and because of the dune configuration their azimuths differ by 4811V38.4U. Time series of shoreline positions measured at all 27 survey lines were employed in the MSSA. The analyzed data set covers 16 years of measurements. Surveys started on 22nd Sept. 1983 and typically have been repeated every 4 weeks (28 days), and the last survey included in this study was carried out on 2nd Sept. 1999. However, the data set includes two longer time periods for which no measurements were taken: there is one gap in the data between 15th Dec. 1992 and 25th Nov. 1993, and another gap between 17th Aug. 1994 and 3rd Jan. 1996. In order
3. Ordinary and multichannel singular spectrum analysis Ordinary singular spectrum analysis (SSA) can be regarded as a particular application of EOF analysis where the data matrix contains values measured at a
N
Fig. 2. Geodetic base at CRS Lubiatowo.
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
specific location lagged in time (Vautard et al., 1992). Therefore, the column vectors of the data matrix contain the measured time series successively delayed in time up to the maximum shift known as the embedding dimension. The eigenvectors of the corresponding covariance matrix yield the dominant patterns in the time series under study. In contrast to traditional EOF analysis, SSA provides information on the dynamics of the underlying system (Von Storch and Navarra, 1995). Actually, it can be shown that the complete dynamics of the underlying system can be reconstructed from the measured data without any knowledge of the evolution equations (Broomhead and King, 1986). The SSA method is based on the Karhunen-Loe´ve expansion theorem, where the lagged covariance matrix of the analyzed time series x i , where 1ViVN, N—series length, is utilized in the expansion. This matrix is often called the Toeplitz matrix as it has the Toeplitz structure with constant diagonals corresponding to equal lags, 2 3 c ð 0Þ c ð 1Þ N c ð M 1Þ 6 c ð 1Þ 7 c ð 0Þ c ð 1Þ : 7 Tx ¼ 6 ð1Þ 4 5 : c ð 1Þ c ð 0Þ : c ð M 1Þ N c ð 1Þ c ð 0Þ where: cð j Þ ¼
N j
1 X xi x Þ xiþj x N j i¼1
ð2Þ
and –x is the mean value of the analyzed time series. The number M is the user defined embedding dimension and its selection is crucial for the analysis results. The larger M is, the better the spectral resolution of oscillatory components in the time series. However, to eliminate statistical errors in the autocovariance function for large j, it is recommended that M should not be greater than (1/3)N with the spectral resolution M 1. The Toeplitz matrix is symmetric, with positive terms on the main diagonal, so all its eigenvalues k k are positive, unless the data are perfectly noise-free and originate from a dynamical system with purely quasi-periodic, deterministic behavior. In such a case, some eigenvalues equal zero and the matrix is singular. The corresponding eigenvectors E k are
83
orthogonal, so the sum of products of the corresponding terms of the (column) eigenvectors E j and E i is equal to the Kronecker delta; M X
Ekj Ekl ¼ djl ; 1VjV M ; 1V lV M
ð3Þ
k¼1
Interestingly, the (column) eigenvectors of the Toeplitz matrix are either symmetric or antisymmetric with respect to (1/2)Mth term in each column. Eq. (3) denotes summation along columns. The same holds true for the sum of products along the jth and lth rows: M X
Ekj Ekl ¼ djl ; 1V jV M ; 1V lV M
ð3aÞ
k¼1
In addition, the spectral decomposition formula is valid as well: M X
kk Ekj Ekl ¼ Tx; jl ¼ cðjj ljÞ; 1VjVM ; 1VlVM
k¼1
ð4Þ The eigenvectors form the time-invariant part of the SSA decomposition, whereas the variability of a given system is contained in the principal components (PCs). The kth PC is a projection coefficient of the original signal onto the kth eigenvector: aki ¼
M X
Xiþj1 Ekj ; 1V iV N M þ 1
ð5Þ
j¼1
According to this formula, we have to take M elements of the original x series from the ith to i+Mth element, compute their products with the corresponding elements of the kth (column) eigenvector and sum them up to obtain the ith element of the kth PC. Hence, the PCs are time series of length NM. Even though the PCs are orthogonal to each other, it does not mean they are independent, because Eq. (5) shows that M consecutive elements of the original time series are needed to compute one term of every PC. Therefore, there are k common elements of this series for the ith term of the rth PC ari and the jth term of the sth PC a sj , such that k=M|ji|N0 (lag|ji|). Hence, the correlation structure of the original series must be imprinted in the sequence of the PC terms, producing nonzero
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
84
correlations for nonzero lags. For natural systems, where the assumption of normally distributed data generated by the system is widely accepted, it means that time series of such data are dependent on each other. The PCs in the SSA decomposition are not convenient, since they do not provide a unique expansion of the signal into a sum of different components. Using the PCs, we can expand our series as: xiþj ¼
M X
aki Ejk ; 1V j V M
ð6Þ
k¼1
There may be up to M subsets of the original time series containing the specific element x i+j . Thus, there are up to M different ways of reconstructing the components of the signal with Eq. (6), which in general do not produce the same results. Furthermore, instead of a series of length N, a series having NM+1 elements is obtained, so it is difficult to apply Eq. (6) for filtering and prediction. It is possible though to construct a least-square optimum series corresponding to a given (sub)set of the eigenelements. If the following expression is minimized, !2 N M M Xþ1 X X HW ð yÞ ¼ Yiþj aki Ekj ð7Þ i¼1
j¼1
kaW
where y is the desired series and W is the subset of eigenelements, the solution y=R W x is given by: ðRW xÞi ¼
M X 1 X ak Ek M j¼1 kaW ijþ1 j
for M V iV N
M þ1
ð8aÞ
ðRW xÞi ¼
i X 1X ak Ek for 1 V iV M 1 i j¼1 kaW ijþ1 j
ð8bÞ
ðRW xÞi ¼
1 N iþ1
M X
X
akijþ1 Ekj
j¼1N þM kaW
for N M þ 2ViVN
ð8cÞ
If W contains only a single index k, the resulting series is called the kth reconstructed component (RC), which will be denoted by x k , the RCs are thus
P k additive, i.e., RW x ¼ kaW x . Hence, the series x can be uniquely expanded as the sum of its RCs: x¼
M X
xk
ð9Þ
k¼1
To explain Eqs. (8a) (8b) (8c) in more detail, let us trace the computations of the ith element of the 1st RC. First, all combinations of the terms of the 1st PC (a1j ) and the 1st eigenvector (E k1) such that j+k=i+1 are multiplied, then their sum is taken and finally its mean value provides the desired term of the 1st RC: other RCs require similar computations on the underlying PCs and eigenvectors. The transformation of the original series x into x k is nonlinear, for the eigenvectors E k depend on x in a nonlinear fashion. Furthermore, the RCs are correlated even at lag 0, so additional correlation analysis is required for the extracted RCs. Consequently, variances of RCs are not cumulative. In some cases measurement time series from more than one location are analyzed simultaneously leading to MSSA. MSSA is mathematically equivalent to extended EOF (EEOF) analysis (Weare and Nasstrom, 1982), although EEOF is used for data sets with many points in space and only few lags in time. Apart from autocovariances of all time series included in an MSSA analysis, their cross-covariances are taken into account as well. It is schematically described by expression (10) that presents a block MSSA system matrix made up from covariances on the main diagonal and cross-covariances below and above. The Ti,i elements represent the autocovariance structure of ith time series and their structure is identical to Tx from expression (1). Consequently, the Ti,j terms describe cross-covariances between ith and jth time series from lag 0 up to the embedding dimension M. If the number of time series in an MSSA study is equal to L, the unfolded block matrix has ML ML terms and the embedding dimension of the whole problem is ML. Like in the case of ordinary SSA analysis, the MSSA matrix is symmetric, cf. Plaut and Vautard (1994), so both eigenstructures are formally identical to each other. Thus, using the rationale expressed through Eqs. (8a) (8b) (8c), we can calculate reconstructed components of an MSSA study. We should remember though that unlike in the ordinary SSA analysis, where the analyzed time series is
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
85
resolved into M reconstructed components, being time series as well, the MSSA reconstructed components are random fields. Thus, as a result of MSSA analysis we obtain ML random fields. This may make its interpretation quite complicated, since the initial field having NL elements is decomposed into MNL numbers, which can rapidly turn into a huge value. 2 3 T1;1 T1;2 : : : T1;L 6 T2;1 T2;2 : : : T2;L 7 7 T¼6 ð10Þ 4::: ::: ::: ::: 5 TL;1 TL;2 : : : TL;L
4. Analysis of shoreline change In order to analyze the temporal and spatial variation in shoreline position with the MSSA method, the measurements from all 27 survey lines at CRS Lubiatowo were employed. All further references to a dlineT or dprofile lineT will mean shoreline positions at that line. The time series from each line encompassed 169 individual measurements in addition to 28 values obtained through linear interpolation for the gaps where data were missing. Table 1 illustrates the mean value and standard deviation for all the studied profile lines, whereas Fig. 3 displays the entire data set, also indicating the
Table 1 Mean value and standard deviation (S.D.) of the analyzed time series of shoreline position at respective profile lines Profile number
Mean value (m)
Standard deviation (m)
Profile number
Mean value (m)
Standard deviation (m)
29 28 27 26 25 24 23 22 21 20 19 18 17 16
69 59.1 49.5 38.7 29.6 21.9 14.9 9.9 7.6 4.6 3.9 3.8 4.7 6.9
19.8 19.5 17.2 14.2 13.3 14 14 13.7 11.1 12.4 14.7 13.6 14.2 13.5
15 14 13 12 11 3 4 5 6 7 8 9 10
10.3 15.4 19.4 22.3 22.5 32.3 40.6 48.2 57.5 67 73.1 77.2 80.7
13.5 12.1 11.7 12.3 12.4 11.9 11.5 9.8 10.6 13.3 13.5 12.6 13.1
Fig. 3. Envelope of shoreline records at CRS Lubiatowo, 1983–1999.
envelope that contains the variation in shoreline position between 1983 and 1999 for all lines. Table 1 demonstrates that the overall behavior of the shoreline quantified in terms of standard deviation is quite similar, except for the two westernmost profiles 29 and 28 (partly 27 also), which show a larger scatter. Also, Fig. 3 shows that the variability is more or less the same all over the study area with a similar difference between the maximum and minimum shoreline position for all lines. The first step in the MSSA encompassed the choice of the maximum time lag M. The value of M=30 was selected, which is considerably less than the maximum length possible. This value is equal to the time lag in the ordinary SSA study, cf. Ro´z˙yn´ski et al. (2001), and was selected mostly to reduce the amount of computations, knowing from that analysis this choice includes all long-term behavior. It should be stressed here that the MSSA analysis would have been hardly possible, had the ordinary SSA study not been executed before. Hence, it highlights a synergistic character of ordinary SSA and MSSA. Fig. 4 presents the 1st RC of the MSSA analysis for all profiles. For lines 29–25 in the westernmost segment of the monitored area the 1st RC is very similar to the results of individual SSA analyses of these lines and varies in aF20 m range about the mean empirical shoreline position. For other pro-
86
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
Fig. 4. Shoreline standing wave with T=several decades and L=1500 m.
files, the 1st RC is not as spectacular and the presence of this pattern could not be found with the SSA method. Nevertheless, the range of variations in that portion of the studied coastal segment falls between +4 and 12 m. The variations have a character of a long standing wave; the lines for profiles 21 and 07 can be regarded as the nodes of that wave, because of their relatively small overall variations. The distance between these two profiles L=1500 m can be interpreted as the standing wavelength. Importantly, there is also a temporal node for all profiles (day 3800, i.e., year 1994). The period can be deduced from the trajectories of profiles 29 and 28 confronted with the lines 20 and
19. They represent the switch from maximum erosion to accretion of the 1st RC (lines 29 and 28) and vice versa (lines 20 and 19). As it occurred during 16 years (1983–1999), we can say the period of this standing wave is about 32 years. However, since this reasoning is entirely based on visual assessment of the 1st RC that embraces only 16 years of observations, it is safer to judge the period of this standing wave is several decades. The explanation of what drives this wave is not possible, because the time span of bathymetric surveys of nearshore bed topography was too short to find a pattern in nearshore bathymetry with a similar period. Nevertheless, it indicates that such a pattern
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
might be detected in the future when more observations of seabed in the nearshore region is available, which points to the direction of future research. Fig. 5 demonstrates trajectories of the 2nd MSSA reconstructed components. They can be interpreted as standing waves with a period of 7–8 years, with nodes at profiles 27–26 and 16–15. The antinodes can be observed at profiles 22–19 and 03–06. Therefore, a wavelength of 1000 to 1400 m can be identified from the locations of nodes and antinodes. The amplitudes are of similar magnitudes (F10 m), from their
87
extremes the period of 7–8 years could be established. Interestingly, this period coincides with the period of the climatic North Atlantic Oscillation (NAO), so it is possible that its influence is imprinted in shoreline variations of an isolated Baltic Sea. What is important, similar pattern was found for a number of profile lines with the SSA routine, where it was the most dominant feature. Like for the 1st RCs, the MSSA showed it appears all over the studied coastal segment as a standing wave. The origin of this wave and its relation to the patterns describing the evolution of nearshore bed
Fig. 5. Shoreline standing wave with T=7–8 years and L=1200–1400 m.
88
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
Fig. 6. Shoreline standing wave with T=20–22 years and L=1400–1600 m.
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
topography is unclear for the time being. Similar period was found with a traditional EOF routine evaluating jointly cross-shore bathymetric patterns at lines 4, 5, 6 and 7, cf. Ro´z˙yn´ski (2003). This pattern describes oscillations of two inner bars that bounce between the shoreline and outer bars. Unfortunately, a coupling of such oscillations with the alongshore standing wave cannot be established conclusively, despite their strikingly similar periods, because separate EOF analyses of individual cross-shore profiles were unable to trace inner bars oscillations, which if related with the shoreline standing wave, should be phase shifted from profile to profile. Nevertheless, the resemblance of periods of this shoreline standing wave and the NAO shows the avenue of future research, which should be aimed at better identification of the influence of large climatic oscillation on apparently isolated Baltic Sea. The 3rd MSSA pattern is presented in Fig. 6. Upon its inspection, this pattern demonstrates yet another standing wave with nodes at profiles 23, 22 and 21 and then at 05 and 06. The antinodes can be observed at profiles 16 and 17 and then at 08, 09 and 10. The temporal nodes are visible at day 1700 (1988) and 5500 (1999). From these trajectories, we can estimate the wavelength as to be between 1400 and 1600 m and the period to be between 20 and 22 years. The magnitude of amplitudes of this standing wave is similar to the amplitudes of the previous two and equals 14 to +10 m at lines in the center of the study area, where it is best pronounced, to F6 m at other lines. Like previously, this pattern was identified with the SSA routine at profiles where it was the most conspicuous (profiles 18–14). As it was with the 2nd MSSA pattern, it was not possible to directly relate this wave to changes in nearshore bathymetry. Being consistent with the results of SSA analyses, we can say that all three patterns represent the forced shoreline response to the wave climate, which consists of standing waves with different wavelengths and periods. The only difference is that in the SSA study it was assumed that long-term shoreline variations represent the response to changes in the wave climate. In light of the MSSA study, the wave climate in general can cause the shoreline to perform standing oscillations. Still, the effects of long-term wave climate fluctuations can be
89
traced as well (2nd pattern in particular). This is a real progress in the understanding of long-term shoreline evolution. Interestingly, locations of nodes and antinodes of these waves appear not to follow any consistent pattern and the evolution at particular spatial points is usually dominated by the most pronounced standing wave at these points, giving the impression that long-term shoreline evolution can vary substantially from one subsegment to another. The MSSA analysis found that in fact the shoreline evolution contains three periodic components at all points. The presented study revealed initially unexpected layout of long-term shoreline response (standing waves), which was found in the described three reconstructed components of shoreline change. The analysis of trajectories of higher RCs detected chaotic behavior, already known from conventional SSA analysis. Their interpretation though brought no new discoveries in relation to that study. For this reason, only three top MSSA reconstructed components were analyzed in detail.
5. Conclusions Multichannel Singular Spectrum Analysis (MSSA) was applied successfully to analyze shoreline evolution at a beach that displays no net long-term erosion or accretion. The MSSA found that the major response consists of three standing wave components with the periods of several decades, 20–22 and 7–8 years for all spatial points. Usually, one of them is much more pronounced than the remaining two at a given point, so conventional SSA routine was unable to extract them. Consequently, it was assumed in the SSA study that the forced component originates from wave climate fluctuations, whereas the MSSA results tend to conclude standing wave components reflect the wave climate in general as well. The present analysis clearly manifests the complementarity of SSA and MSSA methods, illustrated by the fact that interpretation of MSSA patterns is usually very difficult without extra a priori knowledge. This is relatively obvious, when we recall expression (10) containing covariances and cross-covariances among the analyzed time series. The SSA method can
90
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
provide such a priori knowledge by highlighting covariance structures of individual time series. Then, the analysis of the cross-covariance structure of the whole random field in question becomes feasible and more accurate. No relationship was found between the longest standing wave component (L=1500 m, T=several decades) and the evolution of nearshore bathymetry. To achieve that records of nearshore bed topography must be continued until they cover a similar time period. Still, the existence of the shoreline pattern lasting more than three decades indicates we can expect seabed patterns having comparable periods. Consequently, we should try to single out potential drivers of such long-term changes (climatic variations, long-term developments in the Atlantic Ocean, etc.). Thus, the discovery of such a pattern points to where future research should be focused. The coincidence of the period of most frequently encountered 2nd standing wave component (T=7– 8years) with the period of North Atlantic Oscillation (NAO) suggests NAO effects can be imprinted in the shoreline evolution of an isolated sea, such as the Baltic Sea. This is perhaps the most important conclusion, because the existence of periodicity in shoreline change that matches the NAO period points to NAO as a (rather unexpected) driver of morphodynamic evolution in apparently isolated seas. On second thoughts though, it shows that the scale and intensity of NAO cannot be ignored when studying long-term beach processes in such seas. Therefore, possible linkage between NAO and nearshore morphodynamics in the Baltic Sea also points to future research strategies. No relationship was also found for the 3rd standing wave component (T=20–22 years) and nearshore bathymetry. The understanding of its generation, as well as the generation of two other shoreline standing waves, may only become feasible in future when records of nearshore bathymetry are confronted with large enough hydrodynamic records of deep and shallow water waves, nearshore currents and water levels. This should result in better understanding of the feedback between hydrodynamics and nearshore bed topography with multiple bars, which together produce long-term shoreline change in the form of standing waves with different periods.
Notation aki ith element of kth principal component of SSA decomposition c( j) autocovariance of a time series at lag j D 50 mean sediment grain size Ekj jth element of kth eigenvector of Toeplitzmatrix Tx g acceleration of gravity – H, H s, H b, H , H 0 wave heights: measured, significant, at break-ing, average, and deepwater, respectively H W ( y) a least-square minimized quantity to derive reconstructed components h water depth L number of time series in MSSA analysis L0 deepwater wavelength M embedding dimension (window length) m=tana mean slope N number of observations in a time series xi shoreline position at time i, where i=1,. . ., N x¯ average shoreline position xk kth reconstructed component Tx M M system matrix in SSA T ML ML system matrix in MSSA T wave period, x=2p/T y=R W x least-square optimal series corresponding to W d j,l Kronecker delta e surf scaling parameter n surf similarity parameter kk kth eigenvalue of Toeplitz matrix W subset of eigenelements used in derivation of reconstructed components
Acknowledgements The author is indebted to the Commission of the European Community in the framework of the Energy, Environment and Sustainable Development (EESD) Thematic Programme, for its support within the Project: bHuman Interaction with Large Scale Coastal Morphological EvolutionQ HUMOR-EVK3CT-2000-00037. The author is also indebted to Professor M.J.F. Stive and Professor B.G. Ruessink for their in-depth reviews that gave grounds for improvement of the final manuscript in order to meet the journal standards.
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
References Aubrey, D.G., 1979. Seasonal patterns of onshore/offshore sediment movement. Journal of Geophysical Research 84 (C10), 6347 – 6354. Barnett, T.P., Preisendorfer, R., 1987. Origins and levels of monthly and seasonal forecast skill for United States surface air temperatures determined by canonical correlation analysis. Monthly Weather Review 115, 1825 – 1850. Broomhead, D.S., King, P.G., 1986. Extracting qualitative dynamics from experimental data. Physica. D 20, 218 – 236. De Vriend, H.J., 1991. Mathematical modeling and large-scale coastal behavior: Part 1. Physical processes. Journal of Hydraulic Research 29 (6), 727 – 740. Hasselmann, K., 1988. PIPs and POPs: the reduction of complex dynamical systems using principal interaction and oscillation patterns. Journal of Geophysical Research 93 (D9), 11015 – 11021. Larson, M., Kraus, N.C., 1995. Prediction of cross-shore sediment transport at different spatial and temporal scales. Marine Geology 126, 111 – 127. Larson, M., Hanson, H., Kraus, N.C., Newe, J., 1999. Short- and long-term responses of beach fills determined by EOF analysis. Journal of Waterway, Port, Coastal, and Ocean Engineering 125 (6), 285 – 293. Larson, M., Capobianco, M., Jansen, H., Ro´yyn´ski, G., Southgate, H., Stive, M., Wijnberg, K., 2003. Analysis and modeling of field data on coastal morphological evolution over yearly and decadal time scales: Part 1. Background and linear techniques. Journal of Coastal Research 19 (4), 760 – 775. Plaut, G., Vautard, R., 1994. Spells of low-frequency oscillations and weather regimes in the northern hemisphere. Journal of Atmospheric Sciences 51 (2), 210 – 236. Pruszak, Z., 1993. The analysis of beach profile changes using Dean’s method and empirical orthogonal functions. Coastal Engineering 19, 245 – 261. Pruszak, Z., Ro´yyn´ski, G., Zeidler, R.B., 1997. Statistical properties of multiple bars. Coastal Engineering 31, 263 – 280. Pruszak, Z., Ro´yyn´ski, G., Szmytkiewicz, M., Skaja, M., 1999. Quasi-seasonal morphological shore evolution response to
91
variable wave climate. Proc. Coastal Sediments ’99, Conference, ASCE 1081-1093. Ro´yyn´ski, G., 2003. Data-driven modeling of multiple longshore bars and their interactions. Coastal Engineering 48, 151 – 170. Ro´yyn´ski, G., Jansen, H., 2002. Modeling nearshore bed topography with principal oscillation patterns. Journal of Waterway, Port, Coastal, and Ocean Engineering 128 (5), 202 – 215. Ro´yyn´ski, G., Pruszak, Z., Okro´j, T., Zeidler, R.B., 1999. Depth of closure and seabed variability patterns. Proc. ICCE’98, ASCE 2926-2939. Ro´yyn´ski, G., Larson, M., Pruszak, Z., 2001. Forced and selforganized shoreline response for a beach in the southern Baltic Sea determined through singular spectrum analysis. Coastal Engineering 43, 41 – 58. Ruessink, B.G., van Enckevort, I.M.J., Kuriyama, Y., 2004. Nonlinear principal component analysis of nearshore bathymetry. Marine Geology 203, 185 – 197. Southgate, H.N., Wijnberg, K.M., Larson, M., Capobianco, M., Jansen, H., 2003. Analysis of field data of coastal morphological evolution over yearly and decadal time scales: Part 2. Non-linear techniques. Journal of Coastal Research 4 (4), 776 – 789. Stive, M., Aarninkhof, S.J.C., Hamm, L., Hanson, H., Larson, M., Wijnberg, K., Nicholls, R.J., Capobianco, M., 2002. Variability of shore and shoreline evolution. Coastal Engineering 47 (2), 211 – 235. Vautard, R., Yiou, P., Ghil, M., 1992. Singular spectrum analysis: a toolkit for short, noisy, chaotic signals. Physica. D 158, 95 – 126. Von Storch, H., Navarra, A., 1995. Analysis of Climate Variability. Applications of Statistical Techniques. Springer Verlag, Berlin, Germany. Weare, B.C., Nasstrom, J.S., 1982. Examples of extended empirical orthogonal function analysis. Monthly Weather Review 110, 481 – 485. Wijnberg, K.M., Terwindt, J.H.J., 1995. Extracting decadal morphological behaviour from high-resolution, long-term bathymetric surveys along the Holland coast using eigenfunction analysis. Marine Geology 126, 301 – 330. Winant, C.D., Inman, D.L., Nordstrom, C.E., 1975. Description of seasonal beach changes using empirical eigenfunctions. Journal of Geophysical Research 80 (15), 1979 – 1986.
Long-term shoreline response of a nontidal, barred coast Grzegorz Ro´yyn´ski* Institute of Hydroengineering, Polish Academy of Sciences, 7 Kos´cierska, 80-953 Gdan´sk, Poland Received 16 February 2004; received in revised form 22 July 2004; accepted 22 September 2004 Available online 18 November 2004
Abstract Long-term variations of shoreline positions along the southern Baltic coast were investigated using multichannel singular spectrum analysis (MSSA) to determine the most dominant long-term response patterns. The investigated beach is located at Lubiatowo on the Polish Coast and is mildly sloping with multiple bars. Data on coastal morphology have been collected at Lubiatowo including (1) bathymetric surveys since 1987 twice a year, and (2) beach topography surveys since 1983 every 4 weeks on the average, extending from the shoreline to the dune foot. Furthermore, several dedicated field campaigns have been carried out at Lubiatowo, as well as measurements of deep-water wave properties since 1998. MSSA was employed to the whole data set of shoreline position from all survey lines. In summary, three patterns emerged reproducing alongshore standing waves with different periods 7 to 8, 20+ and several decades. They represent long-term shoreline response, such that at some locations the longest wave is most predominant, at other locations the medium cycle predominates, whereas the shortest is the most prominent at yet other locations. However, all three can be detected at every location monitored, eliminating the confusion resulting from ordinary singular spectrum analysis (SSA) analysis, done previously for the same data set. D 2004 Elsevier B.V. All rights reserved. Keywords: Data-driven modeling; Shoreline evolution; Multichannel singular spectrum analysis; Lagged covariances; Reconstructed components
1. Introduction Coastal processes include complex, nonlinear phenomena and involve morphological features that are characterized by a wide range of spatial and temporal scales (De Vriend, 1991, Larson and Kraus, 1995). Despite substantial progress in the research of * Tel.: +48 58 552 2011; fax: +48 58 552 4211. E-mail address: [email protected]. 0378-3839/$ - see front matter D 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.coastaleng.2004.09.007
physical mechanisms that shape these morphodynamic processes, present knowledge is still too limited to properly understand, describe, and forecast the behavior of many coastal systems (Stive et al., 2002). Fortunately, high-quality data sets on coastal morphology become more easily available now, providing an opportunity for extracting valuable information on morphological properties by means of computerintensive analysis and data-driven modeling (Larson et al., 2003).
80
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
In analysis of morphological data, it is common practice to split the state vector (i.e., the morphological data set) into two subvectors: the signal vector which is described by a few characteristic patterns believed to represent the dynamics of the considered process, and the noise vector which contains all processes which seem to be unrelated to the signal vector. For example, linear techniques like empirical orthogonal functions (EOF) are optimal in representing the variance. Therefore, the EOF analysis has been widely used by researchers in coastal morphology to identify patterns in the data (e.g. Winant et al., 1975, Aubrey, 1979, Wijnberg and Terwindt, 1995, Larson et al., 1999). Another linear method, namely canonical correlation analysis (CCA) maximizes the correlation between two simultaneously observed vector time series (Barnett and Preisendorfer, 1987, Ro´z˙yn´ski, 2003). Currently, a lot of efforts is devoted to introduction of nonlinear techniques in order to establish better representation of key morphological processes that are essentially nonlinear. Such methods usually involve the combination of at least two signal processing techniques in which ordinary EOF method usually serves as a data compression tool. For example, Ruessink et al. (2004) used neural network based circular nonlinear principal component analysis to study similar bar behavior at coastal segments in the Netherlands, Japan and USA. Other less known methods, like principal oscillation or interaction patterns, satisfy certain dynamical constraints (Hasselmann, 1988, Ro´z˙yn´ski and Jansen, 2002), which might be identified with some a priori knowledge of a studied system. Ordinary singular spectrum analysis (SSA) was employed to a limited extent for investigating coastal morphology and its evolution. Southgate et al. (2003) utilized SSA to extract long-term trends from measurements of the variation in shoreline position at several beaches around the world with the overall purpose of detecting oscillatory behavior in the filtered signal (also, compare Stive et al., 2002). Ro´ z˙yn´ski et al. (2001) applied this method to determine characteristic patterns in the shoreline response at a sandy beach in the southern Baltic Sea seeking patterns related to forced or self-organized behavior. That study revealed three significant patterns. The most important always represented longterm shoreline change, but its period varied from
location to location, so finally periods of some 9, 16 and around 30 years were identified. The two other patterns were identified to represent deterministic chaos and thus were found to reflect self-organized shoreline behavior. The current study is a direct extension of that SSA analysis and is aimed at a more thorough analysis of the long-term behavior. In particular, the study is focused on understanding why ordinary SSA demonstrated different predominant cycles of long-term shoreline response at different points of the monitored area. Analysis of the whole data set, including its cross-covariance structure was thus intended to come up with a more comprehensive idea of that part of shoreline change. Long-term stability of the studied area provides excellent conditions for application of the MSSA method. First, the field site at Lubiatowo Beach in Poland is described together with the morphological data collected. Subsequently, description of SSA and MSSA is given. Then, MSSA is applied to the whole data set (27 survey lines analyzed jointly) on shoreline position from Lubiatowo and patterns of variation in time and space are established and discussed.
2. Field site and data collected The Coastal Research Station (CRS) of the Institute of Hydroengineering, Polish Academy of Sciences, is located at Lubiatowo on the Polish Coast facing the southern part of the Baltic Sea (see Fig. 1). The beach usually features multiple longshore bars and has a mild mean slope of m=tana=1–1.5% with a median grain size of D 50=0.22 mm. Since the resultant longterm energy flux is oblique to the shore, a predominant west to east littoral drift is produced. Multiple breakers are frequently observed in the surf zone. Long-term records indicate that for average storms the significant wave height outside the surf zone (in water depth h of about 20 m) usually reaches H s=2–2.5 m with a mean period of T=5–7 s. The incoming wave energy is dissipated as the waves propagate onshore, and for h=2–3 m the mean height typically equals – H =0.5–1 m with T=4 –5 s (Pruszak et al., 1999). Closer to the shoreline (h=1 m), the height of the waves during storms reduces to 0.3–0.5 m. Since the Baltic Sea is practically isolated from the Atlantic Ocean, there is no influence from tides. The average
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
81
Fig. 1. Location of Coastal Research Station at Lubiatowo.
surf scaling parameter, defined as e=[(H bx 2)/(gm 2)], where x=2p/T, g the acceleration of gravity, m the average beach slope, and H b the wave height at incipient breaking, in general falls between 1000 and 3000, indicating that the coast at Lubiatowo is a highly dissipative system. Due to the mildly sloping bottom, the values of e correspond to low p values offfi ffiffiffiffiffiffiffiffiffiffiffiffi the usual surf similarity parameter n ¼ m= H0 =L0, which falls below 0.1 for the quantities H 0 and L 0 representing deepwater wave characteristics. Long-term bathymetric surveys have shown that the beach at CRS Lubiatowo normally includes a system of four stable bars (Ro´z˙yn´ski et al., 1999). The first, innermost bar is situated about 120–170 m from the baseline with depth over crest of 1 m, the second 220–300 and 2 m, third 400–500 and 3–3.5 m, and the fourth, outermost bar 600–800 and 4–5 m. Occasionally, an ephemeral bar develops in the vicinity of the shoreline, typically located some 20–80 m from the
baseline. The bars do not migrate, but only oscillate around their average locations (Pruszak et al., 1997, 1999). This observation was confirmed by measurements of nearshore bathymetry commenced in the early 1960s, which gained geodetic compatibility since 1983, when the current baseline was established. The baseline itself is situated 10–50 m from shoreline on a dune crest. Very rarely, it can be destroyed by a series of extreme events. In such instances, exact baseline positions need to be reestablished by measurements referenced to a local state geodetic system. The surveys of the beach topography, from the shoreline to the dune foot, started in 1983. Although the measurements indicate long-term beach equilibrium, the shoreline is not a perfectly stable feature. It was observed that some parts of the shoreline advanced, while others receded or displayed no movement. However, after some time the direction
82
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
to obtain equal sampling interval in the MSSA, data points were derived for these gaps using linear interpolation, which showed little impact on longterm patterns though, which can be explained by the fact that the duration of gaps is much shorter than the identified long-term cycles. As previously pointed out, no net long-term erosion or accretion has been observed at CRS Lubiatowo, although the shoreline regularly moves back and forth over a decade or so (Pruszak, 1993). Thus, these measurements are highly suited for extracting temporal and spatial patterns of shoreline change for a beach where oscillations occur but without any net movement. The available data set covers the period between Sept. 1983 and Sept. 1999, which makes it unique in duration in a sense that there are only few shoreline observations, whose time span covers more than a decade (Dutch Jarkus data set or FRF Duck, NC, USA). The spatial extension of more than 2.5 km is also an important factor. Together they allow for comprehensive analysis in both time and space resulting in extraction of large shoreline standing waves having periods of many years.
of movement for some parts may reverse, whereas others may be stable or commence oscillating back and forth. Over long periods all these combinations of movements are possible for a particular shoreline stretch and the mechanisms governing these phenomena appear to be quite complex. The surveys are attached to a local geodetic system, so all measurements taken at different times and along different lines are fully consistent. In addition, shoreline positions at the profile lines are all related to the mean sea level. In total, the surveyed area covered 2600 m of alongshore distance with 27 cross-shore profiles uniformly spaced at 100 m. Enumeration of these profiles originates from the history of the CRS Lubiatowo geodetic base (cf. Fig. 2). Initially, eight baseline points east of point 3 were fixed (points 3, 4, . . ., 10), after which the baseline was extended westwards to fix points 11, 12,. . . , 29. Thus, the baseline consists of two segments that meet at fix point 3 and because of the dune configuration their azimuths differ by 4811V38.4U. Time series of shoreline positions measured at all 27 survey lines were employed in the MSSA. The analyzed data set covers 16 years of measurements. Surveys started on 22nd Sept. 1983 and typically have been repeated every 4 weeks (28 days), and the last survey included in this study was carried out on 2nd Sept. 1999. However, the data set includes two longer time periods for which no measurements were taken: there is one gap in the data between 15th Dec. 1992 and 25th Nov. 1993, and another gap between 17th Aug. 1994 and 3rd Jan. 1996. In order
3. Ordinary and multichannel singular spectrum analysis Ordinary singular spectrum analysis (SSA) can be regarded as a particular application of EOF analysis where the data matrix contains values measured at a
N
Fig. 2. Geodetic base at CRS Lubiatowo.
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
specific location lagged in time (Vautard et al., 1992). Therefore, the column vectors of the data matrix contain the measured time series successively delayed in time up to the maximum shift known as the embedding dimension. The eigenvectors of the corresponding covariance matrix yield the dominant patterns in the time series under study. In contrast to traditional EOF analysis, SSA provides information on the dynamics of the underlying system (Von Storch and Navarra, 1995). Actually, it can be shown that the complete dynamics of the underlying system can be reconstructed from the measured data without any knowledge of the evolution equations (Broomhead and King, 1986). The SSA method is based on the Karhunen-Loe´ve expansion theorem, where the lagged covariance matrix of the analyzed time series x i , where 1ViVN, N—series length, is utilized in the expansion. This matrix is often called the Toeplitz matrix as it has the Toeplitz structure with constant diagonals corresponding to equal lags, 2 3 c ð 0Þ c ð 1Þ N c ð M 1Þ 6 c ð 1Þ 7 c ð 0Þ c ð 1Þ : 7 Tx ¼ 6 ð1Þ 4 5 : c ð 1Þ c ð 0Þ : c ð M 1Þ N c ð 1Þ c ð 0Þ where: cð j Þ ¼
N j
1 X xi x Þ xiþj x N j i¼1
ð2Þ
and –x is the mean value of the analyzed time series. The number M is the user defined embedding dimension and its selection is crucial for the analysis results. The larger M is, the better the spectral resolution of oscillatory components in the time series. However, to eliminate statistical errors in the autocovariance function for large j, it is recommended that M should not be greater than (1/3)N with the spectral resolution M 1. The Toeplitz matrix is symmetric, with positive terms on the main diagonal, so all its eigenvalues k k are positive, unless the data are perfectly noise-free and originate from a dynamical system with purely quasi-periodic, deterministic behavior. In such a case, some eigenvalues equal zero and the matrix is singular. The corresponding eigenvectors E k are
83
orthogonal, so the sum of products of the corresponding terms of the (column) eigenvectors E j and E i is equal to the Kronecker delta; M X
Ekj Ekl ¼ djl ; 1VjV M ; 1V lV M
ð3Þ
k¼1
Interestingly, the (column) eigenvectors of the Toeplitz matrix are either symmetric or antisymmetric with respect to (1/2)Mth term in each column. Eq. (3) denotes summation along columns. The same holds true for the sum of products along the jth and lth rows: M X
Ekj Ekl ¼ djl ; 1V jV M ; 1V lV M
ð3aÞ
k¼1
In addition, the spectral decomposition formula is valid as well: M X
kk Ekj Ekl ¼ Tx; jl ¼ cðjj ljÞ; 1VjVM ; 1VlVM
k¼1
ð4Þ The eigenvectors form the time-invariant part of the SSA decomposition, whereas the variability of a given system is contained in the principal components (PCs). The kth PC is a projection coefficient of the original signal onto the kth eigenvector: aki ¼
M X
Xiþj1 Ekj ; 1V iV N M þ 1
ð5Þ
j¼1
According to this formula, we have to take M elements of the original x series from the ith to i+Mth element, compute their products with the corresponding elements of the kth (column) eigenvector and sum them up to obtain the ith element of the kth PC. Hence, the PCs are time series of length NM. Even though the PCs are orthogonal to each other, it does not mean they are independent, because Eq. (5) shows that M consecutive elements of the original time series are needed to compute one term of every PC. Therefore, there are k common elements of this series for the ith term of the rth PC ari and the jth term of the sth PC a sj , such that k=M|ji|N0 (lag|ji|). Hence, the correlation structure of the original series must be imprinted in the sequence of the PC terms, producing nonzero
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
84
correlations for nonzero lags. For natural systems, where the assumption of normally distributed data generated by the system is widely accepted, it means that time series of such data are dependent on each other. The PCs in the SSA decomposition are not convenient, since they do not provide a unique expansion of the signal into a sum of different components. Using the PCs, we can expand our series as: xiþj ¼
M X
aki Ejk ; 1V j V M
ð6Þ
k¼1
There may be up to M subsets of the original time series containing the specific element x i+j . Thus, there are up to M different ways of reconstructing the components of the signal with Eq. (6), which in general do not produce the same results. Furthermore, instead of a series of length N, a series having NM+1 elements is obtained, so it is difficult to apply Eq. (6) for filtering and prediction. It is possible though to construct a least-square optimum series corresponding to a given (sub)set of the eigenelements. If the following expression is minimized, !2 N M M Xþ1 X X HW ð yÞ ¼ Yiþj aki Ekj ð7Þ i¼1
j¼1
kaW
where y is the desired series and W is the subset of eigenelements, the solution y=R W x is given by: ðRW xÞi ¼
M X 1 X ak Ek M j¼1 kaW ijþ1 j
for M V iV N
M þ1
ð8aÞ
ðRW xÞi ¼
i X 1X ak Ek for 1 V iV M 1 i j¼1 kaW ijþ1 j
ð8bÞ
ðRW xÞi ¼
1 N iþ1
M X
X
akijþ1 Ekj
j¼1N þM kaW
for N M þ 2ViVN
ð8cÞ
If W contains only a single index k, the resulting series is called the kth reconstructed component (RC), which will be denoted by x k , the RCs are thus
P k additive, i.e., RW x ¼ kaW x . Hence, the series x can be uniquely expanded as the sum of its RCs: x¼
M X
xk
ð9Þ
k¼1
To explain Eqs. (8a) (8b) (8c) in more detail, let us trace the computations of the ith element of the 1st RC. First, all combinations of the terms of the 1st PC (a1j ) and the 1st eigenvector (E k1) such that j+k=i+1 are multiplied, then their sum is taken and finally its mean value provides the desired term of the 1st RC: other RCs require similar computations on the underlying PCs and eigenvectors. The transformation of the original series x into x k is nonlinear, for the eigenvectors E k depend on x in a nonlinear fashion. Furthermore, the RCs are correlated even at lag 0, so additional correlation analysis is required for the extracted RCs. Consequently, variances of RCs are not cumulative. In some cases measurement time series from more than one location are analyzed simultaneously leading to MSSA. MSSA is mathematically equivalent to extended EOF (EEOF) analysis (Weare and Nasstrom, 1982), although EEOF is used for data sets with many points in space and only few lags in time. Apart from autocovariances of all time series included in an MSSA analysis, their cross-covariances are taken into account as well. It is schematically described by expression (10) that presents a block MSSA system matrix made up from covariances on the main diagonal and cross-covariances below and above. The Ti,i elements represent the autocovariance structure of ith time series and their structure is identical to Tx from expression (1). Consequently, the Ti,j terms describe cross-covariances between ith and jth time series from lag 0 up to the embedding dimension M. If the number of time series in an MSSA study is equal to L, the unfolded block matrix has ML ML terms and the embedding dimension of the whole problem is ML. Like in the case of ordinary SSA analysis, the MSSA matrix is symmetric, cf. Plaut and Vautard (1994), so both eigenstructures are formally identical to each other. Thus, using the rationale expressed through Eqs. (8a) (8b) (8c), we can calculate reconstructed components of an MSSA study. We should remember though that unlike in the ordinary SSA analysis, where the analyzed time series is
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
85
resolved into M reconstructed components, being time series as well, the MSSA reconstructed components are random fields. Thus, as a result of MSSA analysis we obtain ML random fields. This may make its interpretation quite complicated, since the initial field having NL elements is decomposed into MNL numbers, which can rapidly turn into a huge value. 2 3 T1;1 T1;2 : : : T1;L 6 T2;1 T2;2 : : : T2;L 7 7 T¼6 ð10Þ 4::: ::: ::: ::: 5 TL;1 TL;2 : : : TL;L
4. Analysis of shoreline change In order to analyze the temporal and spatial variation in shoreline position with the MSSA method, the measurements from all 27 survey lines at CRS Lubiatowo were employed. All further references to a dlineT or dprofile lineT will mean shoreline positions at that line. The time series from each line encompassed 169 individual measurements in addition to 28 values obtained through linear interpolation for the gaps where data were missing. Table 1 illustrates the mean value and standard deviation for all the studied profile lines, whereas Fig. 3 displays the entire data set, also indicating the
Table 1 Mean value and standard deviation (S.D.) of the analyzed time series of shoreline position at respective profile lines Profile number
Mean value (m)
Standard deviation (m)
Profile number
Mean value (m)
Standard deviation (m)
29 28 27 26 25 24 23 22 21 20 19 18 17 16
69 59.1 49.5 38.7 29.6 21.9 14.9 9.9 7.6 4.6 3.9 3.8 4.7 6.9
19.8 19.5 17.2 14.2 13.3 14 14 13.7 11.1 12.4 14.7 13.6 14.2 13.5
15 14 13 12 11 3 4 5 6 7 8 9 10
10.3 15.4 19.4 22.3 22.5 32.3 40.6 48.2 57.5 67 73.1 77.2 80.7
13.5 12.1 11.7 12.3 12.4 11.9 11.5 9.8 10.6 13.3 13.5 12.6 13.1
Fig. 3. Envelope of shoreline records at CRS Lubiatowo, 1983–1999.
envelope that contains the variation in shoreline position between 1983 and 1999 for all lines. Table 1 demonstrates that the overall behavior of the shoreline quantified in terms of standard deviation is quite similar, except for the two westernmost profiles 29 and 28 (partly 27 also), which show a larger scatter. Also, Fig. 3 shows that the variability is more or less the same all over the study area with a similar difference between the maximum and minimum shoreline position for all lines. The first step in the MSSA encompassed the choice of the maximum time lag M. The value of M=30 was selected, which is considerably less than the maximum length possible. This value is equal to the time lag in the ordinary SSA study, cf. Ro´z˙yn´ski et al. (2001), and was selected mostly to reduce the amount of computations, knowing from that analysis this choice includes all long-term behavior. It should be stressed here that the MSSA analysis would have been hardly possible, had the ordinary SSA study not been executed before. Hence, it highlights a synergistic character of ordinary SSA and MSSA. Fig. 4 presents the 1st RC of the MSSA analysis for all profiles. For lines 29–25 in the westernmost segment of the monitored area the 1st RC is very similar to the results of individual SSA analyses of these lines and varies in aF20 m range about the mean empirical shoreline position. For other pro-
86
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
Fig. 4. Shoreline standing wave with T=several decades and L=1500 m.
files, the 1st RC is not as spectacular and the presence of this pattern could not be found with the SSA method. Nevertheless, the range of variations in that portion of the studied coastal segment falls between +4 and 12 m. The variations have a character of a long standing wave; the lines for profiles 21 and 07 can be regarded as the nodes of that wave, because of their relatively small overall variations. The distance between these two profiles L=1500 m can be interpreted as the standing wavelength. Importantly, there is also a temporal node for all profiles (day 3800, i.e., year 1994). The period can be deduced from the trajectories of profiles 29 and 28 confronted with the lines 20 and
19. They represent the switch from maximum erosion to accretion of the 1st RC (lines 29 and 28) and vice versa (lines 20 and 19). As it occurred during 16 years (1983–1999), we can say the period of this standing wave is about 32 years. However, since this reasoning is entirely based on visual assessment of the 1st RC that embraces only 16 years of observations, it is safer to judge the period of this standing wave is several decades. The explanation of what drives this wave is not possible, because the time span of bathymetric surveys of nearshore bed topography was too short to find a pattern in nearshore bathymetry with a similar period. Nevertheless, it indicates that such a pattern
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
might be detected in the future when more observations of seabed in the nearshore region is available, which points to the direction of future research. Fig. 5 demonstrates trajectories of the 2nd MSSA reconstructed components. They can be interpreted as standing waves with a period of 7–8 years, with nodes at profiles 27–26 and 16–15. The antinodes can be observed at profiles 22–19 and 03–06. Therefore, a wavelength of 1000 to 1400 m can be identified from the locations of nodes and antinodes. The amplitudes are of similar magnitudes (F10 m), from their
87
extremes the period of 7–8 years could be established. Interestingly, this period coincides with the period of the climatic North Atlantic Oscillation (NAO), so it is possible that its influence is imprinted in shoreline variations of an isolated Baltic Sea. What is important, similar pattern was found for a number of profile lines with the SSA routine, where it was the most dominant feature. Like for the 1st RCs, the MSSA showed it appears all over the studied coastal segment as a standing wave. The origin of this wave and its relation to the patterns describing the evolution of nearshore bed
Fig. 5. Shoreline standing wave with T=7–8 years and L=1200–1400 m.
88
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
Fig. 6. Shoreline standing wave with T=20–22 years and L=1400–1600 m.
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
topography is unclear for the time being. Similar period was found with a traditional EOF routine evaluating jointly cross-shore bathymetric patterns at lines 4, 5, 6 and 7, cf. Ro´z˙yn´ski (2003). This pattern describes oscillations of two inner bars that bounce between the shoreline and outer bars. Unfortunately, a coupling of such oscillations with the alongshore standing wave cannot be established conclusively, despite their strikingly similar periods, because separate EOF analyses of individual cross-shore profiles were unable to trace inner bars oscillations, which if related with the shoreline standing wave, should be phase shifted from profile to profile. Nevertheless, the resemblance of periods of this shoreline standing wave and the NAO shows the avenue of future research, which should be aimed at better identification of the influence of large climatic oscillation on apparently isolated Baltic Sea. The 3rd MSSA pattern is presented in Fig. 6. Upon its inspection, this pattern demonstrates yet another standing wave with nodes at profiles 23, 22 and 21 and then at 05 and 06. The antinodes can be observed at profiles 16 and 17 and then at 08, 09 and 10. The temporal nodes are visible at day 1700 (1988) and 5500 (1999). From these trajectories, we can estimate the wavelength as to be between 1400 and 1600 m and the period to be between 20 and 22 years. The magnitude of amplitudes of this standing wave is similar to the amplitudes of the previous two and equals 14 to +10 m at lines in the center of the study area, where it is best pronounced, to F6 m at other lines. Like previously, this pattern was identified with the SSA routine at profiles where it was the most conspicuous (profiles 18–14). As it was with the 2nd MSSA pattern, it was not possible to directly relate this wave to changes in nearshore bathymetry. Being consistent with the results of SSA analyses, we can say that all three patterns represent the forced shoreline response to the wave climate, which consists of standing waves with different wavelengths and periods. The only difference is that in the SSA study it was assumed that long-term shoreline variations represent the response to changes in the wave climate. In light of the MSSA study, the wave climate in general can cause the shoreline to perform standing oscillations. Still, the effects of long-term wave climate fluctuations can be
89
traced as well (2nd pattern in particular). This is a real progress in the understanding of long-term shoreline evolution. Interestingly, locations of nodes and antinodes of these waves appear not to follow any consistent pattern and the evolution at particular spatial points is usually dominated by the most pronounced standing wave at these points, giving the impression that long-term shoreline evolution can vary substantially from one subsegment to another. The MSSA analysis found that in fact the shoreline evolution contains three periodic components at all points. The presented study revealed initially unexpected layout of long-term shoreline response (standing waves), which was found in the described three reconstructed components of shoreline change. The analysis of trajectories of higher RCs detected chaotic behavior, already known from conventional SSA analysis. Their interpretation though brought no new discoveries in relation to that study. For this reason, only three top MSSA reconstructed components were analyzed in detail.
5. Conclusions Multichannel Singular Spectrum Analysis (MSSA) was applied successfully to analyze shoreline evolution at a beach that displays no net long-term erosion or accretion. The MSSA found that the major response consists of three standing wave components with the periods of several decades, 20–22 and 7–8 years for all spatial points. Usually, one of them is much more pronounced than the remaining two at a given point, so conventional SSA routine was unable to extract them. Consequently, it was assumed in the SSA study that the forced component originates from wave climate fluctuations, whereas the MSSA results tend to conclude standing wave components reflect the wave climate in general as well. The present analysis clearly manifests the complementarity of SSA and MSSA methods, illustrated by the fact that interpretation of MSSA patterns is usually very difficult without extra a priori knowledge. This is relatively obvious, when we recall expression (10) containing covariances and cross-covariances among the analyzed time series. The SSA method can
90
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
provide such a priori knowledge by highlighting covariance structures of individual time series. Then, the analysis of the cross-covariance structure of the whole random field in question becomes feasible and more accurate. No relationship was found between the longest standing wave component (L=1500 m, T=several decades) and the evolution of nearshore bathymetry. To achieve that records of nearshore bed topography must be continued until they cover a similar time period. Still, the existence of the shoreline pattern lasting more than three decades indicates we can expect seabed patterns having comparable periods. Consequently, we should try to single out potential drivers of such long-term changes (climatic variations, long-term developments in the Atlantic Ocean, etc.). Thus, the discovery of such a pattern points to where future research should be focused. The coincidence of the period of most frequently encountered 2nd standing wave component (T=7– 8years) with the period of North Atlantic Oscillation (NAO) suggests NAO effects can be imprinted in the shoreline evolution of an isolated sea, such as the Baltic Sea. This is perhaps the most important conclusion, because the existence of periodicity in shoreline change that matches the NAO period points to NAO as a (rather unexpected) driver of morphodynamic evolution in apparently isolated seas. On second thoughts though, it shows that the scale and intensity of NAO cannot be ignored when studying long-term beach processes in such seas. Therefore, possible linkage between NAO and nearshore morphodynamics in the Baltic Sea also points to future research strategies. No relationship was also found for the 3rd standing wave component (T=20–22 years) and nearshore bathymetry. The understanding of its generation, as well as the generation of two other shoreline standing waves, may only become feasible in future when records of nearshore bathymetry are confronted with large enough hydrodynamic records of deep and shallow water waves, nearshore currents and water levels. This should result in better understanding of the feedback between hydrodynamics and nearshore bed topography with multiple bars, which together produce long-term shoreline change in the form of standing waves with different periods.
Notation aki ith element of kth principal component of SSA decomposition c( j) autocovariance of a time series at lag j D 50 mean sediment grain size Ekj jth element of kth eigenvector of Toeplitzmatrix Tx g acceleration of gravity – H, H s, H b, H , H 0 wave heights: measured, significant, at break-ing, average, and deepwater, respectively H W ( y) a least-square minimized quantity to derive reconstructed components h water depth L number of time series in MSSA analysis L0 deepwater wavelength M embedding dimension (window length) m=tana mean slope N number of observations in a time series xi shoreline position at time i, where i=1,. . ., N x¯ average shoreline position xk kth reconstructed component Tx M M system matrix in SSA T ML ML system matrix in MSSA T wave period, x=2p/T y=R W x least-square optimal series corresponding to W d j,l Kronecker delta e surf scaling parameter n surf similarity parameter kk kth eigenvalue of Toeplitz matrix W subset of eigenelements used in derivation of reconstructed components
Acknowledgements The author is indebted to the Commission of the European Community in the framework of the Energy, Environment and Sustainable Development (EESD) Thematic Programme, for its support within the Project: bHuman Interaction with Large Scale Coastal Morphological EvolutionQ HUMOR-EVK3CT-2000-00037. The author is also indebted to Professor M.J.F. Stive and Professor B.G. Ruessink for their in-depth reviews that gave grounds for improvement of the final manuscript in order to meet the journal standards.
G. Ro´˙zyn´ski / Coastal Engineering 52 (2005) 79–91
References Aubrey, D.G., 1979. Seasonal patterns of onshore/offshore sediment movement. Journal of Geophysical Research 84 (C10), 6347 – 6354. Barnett, T.P., Preisendorfer, R., 1987. Origins and levels of monthly and seasonal forecast skill for United States surface air temperatures determined by canonical correlation analysis. Monthly Weather Review 115, 1825 – 1850. Broomhead, D.S., King, P.G., 1986. Extracting qualitative dynamics from experimental data. Physica. D 20, 218 – 236. De Vriend, H.J., 1991. Mathematical modeling and large-scale coastal behavior: Part 1. Physical processes. Journal of Hydraulic Research 29 (6), 727 – 740. Hasselmann, K., 1988. PIPs and POPs: the reduction of complex dynamical systems using principal interaction and oscillation patterns. Journal of Geophysical Research 93 (D9), 11015 – 11021. Larson, M., Kraus, N.C., 1995. Prediction of cross-shore sediment transport at different spatial and temporal scales. Marine Geology 126, 111 – 127. Larson, M., Hanson, H., Kraus, N.C., Newe, J., 1999. Short- and long-term responses of beach fills determined by EOF analysis. Journal of Waterway, Port, Coastal, and Ocean Engineering 125 (6), 285 – 293. Larson, M., Capobianco, M., Jansen, H., Ro´yyn´ski, G., Southgate, H., Stive, M., Wijnberg, K., 2003. Analysis and modeling of field data on coastal morphological evolution over yearly and decadal time scales: Part 1. Background and linear techniques. Journal of Coastal Research 19 (4), 760 – 775. Plaut, G., Vautard, R., 1994. Spells of low-frequency oscillations and weather regimes in the northern hemisphere. Journal of Atmospheric Sciences 51 (2), 210 – 236. Pruszak, Z., 1993. The analysis of beach profile changes using Dean’s method and empirical orthogonal functions. Coastal Engineering 19, 245 – 261. Pruszak, Z., Ro´yyn´ski, G., Zeidler, R.B., 1997. Statistical properties of multiple bars. Coastal Engineering 31, 263 – 280. Pruszak, Z., Ro´yyn´ski, G., Szmytkiewicz, M., Skaja, M., 1999. Quasi-seasonal morphological shore evolution response to
91
variable wave climate. Proc. Coastal Sediments ’99, Conference, ASCE 1081-1093. Ro´yyn´ski, G., 2003. Data-driven modeling of multiple longshore bars and their interactions. Coastal Engineering 48, 151 – 170. Ro´yyn´ski, G., Jansen, H., 2002. Modeling nearshore bed topography with principal oscillation patterns. Journal of Waterway, Port, Coastal, and Ocean Engineering 128 (5), 202 – 215. Ro´yyn´ski, G., Pruszak, Z., Okro´j, T., Zeidler, R.B., 1999. Depth of closure and seabed variability patterns. Proc. ICCE’98, ASCE 2926-2939. Ro´yyn´ski, G., Larson, M., Pruszak, Z., 2001. Forced and selforganized shoreline response for a beach in the southern Baltic Sea determined through singular spectrum analysis. Coastal Engineering 43, 41 – 58. Ruessink, B.G., van Enckevort, I.M.J., Kuriyama, Y., 2004. Nonlinear principal component analysis of nearshore bathymetry. Marine Geology 203, 185 – 197. Southgate, H.N., Wijnberg, K.M., Larson, M., Capobianco, M., Jansen, H., 2003. Analysis of field data of coastal morphological evolution over yearly and decadal time scales: Part 2. Non-linear techniques. Journal of Coastal Research 4 (4), 776 – 789. Stive, M., Aarninkhof, S.J.C., Hamm, L., Hanson, H., Larson, M., Wijnberg, K., Nicholls, R.J., Capobianco, M., 2002. Variability of shore and shoreline evolution. Coastal Engineering 47 (2), 211 – 235. Vautard, R., Yiou, P., Ghil, M., 1992. Singular spectrum analysis: a toolkit for short, noisy, chaotic signals. Physica. D 158, 95 – 126. Von Storch, H., Navarra, A., 1995. Analysis of Climate Variability. Applications of Statistical Techniques. Springer Verlag, Berlin, Germany. Weare, B.C., Nasstrom, J.S., 1982. Examples of extended empirical orthogonal function analysis. Monthly Weather Review 110, 481 – 485. Wijnberg, K.M., Terwindt, J.H.J., 1995. Extracting decadal morphological behaviour from high-resolution, long-term bathymetric surveys along the Holland coast using eigenfunction analysis. Marine Geology 126, 301 – 330. Winant, C.D., Inman, D.L., Nordstrom, C.E., 1975. Description of seasonal beach changes using empirical eigenfunctions. Journal of Geophysical Research 80 (15), 1979 – 1986.