Supervised chaotic source separation by a tank of water

Whether listening to overlapping conversations in a crowded room or recording the simultaneous electrical activity of millions of neurons, the natural world abounds with sparse measurements of complex overlapping signals that arise from dynamical processes. While tools that separate mixed signals into linear sources have proven necessary and useful, the underlying equational forms of most natural signals are unknown and nonlinear. Hence, there is a need for a framework that is general enough to extract sources without knowledge of their generating equations and flexible enough to accommodate nonlinear, even chaotic, sources. Here, we provide such a framework, where the sources are chaotic trajectories from independently evolving dynamical systems. We consider the mixture signal as the sum of two chaotic trajectories and propose a supervised learning scheme that extracts the chaotic trajectories from their mixture. Specifically, we recruit a complex dynamical system as an intermediate processor that is constantly driven by the mixture. We then obtain the separated chaotic trajectories based on this intermediate system by training the proper output functions. To demonstrate the generalizability of this framework in silico, we employ a tank of water as the intermediate system and show its success in separating two-part mixtures of various chaotic trajectories. Finally, we relate the underlying mechanism of this method to the state-observer problem. This relation provides a quantitative theory that explains the performance of our method, and why separation is difficult when two source signals are trajectories from the same chaotic system.


I. INTRODUCTION
Blind source separation (BSS) is the separation of source signals from a mixed signal with little or no information regarding the source signals or mixing process. A classic example is the cocktail party problem, where a listener follows any one of many simultaneously occurring conversations at a cocktail party. BSS also has many Chaos ARTICLE scitation.org/journal/cha notable applications in digital signal processing, such as removing artifacts from electroencephalography (EEG) and magnetoencephalography (MEG) recordings. [1][2][3][4][5] When the mixed signal has a lower dimension than the total dimension of the sources, the BSS is called underdetermined. Many methods have been proposed to solve BSS in various scenarios. For example, by assuming various types of statistical independencies or mixing properties of the source signals, unsupervised classical methods such as principal component analysis (PCA), 6,7 independent component analysis (ICA), 8,9 and non-negative matrix factorization (NMF) 10,11 have been proposed. While these methods have dramatically enhanced our ability to parse data from linear and statistical distributions, it has been shown that adaptations of the above methods [12][13][14] as well as many other supervised learning methods, such as the Wiener filter, 15 support vector machines, 16 deep learning networks, 17,18 and recurrent neural networks, 19,20 outperform classical methods when the signals are generated from complex dynamical sources.
In this paper, we focus on a particular type of separation problem: chaotic source separation (CSS). Specifically, we consider the d-dimensional mixed signal to be a superposition of two d-dimensional trajectories, each of which is generated by an autonomous d-dimensional chaotic system. This problem is of particular relevance in the study of high-dimensional biological signals such as those from neural systems, as experimental measurements involve a mixture of electrical activity, correlated artifacts, and hemodynamic response. 2 Hence, it is of interest to study how one can extract a chaotic trajectory of interest from the mixed signal.
Here, we propose to solve this problem with an intermediate dynamical system that is trained by a supervised learning method. Although the dimension of the mixed signal d is only half of the total dimension 2d, the problem can still be solved by a supervised learning framework, where the exact separated trajectories are known during a training period. As a significant extension from previous studies 21,22 that require knowledge of the governing equations of the source chaotic systems, we build on a prior demonstration from the present authors that a recurrent neural network (RNN) can solve the CSS problem in the absence of these equations 19 (a more recent study by Krishnagopal et al. also demonstrated that a reservoir computer can solve the CSS problem 20 ). In this paper, we extend the demonstration by enacting this separation through a dynamically simple intermediate system, which is a simulated tank of water, and provide a quantitative theory explaining why and how such chaotic source separation is solvable. Our theory accurately predicts that separation is harder when the two source signals are generated by the same chaotic system and provides a foundation for the principled study of source separation in nonlinear dynamical data.

A. General scheme
We begin with a simple description of a general scheme of our supervised CSS (Fig. 1). We consider extracting trajectories of two autonomously evolving chaotic systems, s a (t) and s b (t), from their FIG. 1. Schematic plot of the source separation model. Two signals, s a (t) and s b (t), are generated independently by the two sources, which are two chaotic dynamical systems. The two signals are mixed into s + (t) = s a (t) + s b (t). The intermediate system evolves with s + as its driving signal. Supervised learning is applied to find proper readout functions φ a and φ b that extract the separated signals s a and s b from their mixture s + .
CSS is similar to an underdetermined BSS problem in the sense that the dimension of the mixture d is less than the total dimension of the sources, i.e., 2d. As a result, there exist mixed states s + that correspond to multiple distinct pairs of sources. Thus, without utilizing the temporal structure, one cannot find a function that maps the simultaneous state s + (t) to the separated states s a (t) and s b (t).
The essential idea of our method is to implement a highdimensional dynamical system as an intermediate system, which is continuously driven by s + (t). Since the state of the intermediate system, r, incorporates both the immediate value and the history of the mixed signal s + (t), one may obtain the full states of the two sources by training the output functions φ a (·) and φ b (·), as shown in Fig. 1. We assume that the governing equations of the source systems are unknown. However, different from the BSS problem, we do assume that the separated trajectories s a (t) and s b (t) are known for a finite time window. During this time window, we match the recorded state of the intermediate system r(t) with the two separated signals s a (t) and s b (t), and we look for two functions, φ a (·) and φ b (·) that can estimate the separated signals based on the state of where t is within the time window.

B. Intermediate system instantiated by a tank of water
To demonstrate the generalizability of this scheme beyond the RNN used in prior work, 19 we instantiate the intermediate system in silico as a tank of water (Fig. 2). We test the performance of this intermediate system on the CSS problem given the mixed signals of different pairwise sums from six distinct chaotic systems. We show their governing equations in Table I and their attractors in Fig. 3. We notice that trajectories of different chaotic systems have different ranges; to simplify the simulation and ensure an accurate quantification of the error, we deliberately preprocess the chaotic trajectories such that all variables have zero mean and unit variance along the time axis.
We construct the intermediate system as a square tank of water that is constantly perturbed by the mixed signal s + (t). The perturbed water evolves following the nonlinear partial differential equations, where h(x, y, t) is the height of the wave surface; u(x, y, t) and v(x, y, t) are the zonal and meridional speeds, respectively; g = 9.8 is the gravity constant; and b > 0 is the viscous drag coefficient. The tank has a flat 1 × 1 bottom and four vertical hard walls with a reflective boundary condition. When the perturbing term p = 0, Eqs. (3a)-(3c) become the traditional shallow water equations with the presence of a viscous dragging force. 23 Although other forms of perturbation exist, for example, a time-varying bottom, which requires the modification of all three equations, for the simplicity of the demonstration, we drive the water by artificially defining the perturbation term, on the right hand side of Eq. (3a) only. The perturbation term p(x, y, t) can be considered to reflect the speed with which one vertically and inhomogeneously adds or removes water from right above   Fig. 2. To guarantee the conservation of the water volume, we renormalize each filter such that where V = [0, 1] × [0, 1] ⊂ R 2 for each i = 1, 2, . . . , d. Thus, we now have a tank of water that is constantly being perturbed by the input signal s + while preserving its total volume. The simulation of this perturbed shallow water system is done by a modified Lax-Wendroff method. 24 The method preserves the second order spatial and temporal accuracy even with the presence of the three source terms in the partial differential equations (p, buh, and bvh). The viscous coefficient b is empirically set to 0.3. In this finite difference method, we discretize the wave surface into a 128 × 128 grid and integrate it with a time step t = 0.03.
Starting from the initial quiescent wave surface h(x, y, 0) = 1, we drive this dissipative water wave system by the d-dimensional mixed trajectory s + (t). After a transient period (T dump = 600) that is long enough to wash out the effect of the initial condition, we record the water's reaction to the mixed signal. To reduce the amount of data recorded, we sparsely measure the deviation of water elevation from the equilibrium height h = 1 at 2000 randomly selected locations, denoted as h(t) ∈ R 2000 . Then, with the available separated trajectories during the training period, we construct output functions that map h(t) into separated signals [s a (t), s b (t)]. Although many other forms of output functions may also work, we adopt the following nonlinear form with a tanh-type saturation on the quadratic and cubic terms: where W ∈ R 6×6000 is the coefficient matrix of the nonlinear output function. With the recorded h and the available s a and s b during the training phase (T tran = 600 with 20 000 time points), the output weight matrix W is calculated by the least squares method with the Tikhonov regularization, α = 0.001. We note that other output functions that outperform this one should exist. However, the purpose of this simulation is to demonstrate that CSS is indeed solvable by such an intermediate system, rather than to develop an optimal design. In Fig. 2, we show a schematic of the intermediate system separating a mixed signal, which is a summation of a Lorenz trajectory and a Rössler trajectory. Given the six distinct chaotic systems listed in Table I, we train and test the separation performance of a shallow water system driven by 6 2 − 6 * (6 − 1)/2 = 21 mixed signals. Each mixed signal is the sum of two chaotic trajectories, s a (t) and s b (t), that are from the ith and jth chaotic system, respectively, for 1 ≤ i ≤ j ≤ 6. To test the system's performance in separating each mixed signal, we reinitialize the water at quiescence and drive it with a new mixed signal. The separated signal s a from the post-training water system, after a transient period (T = 600), is plotted in Fig. 4. Specifically, the trajectory on row i and column j is the separated s a , where the mixed signal is a summation s + = s a + s b , with s a and s b   FIG. 4. Signal separation. Here, we show the signal s a separated by a posttraining water tank from the mixed signal s + = s a + s b , where s a and s b are the distinct trajectories generated by two of the six chaotic systems listed in Table I, respectively.
from system i and system j. For cases where i = j, we ensure that the two trajectories being mixed are distinct, i.e., s a (t) = s b (t), by choosing different initial conditions. By visually collating the separated signals in Fig. 4 and the six chaotic attractors in Fig. 3, we note that the performance of the separation varies across chaotic systems. Specifically, the Sprott N and Rössler trajectories separated from a mixture with other systems seem to have much higher quality compared with others (see the first and second rows in Fig. 4). To quantify performance, we calculate the mean squared error (MSE) between the actual trajectory and the separated one during the post-training period following T = 600 (see Fig. 5). We do not find a concrete relationship between the separation performances and the Lyapunov dimensions of the chaotic attractors. However, we note that the quality of the separation appears particularly poor when the two source signals are trajectories from the same chaotic system (see the diagonal line in Fig. 4). In Sec. III, we explain the underlying mechanism behind this supervised CSS and give an explanation for the diminished performance when signals are taken from the same attractor.

III. UNDERLYING MECHANISM OF SUPERVISED CHAOTIC SOURCE SEPARATION
We notice that the chaotic source separation (CSS) problem is essentially a nonlinear state-observer problem, and the intermediate system plays the role of the state-observer. To explicitly state this role, we rewrite the dynamical equations of the two source systems by combining them into a single 2d-dimensional autonomous dynamical system, denoted bẏ where is the concatenation of the two source systems. Thus, the simultaneous mixed signal s + (t), a summation of the two source systems, can be considered as an output from the combined system depicted in Eq. (7), i.e., y = g(x), (9) where g(x) ≡ s a + s b = s + ∈ R d . The CSS problem can then be considered as the problem of uncovering x in Eq. (7) through y in Eq. (9), which is the measurement of the combined system. More precisely, a state-observer is a system that estimates the full state variable of the combined system, x, by considering the measured y from the combined system. In our method, the state-observer is the intermediate system driven by y, i.e., the perturbed water tank. Such a dynamical state-observer solves the CSS problem by capitalizing on the power of invertible generalized synchronization. 19,25,26 Specifically, as the intermediate system is driven by the measured y, it evolves nonautonomously according to After a transient period, if the state-observer and the combined system exhibit invertible generalized synchronization, the state of the intermediate system is then uniquely determined by the concurrent state of the combined system x through an invertible map, i.e., r(t) = φ(x(t)). Thus, to estimate the full state of the combined system, one simply needs to train a readout function φ −1 (·) that approximates the inverse of the generalized synchronization function based on the state of the intermediate system r. Notice that invertible generalized synchronization is a property that emerges from the particular choices of both the intermediate system and the combined system. As such, we do not find a general principle of designing an intermediate system that guarantees invertible generalized synchronization with any combined chaotic system.
After elucidating the connection between the CSS problem and the state-observer problem, we emphasize that it is only when the full state is observable that such an invertible generalized synchronization function φ(·) can exist. In other words, the combined system [Eq. (7)] has to be observable through the output function [Eq. (9)]. The classical work of Kalman has discussed the observability of linear dynamical systems (see Ref. 27). In our case, however, the combined dynamical system is nonlinear and autonomous. The necessary and sufficient condition for such a combined system [Eq. (7)] to be observable through the measured output y is discussed by Inouye in Ref. 28. Specifically, the system is observable if and only if the observability mapping is univalent, 28 where the entries are defined as When both f ab (·) and g(·) are analytic functions on R d , the system [Eqs. (7) and (9)] is observable if and only if the equations G k (x) = G k (x ) with k = 1, 2, . . . imply only the trivial solution x = x . With the necessary and sufficient condition for observability, we can now investigate whether the CSS problem can be solved, providing the measured s + when the two source chaotic systems share the same dynamical equation, such that f a (·) = f b (·) = f(·). To answer this question, we define where s = s are distinct trajectories generated byṡ = f(s). We then rewrite Eq. (7) as d dt and rewrite Eq. (9) as By substituting Eqs. (13)-(15) into Eqs. (11) and (12), we observe that for any k, g k (x) = g k (x ) even though x = x , suggesting that the observability mapping is not univalent. Hence, we explain why the CSS performance along the diagonal in Figs. 4 and 5 tends to be worse compared to the off-diagonal counterparts in the same row.

IV. TESTING ROBUSTNESS AND GENERALIZABILITY
A. Robustness to noise To investigate how the CSS performance changes when the source signals are corrupted by observation noise, we modify the simulations in Sec. II. Specifically, we consider that the measured source signals are s a (t) = s a (t) + σ ξ a (t) and s b (t) = s b (t) + σ ξ b (t), and hence, the mixed signal is s + (t) = s a (t) + s b (t), where σ ≥ 0 is the noise strength and ξ a/b (t) are the white noise terms. By comparing Fig. 6 to Fig. 5, we find that the MSE does not significantly increase until the noise strength surpasses 0.1.

B. Generalizability to high-dimensional chaotic signals
Heretofore, we have only tested CSS on source signals that are three-dimensional. We now address the question of whether a tank of water can be trained to deal with high-dimensional chaotic signals. Accordingly, we employ the Kuramoto-Sivashinsky (KS) system and the Lorenz 96 system as the two chaotic source systems.
We obtain two 32-dimensional source time series by (i) integrating the standard Kuramoto-Sivashinsky equation, 29,30 in region 0 ≤ x < L = 22 (discretized into 32 evenly spaced grid points) with a periodic boundary condition and time resolution t = 1/16, and (ii) integrating the Lorenz 96 equations, 31 with time resolution t = 0.001 and a periodic boundary condition, where i = 1, 2, . . . , 32. As in the previous simulations, we preprocess the source signals such that each of their variables has mean zero and unit variance along the time axis. The mixed signal [ Fig. 7 For this simulation, we utilize a spatial discretization of 256 × 256, which is finer than the low-dimensional case. The viscous drag coefficient is set to b = 0.6. While these parameter choices are sufficient for this demonstration, further parameter optimization could lead to better performance for this or other systems. Complementing previous studies on source separation problems, [6][7][8][9][10][11][12][12][13][14][15][16][17][18][19][20] we show that separation of signals from a mixture of chaotic trajectories can be considered as a nonlinear stateobserver problem. With this realization, we propose to solve the problem by employing and training an intermediate system that is continuously driven by the mixed signal. We extend earlier studies where CSS is solved by recurrent neural networks, 19,20 and we show that even a tank of water under this proposed framework can solve the CSS problem. By making the connection between the CSS problem and the nonlinear state-observer problem, we explain the reason why separating two signals generated from the same chaotic system tends to be difficult.
We note that in this paper, we only consider mixed signals that are sums of two chaotic trajectories. Yet, our method can be applied to other mixing equations or mixtures of more than two chaotic trajectories. However, we do expect the method to perform less well when the mixture is more complicated or contains more than two source systems. Future studies could seek principles that guarantee the design of a better intermediate system for different chaotic signals.