## 2 Methods  
The SINDy framework [24] is based on differential polynomial and sparse regression. It is assumed that the equation of motion to be identified,  
\[\dot{\mathbf{x}}(t)=\mathbf{f}(\mathbf{x}(t))=(\mathbf{\theta}(\mathbf{x}^{\mathrm{T}}(t))\mathbf{ \Xi})^{\mathrm{T}}, \tag{1}\]  
has a set of polynomial candidate functions \(\mathbf{\theta}(\mathbf{x}^{\mathrm{T}}(t))\) regarding the state vector \(\mathbf{x}\in\mathbb{R}^{n_{x}}\) and a static matrix of coefficient \(\mathbf{\Xi}\in\mathbb{R}^{n_{\xi}\times n_{x}}\). \(\mathbf{\xi}_{i}\in\mathbb{R}^{n_{\xi}}\), the \(i\)-th column of the matrix \(\mathbf{\Xi}\), corresponds to the \(i\)-th dimension of the derived state vector \(\dot{\mathbf{x}}(t)\). The \(n_{\xi}\)-dimensional vector of candidate functions \(\mathbf{\theta}(\mathbf{x}^{\mathrm{T}}(t))\) includes the constant \(1\) and all powers of the vector \(\mathbf{x}(t)\) up to a certain order [24]. For example, \(\mathbf{\theta}([x_{1},x_{2}])|_{Ord=2}=[1,x_{1},x_{2},x_{1}^{2},x_{1}x_{2},x_{2} ^{2}]\). These signals of \(\dot{\mathbf{x}}\) and \(\mathbf{\theta}(\mathbf{x}^{\mathrm{T}})\) from measurements and calculation can be written in the following matrices:  
\[\dot{\mathbf{X}}= \begin{bmatrix}\dot{\mathbf{x}}^{\mathrm{T}}(t_{1})\\ \dot{\mathbf{x}}^{\mathrm{T}}(t_{2})\\ \vdots\\ \dot{\mathbf{x}}^{\mathrm{T}}(t_{n_{t}})\end{bmatrix} ,\mathbf{\Theta}(\mathbf{X})= \begin{bmatrix} \mathbf{\theta}(\mathbf{x}^{\mathrm{T}}(t_{1}))\\ \mathbf{\theta}(\mathbf{x}^{\mathrm{T}}(t_{2}))\\ \vdots\\ \mathbf{\theta}(\mathbf{x}^{\mathrm{T}}(t_{n_{t}}))\end{bmatrix} . \tag{2}\]  
Thus, the solving process of Equation (1) can be cast into a regression task [24]:  
\[\dot{\mathbf{X}}=\mathbf{\Theta}(\mathbf{X})\mathbf{\Xi}. \tag{3}\]  
The model is established by identifying the coefficients because in most physical systems only a small subset of all candidate functions are relevant. The coefficient matrix \(\mathbf{\Xi}\) is sparse in the high-dimensional search space. In this situation, the Least-Squares (LS) method  
\[\mathbf{\hat{\xi}}_{i} =\operatorname*{argmin}_{\mathbf{\xi}_{i}\in\mathbb{R}^{n_{\xi}}}\{ \|\dot{\mathbf{X}}_{*,i}-\mathbf{\Theta}\mathbf{\xi}_{i}\|_{2}^{2}\} \tag{4}\] \[=(\mathbf{\Theta}^{\mathrm{T}}\mathbf{\Theta})^{-1}\mathbf{\Theta}^{\mathrm{T }}\mathbf{X}_{*,i} \tag{5}\]  
is not suitable for identification in practice because the parameters with the true value zero will be estimated as a non-zero value due to the noise in the signals. At the same time, overfitting occurs [39]. As one method of sparse regression, the LSST approach promotes sparseness by zeroing the estimated coefficients of the LS method below a fixed threshold \(\lambda>0\). The candidate functions corresponding to non-zero coefficients are selected for the next iteration of LS regression [24]. After a certain number of iterations, usually 10 times as suggested by [PERSON] [24], the algorithm converges to a fixed point [40]. This process is presented in Figure 1a. The parameter \(\lambda\) needs to be small enough to let the correct candidate functions be selected and large enough to avoid noise being picked up [40]. However, a suitable threshold cannot be found if the noise is so intensive that a zero coefficient is estimated to be greater than a non-zero coefficient. This situation often appears when the coefficients are of significantly different orders of magnitude. To overcome this issue, another sparse regression method was developed, which is based on LASSO. LASSO is an optimization method for estimating the coefficients \(\xi_{i}\) in the objective function,  
\[\hat{\xi}_{i}=\operatorname*{argmin}_{\xi_{i}\in\mathbb{R}^{n_{\xi}}}\{\frac{1}{2 n_{y}}\|\hat{X}_{*,i}-\mathbf{\Theta}\xi_{i}\|_{2}^{2}+\lambda\|\xi_{i}\|_{1}\}. \tag{6}\]  
Compared with the LS method, LASSO has additionally a penalty term \(\|\xi_{i}\|_{1}\), i.e., the \(l_{1}\) norm of the optimized vector.