Data-driven prediction and analysis of chaotic origami dynamics

Advances in machine learning have revolutionized capabilities in applications ranging from natural language processing to marketing to health care. Here, we demonstrate the efficacy of machine learning in predicting chaotic behavior in complex nonlinear mechanical systems. Specifically, we use quasi-recurrent neural networks to predict extremely chaotic time series data obtained from multistable origami systems. Additionally, while machine learning is often viewed as a"black box", in this study we conduct hidden layer analysis to understand how the neural network can process not only periodic, but also chaotic data in an accurate manner. Also, our approach shows its effectiveness in characterizing and predicting chaotic dynamics in a noisy environment of vibrations without relying on a mathematical model of origami systems. Therefore, our method is fully data-driven and has the potential to be used for complex scenarios, such as the nonlinear dynamics of thin-walled structures and biological membrane systems.

Since chaos is generally defined for deterministic systems as a high sensitivity to initial conditions, it could be considered predictable in a mathematical sense, assuming all relevant information about the system is known. In practical terms, however, it is extremely difficult to accomplish this task hindered by unknown factors such as noise and interactions with the surrounding environment.
Therefore, it remains a formidable challenge to predict chaotic behavior in practice.
Here, we study a data-driven approach to analyze chaotic time series data and predict future response by using machine learning techniques. Recurrent neural networks (RNNs) constitute a powerful machine learning approach for processing and predicting time-series data [5,6] (see Fig. S1a for schematic illustration of a standard RNN). Due to such capabilities, RNNs have been recently employed for physics problems [7][8][9][10][11]. However, RNNs are usually considered to be a "black box" for learning and predicting time series. Thus, interpretation of the neurons' processing a time series, especially for chaotic data, has remained elusive. This is partly because the activation for each neurons in an RNN for the current time step depends on the activation of every other neuron at the previous time step (denoted by a red circle in Fig. S1a). This deep level of coupling interaction makes it challenging to extract meaningful information about the effects of individual neurons.
Recently, quasi-recurrent neural networks (QRNNs) have been developed, particularly for natural language analysis (see Supplementary Note 1 and Fig. S1b) [12]. They exhibit faster processing of time series data and competitive performance compared with other RNNs. Most notably, QRNNs allow access to hidden layers, offering the ability to scope the data set and assign "meaning" to each neuron. The QRNN is composed of convolutional layers which process time-series data in parallel across each time step, and pooling layers in which recurrent relations can be implemented (Fig. S1b). Due to the element-wise calculation in the pooling function, activation of each neuron does not depend on the past outputs of other neurons.
In this study, we demonstrate not only prediction of both periodic and chaotic data, but also analysis of hidden units' distinctive responses to such dynamic conditions by using the QRNN.
To examine effectiveness of our approach, especially in the experimental context, we need to produce unique sets of dynamic data containing periodic, subharmonic, and chaotic trajectories in a controllable manner. One of the examples showing such behavior is a structure with intrinsic bistability [13][14][15]. In this study, we design and fabricate a versatile bistable mechanical system based on origami unit cells, specifically triangulated cylindrical origami (TCO) [16] as shown in Fig. 1a. These TCO cells, in serial connection, can provide highly tunable properties [17][18][19] as well as multi-degree-of-freedom nature, thus providing an ideal playground to examine the effectiveness of our data-driven prediction.
Based on the experimentally measured time-series data, we explore the feasibility of extracting meaningful system information from our data-driven approach. Prior research also attempted to find parameters/functions of the governing equations by using data-driven approach [20][21][22][23].
However, our study differs from these previous work in that no mathematical model of the system or even the nature of dynamics (e.g., definition of chaos) is provided. In our analysis of the hidden layers of the QRNN, we find that the neuron activation of specific hidden units behaves sensitively to the input data, which naturally enables the QRNN to distinguish between chaotic and nonchaotic origami dynamics. This approach has great potential for predicting complex dynamics of structures by allowing access to evolutionary steps and providing more parameters to reach reliable decisions. 3

Experimental demonstration of chaotic behavior of the Triangulated Cylindrical Origami
Origami has been extensively studied recently due to its tailorable static responses, and we show that a unique TCO-based platform can produce rich data sets from its complex dynamics, especially chaos, thus enabling to examine the effectiveness of our data-driven approach. One of the interesting features of the TCO is that its axial and rotational motions are coupled with each other (see Fig. 1a where the folding sequence of the TCO unit cell is depicted, and also Supplementary Movie 1 for folding animation). Figure 1b (Left) shows the flat sheet with crease patterns of the TCO and its folded shape. To describe the initial shape of the TCO unit cell, we define the initial height (h 0 ), initial rotational angle (θ 0 ), and radius of the cross-section (R) as shown in Fig. 1b (Center-Right).
To analyze the potential energy, we fabricate prototypes by using construction paper sheets cut by a laser cutting machine (see the Method; Supplementary Movie 2). Figure 1c shows our fabricated prototype of the TCO unit cell with (h 0 , θ 0 , R) = (50 mm, 70 • , 36 mm). We conduct quasi-static cyclic loading tests on this paper prototype (see Fig. S2) to examine and enhance the repeatability of folding/unfolding behavior. After this cyclic loading, we extract the forcedisplacement relationship from compression tests [18] and obtain the energy curve as shown in Fig. 1d. Here, the energy is normalized by the initial height (h 0 ) and stiffness (K) at the initial unstretched state, i.e., δ = 0 (Fig. S2c). We observe the bistable behavior, such that the TCO unit cell possesses two local minima in its energy landscape where the distance between these two stable states (L b ) is 0.35. This characteristic distance will be used to aid in analyzing the dynamic behavior of the TCO.
Forced dynamic tests of a system of two connected TCO unit cells were used to create the chaotic response data sets. The unit cells had properties (h 0 , θ 0 , R) = (50 mm, ±70 • ,36 mm), with 4 the left-most cross-section attached to a shaker that generated harmonic excitation (see Fig. 2a; the Method). The experiment was conducted for different excitation frequency ranging from 5 Hz to 25 Hz. The folding motion of these two TCO unit cells is measured by two action cameras together with a customized digital image correlation program. In the experiments, we measured displacement (u i ) and rotational angle (ϕ i ) of each cross-section (i = 0, 1, 2, see Fig. 2b for the notation). We use these measured data as well as the velocities,u i andφ i , numerically calculated from u i (t) and ϕ i (t). This measured data is separated into two data sets, the first of which is used for training, and the second for evaluation (Fig. 2b).
A schematic diagram of the dynamic folding behavior is shown in Fig. 2c overlaid on the underlying double-well potential energy landscape. We can define three different regimes [15]; Intrawell, Interwell (Periodic), and Interwell (Chaotic) vibrations as shown in Fig. 2c. The intrawell oscillation means that the system exhibits small oscillations about one of the two local potential minima. If the system overcomes the energy barrier and goes to the other stable state, we observe the interwell vibrations, which is typically either periodic or chaotic. Quasi-periodic responses may also occur in nonlinear systems, though that was not observed in the origami experiments.
Figure 2d (Top) shows the measurement results for an excitation frequency (f ex ) of 12 Hz. The displacement of the left-most section (u 0 ), which is attached to the shaker, shows input sinusoidal waves. However, u 1 indicates chaotic motion. Here we define δ i = u i−1 − u i and plot a phase plane for δ 1 in which blue dots are all measurement data and red dots represent Poincaré map [24] as shown in Fig. 2d (Bottom). The horizontal axis is normalized by the distance between two local minimum points (L b ). Therefore δ i /L b = 1 indicates that the TCO unit cell transits to the second stable regime. It should be noted, however, that it is possible for a forced dynamic response to transit over the unstable potential hilltop, and reverse course without ever reaching the stable potential minimum. Thus, this normalization provides only a nominal indicator of a completed 5 transit between the two stable equilibria. The experiment result for f ex = 12 Hz clearly shows that the vibrations take place not only around the first stable state, but also around the other energy local minimum state aperiodically, which corresponds to the chaotic interwell vibration. This manifests the capability of our TCO system to form chaotic dynamics (see Supplementary Movie 3 for experimental measurements at different excitation frequencies).

Prediction based on quasi-recurrent neural networks
Based on the experimentally measured data, we predict whether a response is chaotic or periodic by employing the QRNN technique [12] (see Supplementary Note 1 for the detail information about the QRNN). This prediction relies solely on the data obtained from the experiment, and therefore, prior knowledge about a mathematical model of the system is not required. To predict the chaotic/periodic folding motion of the TCO structure, we use the QRNN consisting of three hidden layers. Each layer is composed of 352 units. The input data X contains n = 12 components (u i ,u i , ϕ i ,φ i for three different cross-sections), and each component has T = 128 data points, i.e., time steps, which correspond to data length of 0.53 s given the action camera's sampling frequency of 240 fps. Using the input data of X ∈ R 128×12 from t 1 to t 128 , the QRNN predicts all 12 variables for next 32 time steps. The total duration of the measured data contains for the training results). Please note that we train one QRNN system by using all frequency cases, which enables to predict both periodic and chaotic cases, instead of training neural networks for a specific frequency and predict the dynamics at the corresponding frequency case. It is worth noting that, even with a well-characterized chaotic system (i.e., having access to the governing equations), the sensitivity to initial conditions (and numerical rounding errors) means that quantitative deviation is always expected. Hence, qualitative matching is a more realistic goal. This can be seen by the fact that the frequency spectrum obtained from time series of u 1 and u 2 shows qualitatively similar trend as shown in the center insets in Fig. 3a-d. Specifically, the increase of lower frequency components for chaotic responses is successfully captured. Additionally, even though the prediction from the QRNN deviated from the experimental results as the number of excitation cycle increases, the QRNN outputs surprisingly reasonable peak deformations of both TCO unit cells. For example, the given TCO configuration with θ 0 = +70 • follows the counter-clockwise rotation (ϕ > 0) under compression and vice versa (i.e., clockwise rotation under tenstion), which is also predicted by the QRNN (see Fig. S4 for the configuration space of the first and second TCO units as a function of u and ϕ).
The phase portraits, plotted as a function of displacement and velocity, are shown in the rightside insets in Fig. 3a-d Hz; left inset in Fig. 4c) and chaotic (f ex = 12 Hz; center inset) cases. In addition, in the case of f ex = 20 Hz (chaotic), we observe notably different behavior between the 71st and 313th hidden units in response to the same input data.
To analyze the neuron activation for different periodic/chaotic regimes, we compare FFT analysis on the experimentally measured data u 1 (t) (Fig. 4d) with the Lyapunov exponent calculation shows drastically different behaviors depending on the regimes, i.e., positive (negative) values for periodic (chaotic) cases. In addition, it is interesting that 313th hidden unit exhibits weaker activation only for the first two chaotic regimes. We confirm this unique response of hidden units from 10 different training runs (see Fig. S5 for different patterns of hidden unit responses, and also the Method for how we find such a hidden unit with unique response).
We observe that the QRNN allows not only the prediction of origami folding behavior, but also the classification of chaotic/periodic input data of the TCO systems by simply monitoring the neuron activation of the hidden units in the final hidden layer. The hidden layer analysis can be 9 useful to understand how the QRNN responses to the input data. This simple approach has great potential to characterize and determine chaotic vibrations, providing a reason behind the decision.
Note that the concept of chaos is not implemented specifically during training process. Thus, the QRNN itself acquires prediction capability in the process of the training based on the hidden units' distinctive responses to input data between periodic and chaos cases.

Discussion
In this study, we have demonstrated a data-driven approach to predict and analyze This uncertainty in the quality of each prototype significantly influences the repeatability and consistency of the folding/unfolding motion. To avoid this, each unit cell underwent 200 cycles with a controlled displacement from -3 mm (tension) to 15 mm (compression) at 6 mm/s as a preconditioning process (detailed can be found in Ref. [19]).

Dynamic test
We conducted the dynamic test on a chain composed of two TCO unit cells with (h 0 , θ 0 , R) = (50 mm, ±70 • ,36 mm). The left end of the two-TCO-unit system was connected to the shaker (LDS V406 M4-CE, Brüel & Kjaer). The shaker was excited to apply single-frequency harmonic excitation to the system. To track the folding/unfolding motion of the TCO unit cells during the dynamic testing, we used two action cameras (GoPro Hero4) and customized Python codes for non-contact digital image correlation. We captured images from these two cameras at 240 fps, and identified the coordinates of the spherical markers attached to the corner of the interfacial polygons (see Fig. 2a) based on the triangulation method [19].
Prediction based on quasi-recurrent neural networks The QRNN used for the prediction consists of three hidden layers, and the filter width (k) of 6 is used for each convolution layer (see Fig. 4a). Training was carried out on minibatches of 50 data sets by using the Adam optimization algorithm [25] with learning rate of 0.001 and decay rate of 0.8. To create each data set, we randomly chose excitation frequency out of 21 different frequencies (f ex = 5 Hz ∼ 25 Hz), and from the selected frequency data with 5600 time steps, we randomly selected input data with 128 time steps and next 32 times step data to calculate error of prediction. We implemented this QRNN in custom codes written in Python by using the machine learning library, TensorFlow, and we performed the training/prediction.
Analysis on the hidden units To compare the activation of the hidden units for different frequencies, we calculated the average value over the duration of the input data (0.53 s). Let c ave,n (f ex ) be the average activation of n th hidden unit in the final hidden layer for the excitation frequency of f ex , we constructed a vector composed of C ave,n (f ex ) for different frequencies as C ave,n = [c ave,n (5 Hz), c ave,n (6 Hz), · · · , c ave,n (25 Hz)], which is visualized as a bar graph as shown in Fig. 4g. To find unique response of a hidden unit, we examined every hidden units by calculating the difference between C ave,n and a desired response pattern. For example, if we look for a hidden unit showing negative activation for chaotic data, we set a vector C target = [c 5 , c 6 , · · · , c i , · · · , c 25 ] where Then, we selected the hidden unit which has the minimum value of the norm C ave,n − C target .
By changing C target , we can explore the different types of hidden unit activation patterns (see Fig.   S5).

Data Availability
Data supporting the findings of this study are available from the corresponding author on request.