Quasi-chaotic behaviors of narrow-band response of a non-deterministic resonant system: application to analysis of ship motion in irregular seas

This paper describes similarities between the narrow-band response of a resonant system excited by random inputs and low-dimensional chaos, with particular emphasis on the geometric characteristics of trajectories reconstructed in m-dimensional phase space from measured scalar time series data with a time-delay coordinate system. In this study, the time series data of ship roll angle in irregular waves were analyzed as an example of the narrow-band response of a resonant system. These time series data were measured by one of the authors in Tokyo Bay. The similarities between the narrow-band response of a resonant system excited by random inputs and low-dimensional chaos are verified by numerical simulation data.


Introduction
It is necessary to distinguish low-dimensional dynamics and randomness in measured time series. Because, it is necessary for avoiding the mistake of concluding non-chaotic data as chaotic data. If the decision is wrong, there is a possibility of wrong time series analysis or construction of a wrong model. Therefore, detecting determinism in a time series is very important, and numerous studies have focused on this problem [4,5,10,12,18,24,29,33,40,42,[45][46][47]51,54,55,60]. The purpose of the present research is closely related to this problem. Obtaining quantitative indicators from measured scalar time series data x(t) requires reconstructing the appropriate state vector x(t) in m-dimensional phase space with a time-delay coordinate system K. Ueno  where τ is the time delay and m is the embedding dimension. In this paper, the trajectory of the vector x(t), a function of time t, is called the reconstructed trajectory. The analytic methods used in this study are based on Takens' embedding theorem [52], as extended to both forced and stochastic systems by Stark et al. [48][49][50]. This study particularly focuses on the geometric characteristics of the reconstructed trajectories in m-dimensional phase space. In most oscillation phenomena, the system has at least one resonant frequency at which large amplitudes can be generated by small inputs. At other frequencies, transmission is reduced and, at very high frequencies, the effective mass may be so high that the output is not measurable. Furthermore, since the output spectrum is confined to a narrow band of frequencies in the vicinity of the resonant frequency, the response is narrow-banded and the typical time series of the output resembles a sine wave in which amplitudes and phases vary [39]. Namely, the existence of the resonant frequency of the system acts as a kind of filter (see Fig. 1).
In this study, the time series data of ship roll angle in irregular waves were analyzed as an example of a narrow-band response of a resonant system excited by random inputs. The roll motion is the oscillatory motion around the x-axis in Fig. 2.
There are a number of similarities between the narrow-band response of a resonant system excited by random inputs and low-dimensional chaos.
The first is that the correlation dimension of the time series estimated by the Grassberger-Procaccia algorithm [14,15] converges at a low embedding dimension m [56][57][58]. This convergence is consistent with the result obtained by Rapp et al. [43], who proved that filtered noise can mimic low-dimensional chaotic attractors when examined with the Grassberger-Procaccia algorithm alone. Our study verified the Rapp et al. result with real field measurement data. For the first time, Kawashima [25] estimated the correlation dimension for the time series data of ship roll angle in irregular waves.
The second similarity is that the parallelness degree of the adjacent trajectories reconstructed in m-dimensional phase space is high [59]. Methods that can distinguish deterministic time series and stochastic time series by estimating parallelness degrees of adjacent reconstructed trajectories were proposed by Kaplan and Glass [24], Wayland et al. [60], and Fujimoto et al. [10]. In this study we propose a simple method to evaluate parallelness degrees of the adjacent trajectories quantitatively. The method proposed in this paper has a key advantage of involving a simple calculation and can be used to discriminate determinism and stochastic properties in a time series.
The third similarity is that short-term prediction is possible [56][57][58]. The Jacobian matrix estimation method is used as a prediction method in this study. The algorithm of this method was proposed by Sano and Sawada [44] and Eckmann et al. [7] to measure the Lyapunov spectrum at first. Farmer and Sidorowitch [9] used this algorithm for prediction of chaotic time series data. For predicting a chaotic time series, in addition to the above, other methods with reconstructed trajectories have been proposed; these include the method of analogs by Lorenz [31], the simplex projection method by Sugihara and May [51], the Voronoi tesselation by Mees [33], the local optimal linear-reconstruction method by Jiménez et al. [22], the local fuzzy reconstruction method by Iokibe et al. [21], the improved method of analogs by Ikeguchi and Aihara [17], and the regularized local linear method by Kugiumtzis [26].
The fourth resemblance is that not only is the conventional correlation coefficient of the predicted value for one step ahead and measured value high but so also is the difference correlation coefficient. The difference correlation coefficient was defined by Ikeguchi and Aihara [18,19].
The time series data of ship roll angle were measured in Tokyo Bay by one of the authors. We used a 19GT boat (L pp = 16.55 [m], B = 4.50 [m] and D = 1.55 [m]) to measure roll angle. The eigen frequency of roll motion was 0.286 Hz, and the representative frequency of waves at the time of measurement was 0.233 Hz. To dislodge the noise of measurement, the time series data of ship roll angle were passed through a bandpass filter (0.05-0.5 Hz). The similarities between the narrow-band response of a resonant system excited by random inputs and low-dimensional chaos are verified by numerical simulation data.
Considerable research has been conducted on the topics of chaotic resonators and chaotic oscillators (for instance, Reference [2,6,8,11,13,16,20,27,28,32,34,35,37,38,41,53,61,62]). However, few studies have considered the narrow-band response of a non-deterministic resonant system, as in the case of our research. In other words, very little research has been done in the cases where a forced moment in the equation of motion is non-deterministic. Our research, on the other hand, discusses the similarities between the narrow-band response of a resonant system and low-dimensional chaos. Therefore, our research differs in content from prevailing studies, which dealt with chaotic resonators and chaotic oscillators.

Correlation dimension
2.1 Estimation of the correlation dimension by the Grassberger-Procaccia algorithm [14,15] The following correlation integral is used to evaluate the correlation dimension D 2 of the reconstructed trajectory in m-dimensional phase space in practice: where, N indicates the total number of vectors x(t i ) reconstructed in m-dimensional phase space, ε indicates the radius of the hypersphere whose center is x(t j ), (·) indicates the Heaviside function, and x(t i ) − x(t j ) is the distance between x(t i ) and x(t j ). In this paper, the distance is defined by the Euclidean norm. Next, we consider the number N . If it is too small, the result will be unreliable, whereas if it is too large, computations will take too long. An expedient is to draw M samples (x(t n( j) ), j = 1, 2, . . . , M) out of N and evaluate the following [36]: In the actual estimation, the correlation dimension D 2 is obtained from the slope of the regression line of the rectilinear part of log C(ε) as a function of log ε (see Fig. 4). Figure 3 shows an example of a reconstructed trajectory of the time series data of the ship roll angle in irregular waves. In this case, the embedding dimension m is three.

Estimated result of the correlation dimension
(The time series data of the ship roll angle in irregular waves is shown in Fig. 6a.) It is difficult to define 'irregular wave' and 'irregular seas' appropriately. However, these terms are used quite often in marine engineering and naval architecture. In the present research, waves that do not have a constant period and amplitude are treated as irregular waves. Seas are considered to be irregular when their condition gives rise to irregular waves. The sampling period of the time series data, δt, is 0.1 [s]. The estimated result of the correlation dimension of the time series data of the ship roll angle in irregular waves is shown in Fig. 4. Here, N = 8,000 and M = 4,000. Figure 4a shows the relation between ε and C(ε) obtained from (2.2) by a double logarithm. In Fig. 4a, the value of the slope of the regression line applied to the linear part is the correlation dimension D 2 . The relation between the embedding dimension m and the correlation dimension D 2 is depicted in Fig. 4b, which shows that the correlation dimension D 2 converges at low embedding dimension m. This result is the first similarity between the narrow-band response of a resonant system (in this case, time series data of ship roll angle in irregular waves) and low-dimensional chaos. In addition, this characteristic  that the correlation dimension converges at low embedding dimension is a necessary condition for making short-term prediction of time series possible. If the correlation dimension is small, this suggests that the motion may be described by a small number of variables, even when the original system has many degrees of freedom [36].  proposed by Kaplan and Glass [24], Wayland et al. [60], and Fujimoto et al. [10]. We propose a simple method to evaluate parallelness degrees of the adjacent trajectories quantitatively. In the method proposed here, a simple calculation is one of the advantages. This method can be used to discriminate determinism and stochastic properties in a time series.
Let N Data be the total number of the data points in the time series. Select randomly N 1 vectors from the trajectories reconstructed in m-dimensional phase space with a time-delay coordinate system. The distance between two position vectors is defined by the Euclidean norm in this space. Select the N 2 neighboring vectors of the vector x(t i,1 ) as the following condition is satisfied: After evolution of a time interval δt, these vectors will proceed from Define α i ,ᾱ, and β as follows: If the adjacent trajectories that are reconstructed in m-dimensional phase space are nearly parallel, the value of the inner product of the unit vectors, is nearly one. Moreover, if the parallelness degree of the adjacent trajectories that are reconstructed in m-dimensional phase space is high en masse, the value ofᾱ is nearly one, too. In this case, the value of the indicator β must be nearly zero.

Estimation result of indicator β
In this section, three kinds of time series data (see Fig. 6) are analyzed. The first is the time series of ship roll angle in irregular waves (see Fig. 6a). The second one is the x component of the Lorenz model [30] (see Fig. 6b). Under the initial condition (x(0) = −1.81, y(0) = 0.01, and z(0) = 0.01), the fourth-order Runge-Kutta method was used to solve the following equations numerically with a time step of 0.01 to obtain the numerical solutions of the Lorenz model (Eq. 3.6): This is the second similarity between the narrow-band response of a resonant system and low-dimensional chaos. Furthermore, this characteristic that parallelness degree of adjacent reconstructed trajectories is high (i.e., β is nearly zero) is a necessary condition for making short-term prediction of time series possible by using the deterministic method.
3.3 Distinction between the narrow-band response of a resonant system and deterministic chaos by phase-randomized surrogate data If the function obtained by Fourier transforming time series data x(t) is X ( f ), the power spectral density function P( f ) is defined by the following equation in this study: where f is the frequency. In the actual calculation, T = N Data δt, where N Data is the total number of the data points in the time series and δt is the time interval. The phase-randomized surrogate data [54] x s (t) of the original time series data x(t) has the following properties. The power spectrum density function of the surrogate data x s (t) is equal to that of the original data x(t). However, it is randomized so that the initial phase of each frequency component constituting the time series data is different from the original data. Stam et al. [47] point out problems seen in the case of analyzing phase-randomized surrogate data of time series data with strong periodicity. In this study, we focused on these problems. In the time series data of ship roll angle in irregular waves used in this study, all the amplitudes are less than 20 [deg]. Therefore, the influence of the nonlinear restoring moment term of the equation of ship roll motion can be ignored. The power spectral density functions P( f ) of each time series data are shown in Fig. 8.
Surrogate data of time series data of ship roll angle in irregular waves, which represent the narrow-band response of a resonant system, and those of the x component of the Lorenz model, which is used for deterministic chaos, are shown in Fig. 9. For creation of surrogate data for time series data of ship roll angle in irregular waves, 8,000 data points taken at a time interval δt = 0.1 [s] were used. For creation of surrogate data of the x component of the Lorenz model, 10,000 data points taken at a time interval δt = 0.01 were used. For the time series data of ship roll angle in irregular waves, which represent the narrow-band response of a resonant system, values of indicator β of the original data and those of surrogate data are almost equal (see Fig. 10). This implies that time series data of ship roll angle in irregular waves are stochastic and a linear superposition. In other words, the complexity of time series data of ship roll angle in irregular waves is independent of the initial phase. However, for the x component of the Lorenz model, values of the indicator β of surrogate data greatly differ from those of the original data (see Fig. 10). This implies that the Lorenz model gives deterministic chaos, which subtly depends on initial values. This result agrees with that of Wayland et al. [60].

Prediction method
The Jacobian matrix estimation method (JMEM) is used as a prediction method in this study. Its algorithm was first proposed by Sano and Sawada [44] and Eckmann et al.  Embedding dimension m β [7] to measure the Lyapunov spectrum. Farmer and Sidorowitch used this algorithm for prediction of chaotic time series data [9].
) denote the newest observed vector that reconstructs in m-dimensional phase space. The distance between two position vectors is defined by the Euclidean norm.
We select L neighboring vectors x ε (t p(i) )(i = 1, 2, . . . , L) of the newest observed vector x(t j ) from the past observed vectors. Namely, we select the neighboring vectors in the hypersphere whose radius ε is centered at the newest observed vector x(t j ). The following condition is satisfied as selecting neighboring vectors (see Fig. 11): 1, 2, . . . , L).
The Newest x(t j )

Fig. 11 Jacobian matrix estimation method
According to this condition (4.1), the vector x ε (t p (1) ) is nearest to the newest observed vector x(t j ) in these neighboring vectors. When the vector x ε (t p (1) ) is the source, the displacement vectors for other neighboring vectors are the following: After evolution of a time interval δt, the neighboring vectors will proceed from We consider the following map f, which transforms from y i to z i : Here, 0 is the null vector. The following is obtained by a Maclaurin expansion: . By ignoring terms that are greater than or equal to the second degree, the following approximate equation is obtained: Therefore, we can express this as Here, we use the least-square algorithm to estimate the Jacobian matrix A, which minimizes the summation of the squared error norm between z i and Ay i as follows: By denoting the (k, l) component of the Jacobian matrix A by a k,l (= ∂ f k ∂ y l ) and applying condition (4.9), m × m equations are obtained with ∂ S ∂a k,l = 0. Here, we use the following expression for A: where V and C are m × m matrices, v k,l and c k,l are the (k, l) components of matrices V and C, respectively, and A = C V −1 .
Using matrix A, the following expression for the prediction is obtained: where the first component of the vectorx(t j +δt) will be the prediction value. Predicted values of n steps ahead are obtained by repeating the above operation.

Prediction result
The prediction results of the time series data of ship roll angle in irregular waves are shown in Fig. 12. Predicted values for 30 steps ahead, with data of the past 5,000 [sec] [sec] [deg] [800] x(t) = a 1 x(t − δt) + a 2 x(t − 2δt) + · · · + a n y(t − nδt) + ν, (4.12) AIC(m) = N Data log(2πσ 2 n ) + N Data + 2(n + 1), (4.13) where n is the degree of the self-regression formula, δt = 0.1 [s] is the time interval, a i (i = 1, . . . , n) is the self-regression coefficient, and ν is white noise, which is independent from the past of x(t) with the condition that its average is zero and its variance is σ 2 n . N Data is the total number of data points in the time series. Figure 13 gives scatter diagrams corresponding to Fig. 12. In Fig. 13, r 1 indicates a conventional correlation coefficient and rmse indicates root-mean-square error. For time series data of ship roll angle in irregular waves, which represent the narrow-band response of the resonant system, it has been confirmed that prediction results by JMEM, which was devised to predict deterministic chaos, are more precise than those obtained by the linear AR model, which is the conventional method. This is the third similarity between the narrow-band response of a resonant system and low-dimensional chaos.
When we predict time series data of the narrow-band response of a resonant system using JMEM, we should be careful about the number of neighboring vectors, L. If it is too small, the prediction result will be unreliable. This is a difference between the narrow-band response of a resonant system and low-dimensional chaos. However, if it is too large, there is no remarkable difference between the prediction result using JMEM and that from using the linear AR model.

Influence of high-frequency noise
In this study, the time series data of ship roll angle were passed through a bandpass filter (0.05-0.5 Hz) to dislodge the trend and the high-frequency noise. It is difficult to judge only by seeing raw data whether measured time series data are polluted by minute high-frequency noise. The raw data x o (t) of the time series data of ship roll motion and the filtered data x(t) are shown in Fig. 14a. For examining influences of high-frequency noise, it is effective to amplify the high-frequency component by numerical differentiation [58]. The third derivatives of the raw data and of the filtered data are shown in Fig. 14b. The graph of the third derivative clearly shows that the raw data are polluted by high-frequency noise. Next, using the indicator β defined in and dx(t) dt ). Black squares correspond to the second derivative Sect. 3.1 of this paper, we examined the influence of high-frequency noise. For the raw data and the filtered data, first, second, and third derivatives were obtained. For each, values of the indicator β for embedding dimension m are shown in Fig. 15. The second and third derivatives of the raw data for which the high-frequency noise was amplified have randomness close to the white noise shown in Sect. 3.2 of this paper. The influence of the high-frequency noise clearly appeared even in the prediction for 10 steps ahead. The prediction results for the raw data and the filtered data by JMEM are shown in Fig. 16. Here, L = 250, m = 11 and τ = 0.7 [s]. For the raw data, predictability deteriorates with the influence of high-frequency noise. As above, it is necessary to pay attention to high-frequency noise in an analysis of time series data of the oscillatory system narrow-band response actually measured.

The difference correlation
Ikeguchi and Aihara [18,19] showed that the difference correlation can distinguish deterministic chaos from 1/ f α -type colored noise. The difference coefficient of correlation r 2 is the coefficient of correlation between the first-difference time series where x(t i ) andx(t i ) are actual and predicted time series, x and x are the averages, and N q is the number of data points in the time series x(t i ). x(t) rmse=0.00749 rmse=0.00750

Fig. 17
The conventional coefficient of correlation r 1 and the difference coefficient of correlation r 2 (one step ahead) Both chaos and 1/ f α -type colored noise are predictable over the short term and the conventional correlation coefficient r 1 is high. However, Ikeguchi and Aihara showed, using the difference correlation coefficient r 2 defined by (4.14), that chaos takes a high value, but 1/ f α -type colored noise takes a low value [18,19]. Figure 17 shows the conventional correlation coefficient r 1 and the difference correlation coefficient r 2 for time series data of the ship roll motion, which represent the narrow-band response of the resonant system. This is the result of forecasting one step ahead by JMEM. Here, L = 250, m = 11 and τ = 0.7 [s]. For the results for time series data of ship roll, which represent the narrow-band response of the oscillatory system, not only the conventional correlation coefficient r 1 but also the difference correlation coefficient r 2 take on high values, similar to chaos. This is the fourth resemblance of the narrow-band response of a resonant system and low-dimensional chaos.

Power spectral density function
Distinguishing low-dimensional chaos from the narrow-band response of the resonant system by checking for a peak in the power spectrum density function poses a potential problem. This is because a characteristic peak can appear when the power spectrum density function is obtained from a short part of the time series data in low-dimensional chaos. The power spectrum density function for 4 ≤ t ≤ 8 (see Fig. 6b) of the x component of the Lorenz model [30] is shown in Fig. 18.

Numerical simulation
In this section, we verify the similarities between the narrow-band response of a resonant system excited by random inputs and low-dimensional chaos by numerical simulation data. Under the given initial conditions, a fourth-order Runge-Kutta method was used to solve the following equation numerically at a time step of 0.1 to obtain the simulated data: where b 1 = 0.48, b 2 = 1.25, c 1 = 1.25, and p(t) is an external force. The total number of the data points in the time series was N Data = 8, 000. Figure 19a shows part of the time series data of p(t). Figure 20a shows the power spectrum density function of p(t). The features of the time series data of p(t) are shown in the graph of spectral density function (Fig. 21). The spectral density function of p(t) has the properties of Pierson-Moskowitz spectrum. Figure 19b shows part of the time series data of φ(t). Figure 20b shows the power spectrum density function of the simulated data φ(t).  Figure 21 shows the estimated results of the correlation dimension of the simulated data φ(t), indicating that the correlation dimension D 2 converges at the low embedding dimension m in this example too (see Fig. 21).
Values of indicator β that correspond to embedding dimension m are shown in Fig. 22. The values of the indicator β of simulated data φ(t), which represent the narrow-band response of a resonant system, are near the values of β of the x component of the Lorenz model, which is used for deterministic chaos, rather than the values of β of the white noise.
The prediction results of the simulated data φ(t) are shown in Fig. 23. The first 1,000 data points in the time series were excluded to guarantee that a steady state was reached in the prediction. Prediction values for 10 steps ahead, with the data of the past 5,000 points as reference value at the beginning, are shown. Here, L = 250, m = 8 and τ = 0.2. For the time series data of the narrow-band response of a resonant system like φ(t), it has been confirmed that the prediction results by JMEM, which was devised to predict deterministic chaos, are more precise than those obtained by the linear AR model, which is the conventional method.  Figure 24 shows the conventional correlation coefficient r 1 and the difference correlation coefficient r 2 for the simulated data φ(t), which represent the narrow-band response of the resonant system. For the simulated data φ(t), which represent the narrow-band response of the oscillatory system, both the conventional correlation coefficient r 1 and also the difference correlation coefficient r 2 take on high values, similar to chaos.

Conclusion
Similarities between the narrow-band response of a resonant system and low-dimensional chaos are verified by using the data obtained by measurement on the actual boat and that obtained by performing a simulation. The results of the analysis of the data obtained by simulation of the model expressed by Eq. (6.1) are similar to those obtained by the analysis of the data measured on the actual boat described in Sects. 2-4. Therefore, it can be concluded that the model expressed by Eq. (6.1) is valid. Further, it can also be concluded that the model expressed by Eq. (6.1) is an example of the narrow-band response of a resonant system excited by random external forces. Further, the graph of spectral density function shown in Fig. 20 also indicates that this model is an example of the narrow-band response of a resonant system.
The authors do not consider ship roll motion in irregular waves as a completely deterministic phenomenon. A number of similarities between ship roll motion in irregular waves and low-dimensional chaos have been mentioned. Thus, it has been indicated that for the measured time series data, it is difficult to differentiate the narrow-band response of a resonant system from chaos by using conventional methods like correlation dimension. This is because it may lead to a wrong conclusion that the narrow-band response of a resonant system is equivalent to chaos. In addition, techniques that can be used for differentiating between the narrow-band response of a resonant system and low-dimensional chaos are discussed in this paper. Evaluating phase-randomized surrogate data and the indicator β is effective for distinguishing low-dimensional chaos from the narrow-band response of a resonant system.
It is contemplated that the method of analysis employed in the present research is applicable not only for ship roll motion in irregular waves but also for the other narrow-band response of a resonant system.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited.