Machine-learning construction of a model for a macroscopic fluid variable using the delay-coordinate of a scalar observable

TL;DR

This paper presents a machine-learning approach using reservoir computing to model a macroscopic fluid variable from scalar time-series data without prior physical knowledge, successfully capturing chaotic dynamics.

Contribution

It introduces a data-driven reservoir computing model that effectively reconstructs fluid dynamics from scalar observations, emphasizing the importance of delay-coordinate parameters.

Findings

The model accurately approximates the actual time-series over various intervals.

It captures key characteristics of the chaotic invariant set.

Proper delay-coordinate selection enhances model complexity and accuracy.

Abstract

We construct a data-driven dynamical system model for a macroscopic variable the Reynolds number of a high-dimensionally chaotic fluid flow by training its scalar time-series data. We use a machine-learning approach, the reservoir computing for the construction of the model, and do not use the knowledge of a physical process of fluid dynamics in its procedure. It is confirmed that an inferred time-series obtained from the model approximates the actual one in each of various time-intervals, and that some characteristics of the chaotic invariant set mimic the actual ones. We investigate the appropriate choice of the delay-coordinate, especially the delay-time and the dimension, which enables us to construct a model having a relatively high-dimensional attractor easily.

Tables

Table 1. Table 1. The list of variables and matrices in the reservoir computing.

variable
$𝐮\left(\in {𝐑}^M\right)$	input variable
$𝐫\left(\in {𝐑}^N\right)$	reservoir state vector
$𝐬\left(\in {𝐑}^M\right)$	actual output variable obtained from Navier-Stokes equation
$\widehat{𝐬}\left(\in {𝐑}^M\right)$	inferred output variable obtained from reservoir computing
$𝐀\left(\in {𝐑}^{N\times N}\right)$	weighted adjacency matrix
${𝐖}_{\text{in}}\left(\in {𝐑}^{M\times N}\right)$	linear input weight
${𝐖}_{\text{out}}\left(\in {𝐑}^{N\times M}\right)$	matrix used for translation from $𝐫$ to output variable $\widehat{𝐬}$
$𝐜\left(\in {𝐑}^M\right)$	vector used for translation from $𝐫$ to output variable $\widehat{𝐬}$
$\tilde{x}$	normalized variable of $x$

Table 2. Table 2. The list of parameters and their values used in the reservoir computing in each section.

parameter		Sec. 4	Sec. 5
$M$	dimension of input and output variables	14	Table. 3
$N$	dimension of reservoir state vector	3000	2000
$D$	parameter of determining $𝐀$	120	80
$\Delta t$	time step for reservoir dynamics	0.5
$T_0$	transient time for $𝐫$ to be converged	3750
$T$	training time	40000
$L_0$$\left(=T_0/\Delta t\right)$	number of iterations for the transient	7500
$L$ $\left(=T/\Delta t\right)$	number of iterations for the training	80000
$\rho$	maximal eigenvalue of $𝐀$	0.7
$\sigma$	scale of input weights in ${𝐖}_{\text{in}}$	0.5
$\alpha$	nonlinearity degree of reservoir dynamics	0.6
$\beta$	regularization parameter	0.1

Equations

$$\frac{d\varphi }{dt}=\mathbf{f}(\varphi ),$$

$$\mathbf{u}={\mathbf{h}}_1(\varphi )\in {\mathbb{R}}^M\text{and}\mathbf{s}={\mathbf{h}}_2(\varphi )\in {\mathbb{R}}^M.$$

$$\mathbf{r}\in {\mathbb{R}}^N(N\gg M),$$

$$\mathbf{r}(t+\Delta t)=(1-\alpha )\mathbf{r}(t)+\alpha \tanh (\mathbf{Ar}(t)+{\mathbf{W}}_{\text{in}}\mathbf{u}(t)),$$

$$\tanh (\mathbf{q})=(\tanh (q_1),\tanh (q_2),\cdots ,\tanh (q_N){)}^{\text{T}},$$

$$\hat{\mathbf{s}}(t)={\mathbf{W}}_{\text{out}}\mathbf{r}(t)+\mathbf{c}.$$

$$\sum _{l=1}^L\Vert ({\mathbf{W}}_{\text{out}}\mathbf{r}(l\Delta t)+\mathbf{c})-\mathbf{s}(l\Delta t){\Vert }^2+\beta [Tr({\mathbf{W}}_{\text{out}}{\mathbf{W}}_{\text{out}}^{\text{T}})],$$

$$\hat{\mathbf{s}}(t)={\mathbf{W}}_{\text{out}}^{\ast }\mathbf{r}(t)+{\mathbf{c}}^{\ast },$$

$${\mathbf{W}}_{\text{out}}^{\ast }$$

$${\mathbf{c}}^{\ast }$$

$$\mathbf{r}(t+\Delta t)=(1-\alpha )\mathbf{r}(t)+\alpha \tanh (\mathbf{Ar}(t)+{\mathbf{W}}_{\text{in}}\hat{\mathbf{s}}(t)),$$

$$\tilde{x}(t)=[x(t)-X_1]/X_2,$$

$$\displaystyle\begin{cases}\partial_{t}v-\nu\Delta v+(v\cdot\nabla)v+\nabla\pi=f,~{}\nabla\cdot v=0,~{}\mathbb{T}^{3}\times(0,\infty),\\ v\big{|}_{t=0}=v_{0}\quad\text{with $\nabla\cdot v_{0}=0$},~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\mathbb{T}^{3},\end{cases}$$

$$E(t)=\sum _{\kappa \in D}\sum _{\zeta =1}^3{\left({\mathcal{F}}_{[v_{\zeta }]}(\kappa ,t)\right)}^2,$$

$${\mathcal{F}}_{[v_{\zeta }]}(\kappa ,t):=\frac{1}{(2\pi {)}^3}{\int }_{{\mathbb{T}}^3}{v}_{\zeta }(x,t)e^{-i(\kappa \cdot x)}dx(\zeta =1,2,3),$$

$${\check{R}}_{\lambda }(t):=\frac{\sqrt{(2/3)E(t)}\lambda }{\nu }=\sqrt{\frac{20E(t{)}^2}{3\nu ϵ(t)}},$$

$$ϵ(t)=2\nu \sum _{\kappa \in D}\sum _{\zeta =1}^3{\left|\kappa \right|}^2{\left({\mathcal{F}}_{[v_{\zeta }]}(\kappa ,t)\right)}^2,$$

$$\lambda ={\left(\frac{15\nu (2/3)E(t)}{ϵ(t)}\right)}^{1/2},$$

$$R_{\lambda }(t)=\sum _{l=99}^0{\check{R}}_{\lambda }(t-l\Delta t^{\ast })/100,$$

$$\mathbf{u}(t)$$

$$\mathbf{s}(t)$$

$$C(x)=\frac{\frac{1}{J}\sum _{j=0}^{J-1}(R_{\lambda }(t_0+j\Delta t^{\ast })-{\bar{R}}_{\lambda })(R_{\lambda }(t_0+j\Delta t^{\ast }+x)-{\bar{R}}_{\lambda })}{\sqrt{\frac{1}{J}\sum _{j=0}^{J-1}(R_{\lambda }(t_0+j\Delta t^{\ast })-{\bar{R}}_{\lambda }{)}^2}\sqrt{\frac{1}{J}\sum _{j=0}^{J-1}(R_{\lambda }(t_0+j\Delta t^{\ast }+x)-{\bar{R}}_{\lambda }{)}^2}},$$

$$\{\begin{array}{l}\text{(i) the time average along }|{\hat{s}}_1(t^{\prime })|<3\text{ for }{t}^{\prime }\le 3000,\\ \text{(ii) the error }{\epsilon }_2(t^{\prime })=|s_1(t^{\prime })-{\hat{s}}_1(t^{\prime })|=|{\tilde{R}}_{\lambda }(t^{\prime })-{\hat{\tilde{R}}}_{\lambda }(t^{\prime })|<e_{60}\text{ for all }{t}^{\prime }\le 60,\\ \text{(iii) the error }{\epsilon }_2(t^{\prime })<e_{90}\text{ for all }{t}^{\prime }\le 90,\end{array}$$

Full text

Kengo Nakai

Graduate School of Mathematical Sciences, The University of Tokyo

3-8-1 Komaba, Tokyo 153-0041, Japan

Yoshitaka Saiki

Graduate School of Business Administration, Hitotsubashi University

2-1 Naka, Kunitachi, Tokyo 186-8601, Japan

and

JST PRESTO, 4-1-8 Honcho, Kawaguchi-shi, Saitama 332-0012, Japan

and

Institute for Physical Science and Technology, University of Maryland

College Park, MD 20742, USA

1. Introduction

Reservoir computing is a brain-inspired machine-learning technique that employs a data-driven dynamical system. The framework was proposed as Echo-State Networks [Jaeger, 2001, Jaeger and Haas, 2004] and Liquid-State Machines [Maass et al., 2002], and it has been found to be effective in the inference of a future such as time-series, frequency spectra and the Lyapunov spectra [Antonik et al., 2018, Ibáñez-Soria et al., 2018, Inubushi and Yoshimura, 2017, Lu et al., 2017, Pathak et al., 2018, Pathak et al., 2017, Verstraeten et al., 2007].

A reservoir is a recurrent neural network whose internal parameters are not adjusted to fit the data in the training process. Only an output layer is trained. Therefore, the total computational costs are relatively low in comparison with many other machine learning techniques having the same dimensional neural networks. Many physical phenomena including a fluid flow are deterministic, and thus can be described by a high-dimensional dynamical system, even thought they have a complex behavior. That is why the reservoir computing with a high-dimensional neural networks can be useful for the construction of a model for such a phenomenon.

In our previous paper [Nakai and Saiki, 2018], we infer both microscopic and macroscopic behaviors of a three-dimensional chaotic fluid flow using reservoir computing. We presented two ways of inference of the complex behavior: the first, called partial inference, requires continued knowledge of partial time-series data during the inference as well as past time-series data, while the second, called full inference, requires only past time-series data as training data. For the first case, we are able to infer long-time motion of microscopic fluid variables. For the second case, we showed that the reservoir dynamics constructed from only past data of energy functions can infer the future behavior of energy functions and reproduced the energy spectrum.

In various experiments and observations of high-dimensional complex phenomena, there are usually much smaller number of measurements than the Lyapunov dimensions of the attractor. Even in such cases, we can efficiently construct a dynamical model by generating high-dimensional input data $\mathbf{u}$ for the reservoir computing by using the delay-coordinate [Nakai and Saiki, 2018, Sauer et al., 1991, Takens, 1981].

The current paper focuses on the model construction and the full-inference of a macroscopic variable, the Taylor microscale Reynolds number, when the scalar time-series is accessible as measurements. We evaluate the model in many ways, and discuss details of the appropriate choice of the delay-coordinate created from the single observable. This will be useful for readers who wish to construct a reservoir model by themselves.

After reviewing the procedure of the reservoir computing in Sec. 2 and the generation of time-series data of a fluid flow in Sec. 3, we show that the constructed reservoir model recovers various properties of a fluid flow obtained from the Navier-Stokes equation in Sec. 4. We investigate the effective choice of delay-coordinate in order to construct a model in Sec. 5. We summarize our results in Sec. 6.

2. Reservoir computing

Reservoir computing is recently used in the inference of complex dynamics [Ibáñez-Soria et al., 2018, Lu et al., 2018, Lu et al., 2017, Pathak et al., 2018, Pathak et al., 2017]. It focuses on the determination of a linear function from the reservoir state vector to variables to be inferred (see eq. (5)). Here we review the outline of the method [Jaeger and Haas, 2004, Lu et al., 2017]. In this paper, we construct a model dealing with so called full-inference, in which there is no observable data in the inference phase [Nakai and Saiki, 2018].

We consider a dynamical system

[TABLE]

together with a pair of $\phi$-dependent, vector valued variables

[TABLE]

We seek a method for using the knowledge of $\mathbf{u}$ to determine an estimate $\hat{\mathbf{s}}$ of $\mathbf{s}$ as a function of time when direct measurement of $\mathbf{s}$ is not available. We have a knowledge $\mathbf{u}$ and $\mathbf{s}$ during the training phase for $t\leq T$, $\mathbf{u}$ and $\mathbf{s}$ are unknown during the inference phase for $t>T$. Therefore, $\mathbf{u}$ during the inference phase is replaced by $\hat{\mathbf{s}}$ in the previous step. See eq. (8) for the detail.

The dynamics of the reservoir state vector

[TABLE]

is defined by the neural network

[TABLE]

where $\Delta t$ is a relatively short time step, and

[TABLE]

for a vector $\mathbf{q}=(q_{1},q_{2},\cdots,q_{N})^{\text{T}}$. Here, T represents the transpose of a matrix. The matrix $\mathbf{A}$ is a weighted adjacency matrix, and the $M$-dimensional input $\mathbf{u}$ is fed in to the $N$ reservoir nodes via a linear input weight matrix denoted by $\mathbf{W}_{\text{in}}$. The parameter $\alpha$ ($0<\alpha\leq 1$) adjusts the nonlinearity of the dynamics of $\mathbf{r}$, and is chosen depending upon the complexity of the dynamics of measurements and the time step $\Delta t$.

Each row of $\mathbf{W}_{\text{in}}$ has one nonzero element, chosen from a uniform distribution on $[-\sigma,\sigma]$. The matrix $\mathbf{A}$ is chosen from a sparse random matrix in which the fraction of nonzero matrix elements is $D/N$, so that the average degree of a reservoir node is $D$. The $D$ non-zero components are chosen from a uniform distribution on $[-1,1]$. Then we uniformly rescale all the elements of $\mathbf{A}$ so that the largest value of the magnitudes of its eigenvalues becomes $\rho$.

The output, which is a $M$-dimensional vector, is taken to be a linear function of the reservoir state vector $\mathbf{r}$:

[TABLE]

The reservoir state vector $\mathbf{r}$ evolves following eq. (2) with input $\mathbf{u}(t)$, starting from random initial state $\mathbf{r}(-T_{0})$ whose elements are chosen from $(0,1]$ in order not to diverge, where $T_{0}=L_{0}\Delta t~{}(\gg 1)$ is the transient time for $\mathbf{r}(t)$ ($t>0$) to be on the attractor. We obtain $L=T/\Delta t$ steps of reservoir state vectors $\{\mathbf{r}(l\Delta t)\}_{l=1}^{L}$ by iterating eq. (2), while we record the variables $\{\mathbf{s}(l\Delta t)\}_{l=1}^{L}$ by using the actual measurements from eq. (1) for the training phase.

Determination of $\mathbf{W}_{\text{out}}$ and $\mathbf{c}$. We determine $\mathbf{W}_{\text{out}}$ and $\mathbf{c}$ so that the reservoir output $\hat{\mathbf{s}}$ (eq. (3)) approximates the measurement $\mathbf{s}$ for $0<t\leq T$ (training phase), which is a training process in the reservoir computing. We determine them by minimizing the following quadratic form with respect to $\mathbf{W}_{\text{out}}$ and $\mathbf{c}$:

[TABLE]

where $\|\mathbf{q}\|^{2}=\mathbf{q}^{\text{T}}\mathbf{q}$ for a vector $\mathbf{q}$, and the second term is a regularization term introduced to avoid overfitting $\mathbf{W}_{\text{out}}$ for $\beta\geq 0$. When the training is successful, $\hat{\mathbf{s}}(t)$ should approximate the desired unmeasured quantity $\mathbf{s}(t)$ for $t>T$ (inference phase). Following eq. (3), we obtain

[TABLE]

where $\mathbf{W}^{*}_{\text{out}}$ and $\mathbf{c}^{*}$ denote the solution for the minimizers of the quadratic form (4):

[TABLE]

where $\overline{\mathbf{r}}=\sum^{L}_{l=1}\mathbf{r}(l\Delta t)/L$, $\overline{\mathbf{s}}=\sum^{L}_{l=1}\mathbf{s}(l\Delta t)/L$, and $\mathbf{I}$ is the $N\times N$ identity matrix, $\delta\mathbf{R}$ (respectively, $\delta\mathbf{S}$) is the matrix whose $l$-th column is $\mathbf{r}(l\Delta t)-\overline{\mathbf{r}}$ (respectively, $\mathbf{s}(l\Delta t)-\overline{\mathbf{s}}$) (see [Lukosevivcius and Jaeger, 2009] P.140 and [Tikhonov and Arsenin, 1977] Chapter 1 for details).

In the inference phase for $t>T$, eq.(2) is written as

[TABLE]

by setting $\mathbf{u}(t)$ as $\hat{\mathbf{s}}(t)$ obtained from eq. (5).

We define a reservoir model by eqs. (5) and (8) under the values determined by eqs. (6) and (7) through the training data in a time-interval $[0,T]$. The main variables and matrices in the reservoir computing are summarized in Table 1.

Normalization of a variable. In order to consider the effect of all the variables equally, we take the normalized value $\tilde{x}(t)$ for each variable $x(t)$, which will be used in the procedure of our reservoir computing:

[TABLE]

where $X_{1}$ is the mean value and $X_{2}$ is the variance. When we reconstruct $x(t)$ in the inference phase from $\tilde{x}(t)$, we employ $X_{1}$ and $X_{2}$ obtained in the training phase. Due to the normalization we can avoid adjustments of $\sigma$.

Parameter choice. We apply a method of reservoir computing described above in order to construct a model. The sets of parameter values used are shown in Table 2.

3. Generation of a fluid flow data

Modelling and inference of a fluid flow are important problems in many areas [Di Leoni et al., 2018, Nakai and Saiki, 2018]. In this paper, we construct a model for a macroscopic variable of a fluid flow, especially the time-dependent “Taylor microscale Reynolds number” which reflects the degree of complexity in the fluid flow. We generate training data by the direct numerical simulation of the Navier-Stokes equation, which is also used for the reference data in the inference phase in order to evaluate the constructed reservoir model. It should be remarked that the Navier-Stokes equation and its physical property are not considered at all when constructing a reservoir model.

Generation of training data. In order to generate measurements of the reservoir computing, we employ the direct numerical simulation of the incompressible three-dimensional Navier-Stokes equation under periodic boundary conditions:

[TABLE]

where ${\mathbb{T}}=[0,1)$, $\nu>0$ is a viscosity parameter, $\pi(x,t)$ is pressure, and $v(x,t)=(v_{1}(x,t),v_{2}(x,t),v_{3}(x,t))$ is velocity. Throughout this paper, we set $\nu=0.058$, under which the fluid flow shows an intermittent behavior between laminar and bursting states. See such a behavior in the bottom panel of Fig. 1. We use the Fourier spectral method [Ishioka, 1999] with $N_{0}(=9)$ modes in each of three directions, meaning that the system is approximated by $2(2N_{0}+1)^{3}~{}(=13718)$-dimensional ordinary differential equations (ODEs). The ODEs are integrated by the 4th-order Runge–Kutta scheme, and the forcing is input into the low-frequency variables at each time step so as to preserve the energy of the low-frequency part. See [Ishioka, 1999, Nakai and Saiki, 2018] for the details.

Reynolds number $R_{\lambda}$. We focus on the time-series of the Taylor microscale Reynolds number, a macroscopic variable representing the degree of complexity of a fluid flow. The total energy $E(t)$ is defined by

[TABLE]

where

[TABLE]

and $D=\{(\kappa_{1},\kappa_{2},\kappa_{3})\in\mathbb{Z}^{3}\mid\kappa_{1},\kappa_{2},\kappa_{3}\in[-9,9]\}$. The Taylor microscale Reynolds number $\check{R}_{\lambda}(t)~{}$[Ishihara and Kaneda, 2003] is defined as follows:

[TABLE]

where

[TABLE]

is the average rate of energy dissipation per unit mass and

[TABLE]

is the characteristic length of a turbulent fluid flow. The length roughly corresponds to that of an energy input in this study.

In order to get rid of the high-frequency fluctuation, we take the short-time average

[TABLE]

where $\Delta t^{*}=0.05$ is the time step of the integration of the Navier-Stokes equation. This helps us to obtain essential low-frequency dynamics of a Reynolds number and construct a model with less computational costs with lower dimension $N$ of the reservoir state vectors. The averaged Reynolds number ${R_{\lambda}}$ will be called the Reynolds number, and the time-series generated by the direct numerical simulation in the inference phase will be called the “actual data”.

4. Construction of a model for a macroscopic variable: Reynolds Number

Using the reservoir computing discussed in Sec. 2, we construct a model by training a time-series data of the Reynolds number ${R_{\lambda}}$ (see Sec. 3) that shows an intermittent behavior between laminar and bursting states. For its purpose a delay-coordinate vector created from a scalar observable is introduced to the input and output variables.

4.1. Construction

Delay-coordinate. The choice of variables for the reservoir model is significant. Here, we introduce an $M$-dimensional delay-coordinate vector of the Reynolds number with a delay-time $\Delta\tau$ as input and output variables $\mathbf{u}(t)=(u_{1}(t),u_{2}(t),\cdots,u_{M}(t))^{\text{T}}$ and $\mathbf{s}(t)=(s_{1}(t),s_{2}(t),\cdots,s_{M}(t))^{\text{T}}$ in eq. (1), that is,

[TABLE]

The appropriate choice of the dimension $M$ and the delay-time $\Delta\tau$ of the delay-coordinate will be discussed in Sec. 5.

Determination of a model. Under the parameters listed in Table 2 and randomly chosen matrices $\mathbf{A}$ and $\mathbf{W}_{\text{in}}$, we find a candidate of a reservoir model by fixing $\mathbf{W}^{*}_{\text{out}}$ and $\mathbf{c}^{*}$ following the procedure explained in Sec 2. If the candidate passes a certain criteria concerning the short time inference, the candidate is considered as a model. See Sec. 5 for the details of the criteria. Remark that although we can use a training data as some delay components of input data when $t-(M-1)\Delta\tau<T$, we do not use any training data in the inference phase. Hereafter throughout this section, we choose one of the models, and fix the corresponding set of values $\mathbf{A}$, $\mathbf{W}_{\text{in}}$, $\mathbf{W}^{*}_{\text{out}}$ and $\mathbf{c}^{*}$.

4.2. Evaluation of the model

We evaluate the constructed reservoir model for the Reynolds number from several points of view by comparing its property with that of the actual data obtained from the direct numerical simulation of Navier-Stokes equation.

Time-series. We confirm that an inference of a time-series of the Reynolds number $s_{1}=\tilde{R}_{\lambda}$ is successful for some time after finishing the training phase. The time-series of the inferred variable $\hat{s}_{1}=\hat{\tilde{R}}_{\lambda}(t)$ ($t>T$) is shown with the actual data $\hat{s}_{1}={\tilde{R}}_{\lambda}(t)$ obtained from the direct numerical simulation of Navier-Stokes equation in the top left panel of Fig. 1. The failure in the long-term time-series inference is inevitable just due to the sensitive dependence on the initial condition of a chaotic property of the fluid flow. The two types of errors between the inferred value and the actual one are shown in the top right panel of Fig. 1. Moreover, the long-time behavior of $\hat{\tilde{R}}_{\lambda}(t)$ is shown in the bottom panel, which has qualitatively similar intermittent behaviors to the actual one, intermittent switching between the state of low amplitude fluctuations (laminar state) and the state of high amplitude fluctuations (burst state).

Delay Property. As we employ the delay-coordinate vector for input and output variables of the reservoir computing (eqs. (9),(10)), the relation $s_{1}(t)=s_{m}(t+(m-1)\Delta\tau)$ holds for any $m~{}(m=2,\cdots,M)$ during the training phase. The corresponding relation should also be satisfied in the inference phase. We show the time-series $\hat{s}_{1}(t)$ and $\hat{s}_{14}(t+13\Delta\tau)$ in Fig. 2, which satisfies the relation $\hat{s}_{1}(t)\approx\hat{s}_{14}(t+13\Delta\tau)$. We can confirm that for almost all $t$ the relation $\hat{s}_{1}(t)\approx\hat{s}_{m}(t+(m-1)\Delta\tau)$ is satisfied for any $m$. The results imply that our reservoir computing successfully learns the delay property only through training such data.

Poincaré plane. We investigate the chaotic set computed from a model trajectory to see whether the inferred chaotic set mimics the actual one. For its purpose we describe the Poincaré plane in comparison with that computed from a trajectory of the direct numerical simulation of the Navier-Stokes equation with the same time length in the Fig.  3. Although the sections are similar to each other, but they are not very close to each other. This may be because the length of the intermittent trajectory is not enough to cover various regions especially in the bursting state.

Distribution. Density distributions computed from two inferred trajectories $\{\hat{s}_{1}(t)\}$ and those from two actual trajectories $\{{s}_{1}(t)\}$ are shown in Fig.  4. We can observe that the distributions computed from trajectories of time lengths 5000 are fluctuating, but the inferred distributions seem to have similar properties to the actual distributions. Relatively large fluctuations in distributions for $|\hat{s}_{1}(t)|>1$ should be due to the intermittency.

The reservoir model can be used to infer time-series of another time-interval. We obtained a model just by training the data and it enables us to infer short-time behavior, the shape of an attractor and the density distribution. Here we confirm that the model constructed using a certain training data has the ability to infer a short-time behavior of the Reynolds number for the totally different time-interval. In Fig. 5, the inferred time-series is shown in comparison with the actual one. For this inference we use the same reservoir model as is used in Fig. 1, and only change the initial condition. This figure supports the accuracy of the constructed reservoir model.

In Fig. 6, by using the same model the inference of time-series of the Reynolds number in many different time intervals are shown. For each time-interval, we confirm that the short time inference is successful. This implies that the obtained model can describe the dynamics of the Reynolds number. Remark that the top middle panel in Fig. 6 corresponds to Fig. 5.

5. Choice of delay-coordinate.

We use an $M$-dimensional delay-coordinate vector with a delay-time $\Delta\tau$ (eqs. (9),(10)) as input and output variables $\mathbf{u}$ and $\mathbf{s}$ in eq. (1). In this section we investigate the appropriate choice of time-delay $\Delta\tau$ and the dimension $M$.

Time-correlation. The auto-correlation function $C({x})$ along a trajectory $\{R_{\lambda}(t))\}$ with respect to the time-difference $x$ is computed by

[TABLE]

where $\bar{R}_{\lambda}$ is the time average of $R_{\lambda}(t)$, $\Delta t^{*}$ is the time step of the discrete trajectory, and $t_{0}$ is an initial time of a trajectory. In Fig. 7, we show the auto-correlation function $C(x)$ for a trajectory $\{R_{\lambda}(t)\}$ with respect to the delay-time $x$.

It is observed from Fig. 7 (right) that as $x$ increases from 0, $C(x)$ goes below 0.7 and 0.3 when $x\approx 3.0$ and $5.0$, respectively.

The observation suggests that the value of the delay-time $\Delta\tau$ is to be chosen around $3.0$-$5.0$. If $\Delta\tau<3.0$, the consecutive two components of a delay-coordinate vector in (9), $R_{\lambda}(t)$ and $R_{\lambda}(t-\Delta\tau)$ behave too similarly, and if $\Delta\tau>5.0$, the consecutive two components behave too differently, and some dynamics to be captured may be missing.

Delay-time and dimensions. Based on the above implication about the auto-correlation function in Fig. 7, we investigate the effective delay-time $\Delta\tau$ and dimensions $M$. We focus on the delay-time $\Delta\tau\approx 3.0$-$5.0$ in Table 3. We infer time-series of the Reynolds number $R_{\lambda}(t)$ (actually its normalized value $\tilde{R}_{\lambda}(t)$) using the procedure in Sec. 2 by employing the delay-coordinate in eq. (5). We tried 8160 cases for each set of parameters $(\Delta\tau,M)$, for which matrices $\mathbf{A}$ and $\mathbf{W}_{\text{in}}$ are chosen randomly, and the number of successful cases are counted in TABLE 3. We say that the inference of $s_{1}(t^{\prime})$ ($t^{\prime}=t-T>0$) is successful if the conditions

[TABLE]

hold, where the criteria $(e_{60},e_{90})$ are set as (a)$(0.14,0.30)$ and (b)$(0.13,0.17)$. Remark that the condition (i) is given so as to get rid of a candidate which diverges within a short time, as $|s_{1}(t^{\prime})|<3$ for almost all $t$ even in the bursting region. For each case we use the same training data as in Fig. 1.

It is observed that the delay-time $\Delta\tau$ and the dimension $M$ of the delay-coordinate are chosen so that $\Delta\tau\approx 4.0$-$4.5$, and $M\Delta\tau\approx 55$-$60$, which correspond to $C(\Delta\tau)\approx 0.45$-$0.55$ and its envelope $C_{e}(M\Delta\tau)\approx 0.35$-$0.40$, respectively (see the left panel of Fig. 7 for the envelope $C_{e}$). For $\Delta\tau=4.0$ and $M=14,15$ by computing 16 times more cases, we confirmed that the rate of successful trials does not change much. In addition, even when we change the value of $N$ such as $1000$ or $3000$, we obtain almost the same results.

6. Summary and discussion

By training a time-series data of a macroscopic quantity the Reynolds number of a fluid flow, we construct a closed form system describing its intermittent behavior between laminar and bursting states. For the model construction, we do not use the knowledge of a physical process. We evaluate the obtained model in many ways. In particular, the model is confirmed to have a time-series predictability in many time intervals.

In order to construct a model from a scalar time-series data, we introduce a time-delay coordinate. From our investigations, the time-delay should be chosen to be the lowest value $\Delta\tau(>0)$ so that the auto-correlation function $C$ is $0.45<C(\Delta\tau)<0.55$ at the first time, and that the dimension $M$ of the delay-coordinate should be chosen so that the envelope $C_{e}$ of the auto-correlation function $C$ is $0.35<C_{e}(M\Delta\tau)<0.40$.

It should be remarked that the obtained reservoir model has a chaotic set on which a trajectory approximates the actual one, but the set is not an attractor. This may be due to the lack of a training data, especially in the bursting state. The clarification is remained as a future study.

Acknowledgements

KN was supported by the Leading Graduate Course for Frontiers of Mathematical Sciences and Physics (FMSP) at the University of Tokyo. YS was supported by the JSPS KAKENHI Grant No.17K05360 and JST PRESTO JPMJPR16E5. Part of the computation was supported by the Collaborative Research Program for Young $\cdot$ Women Scientists of ACCMS and IIMC, Kyoto University.