We present a method to learn mean residence time and escape probability from data modeled by stochastic differential equations. This method is a combination of machine learning from data (to extract stochastic differential equations as models) and stochastic dynamics (to quantify dynamical behaviors with deterministic tools). The goal is to learn and understand stochastic dynamics based on data. This method is applicable to sample path data collected from complex systems, as long as these systems can be modeled as stochastic differential equations.

Stochastic dynamical systems are appropriate models for randomly influenced systems. Understanding the complex dynamical behaviors of these systems is a challenge in diverse areas of science and engineering. In deterministic dynamical systems, invariant manifolds and other invariant structures provide global information for dynamical evolution. For stochastic dynamical systems, better quantitative analysis and understanding is needed because of limitations in our current analytical skills or computation capability. Fortunately, researchers are increasingly using data-driven methods for system identification and the discovery of dynamics. In this work, we propose a new approach to determine some computable dynamical quantities from data, such as the mean residence time and escape probability, which offer insights into global dynamics under uncertainty. We demonstrate the algorithm to be effective and robust by reproducing known dynamics and evaluating errors for several prototypical stochastic dynamical systems with Brownian motions.

## I. INTRODUCTION

Stochastic dynamical systems arise in modeling molecular dynamics, mechanical and electrical engineering, climate dynamics, geophysical and environmental systems, among others. Advances in machine learning and data science are leading to new progresses in the analysis and understanding of complex dynamics for systems with massive observation data sets. Despite the rapid development of tools to extract governing equations from data, there has been slow progress in distilling quantities that may be used to explore stochastic dynamics. It is desirable to extract deterministic quantities that carry dynamical information of stochastic differential equations (SDEs). These deterministic quantities include moments for solution paths, probability density functions for solution paths, mean residence time, and escape probability.^{7,9} These concepts help us to understand various phenomena in complex systems under uncertainty.^{1}

In this present paper, we will present a new approach to extract the underlying deterministic quantities, mean residence time, and escape probability that describe certain aspects of stochastic dynamics. In fact, mean residence time for a stochastic dynamical system quantifies how long the system stays in a region, and escape probability describes the likelihood of a system transition from one regime to another. Fortunately, these deterministic quantities can be determined by solving an elliptic partial differential equation as in Duan,^{9} once the underlying stochastic differential equation model for the system evolution is discovered from data.

We assume that a data set is composed of sample paths governed by a stochastic differential equation (SDE) in $Rn$,

where $b$ is an n-dimensional vector function, $\sigma $ is an $n\xd7m$ matrix function, and $Bt$ is an m-dimensional Brownian motion. This is the customary, probabilistic way of writing the equation

Often, $b$ is called “drift” and $\sigma $ is called “diffusion.” Assume that $b$ and $\sigma $ satisfy an appropriate local Lipschitz condition as follows:

for $\Vert x\Vert \u2264N$, $\Vert y\Vert \u2264N$ and $N>0$. Here, the Lipschitz constant $KN$ depends on the positive number $N$. The generator for this SDE system is a linear second-order differential operator:

where $H$ denotes the Hessian matrix of a multivariate function and $Tr$ denotes the trace of a matrix. We view $A$ as a linear differential operator in Hilbert space $L2(Rn)$ with domain of definition $D(A)=H02(Rn)$. Furthermore, we assume that the generator $A$ is uniformly elliptic. That is for $x\u2208D$ and all $\xi \u2208Rn$, there exists a positive constant $C$ such that

We discuss mean residence time and escape probability as deterministic quantities that carry dynamical information for solution orbits of (1). For a bounded domain $D\u2282Rn$ (with boundary $\u2202D$), the first exit time for a solution orbit starting at $x\u2208D$ is a stopping time defined as

The mean residence time is defined as $u(x)\u225cE\tau D(\omega )$. It is the mean residence time of a particle initially at $x$ inside $D$ until the particle first hits the boundary $\u2202D$ or escapes from $D$. It turns out that the mean residence time $u$ of stochastic system (1) can be determined by solving a deterministic partial differential equation as follows:

where $A$ is the generator defined as (2). To prove the existence and uniqueness of $u$, we recall the Dynkin’s formula (Theorem 7.4.1 in Oksendal^{17}),

for $f$ in the domain of definition of the generator $A$, $Ex$ is the expectation with respect to the probability law $Qx$ induced by a solution process $Xt$ starting at $x\u2208Rn$. We take the continuous bounded boundary value $\varphi =0$ and the continuous inhomogeneous term $g=1$. By Theorem 9.3.3 in Oksendal,^{17} which is a consequence of Dynkin’s formula, we know that the linear expectation $Ex[\varphi (X\tau D)]+Ex[\u222b0\tau Dg(Xs)ds]$, which is just $E\tau D=u(x)$, solves $Au=\u22121$ with the boundary condition $u|\u2202D=0$. Thus, we can numerically compute mean residence time $u$ by solving elliptic partial differential equations, so we know how long the stochastic system (1) stays in the region $D$. For more details, see Duan.^{9}

The escape probability is the likelihood that an orbit starting inside a domain $D$ exits from this domain first through a specific part $\Gamma $ of the boundary $\u2202D$. Let $\Gamma $ be a subset of the boundary $\u2202D$. We define the escape probability $p(x)$ from $D$ through $\Gamma $ as the likelihood that $Xt$ starting at $x$ exits from $D$ first through $p(x)=P{X\tau \u2202D\u2208\Gamma}$. We will show that the escape probability $p(x)$ solves a linear elliptic partial differential equation, with a specifically chosen Dirichlet condition as follows:

where $A$ is the generator defined as (2). Taking

we have

This means that $E[\varphi (X\tau \u2202D(x))]$ is the escape probability $p(x)$ that we are looking for. We know $Ap=0$ together with $p|\Gamma =1$ and $p|\u2202D\u2216\Gamma =0$ by Theorem 9.2.14 in Okendal.^{17} For more details, see Duan.^{9} Moreover, we suppose that solution orbits, i.e., “particles” are initially uniformly distributed in $D$, then the average escape probability $P$ that a trajectory will leave $D$ through $\Gamma $ is $P=1D\u222bDp(x)dx$.

To obtain the mean residence time $u$ and escape probability $p$ from sample path data $XTrue$ [with the underlying model (1)], we propose the following machine learning algorithm: First, we collect sample path data $XTrue$ by sample-wisely simulating (1) via the Euler method (and treat these data as our observation data), then we try to learn the stochastic dynamical system model from these data with the following model ansatz:

where the basis $\Theta $ consists of polynomial functions ${1,x,y,\u2026,xNyN}$ together with the (generalized) time derivative of Brownian motion $dBt/dt$, and $\Xi $ is the coefficient (or coordinate, or weight) under this basis. Here, $\u22c5$ denotes scalar product. Each column $\xi k=[\xi 1k,\xi 2k,\u2026,\xi Nk]$ of $\Xi $ is a vector of coefficients, determining which terms are active in the right-hand side for one of the row equations in (1). Brunton *et al.*^{20} considered such a learned model, with a basis not containing noise (i.e., $dBt/dt$).

Furthermore, we set up a regression problem to determine coefficient $\Xi =[\xi 1,\u2026,\xi n]$ by minimizing the mean-square discrepancy (other metrics are possible) between the data $XTrue$ (many samples) and the solution $XLearn$ (many samples) for the learned (i.e., extracted) governing model (9). This will provide us the “learned” drift $b$ and “learned” diffusion $\sigma $, and thus, we also have the extracted stochastic model (9), together with the “learned” generator $A$.

Finally, with this learned generator $A$, we compute mean residence time $u$ and escape probability $p$, by solving the deterministic partial differential Eqs. (4)–(5) and (6)–(8), respectively. We verify that our algorithm is effective by estimating the error of the maximal mean residence time

and the error of the average escape probability

This paper is organized as follows. In Sec. II, we apply our method to learn a two-dimensional quasigeostrophic meandering jet model with additive noise and multiplicative noise and compute the mean residence time and escape probability. In Sec. III, we illustrate our method to learn two three-dimensional systems (a linear damped oscillator and the well-known Lorenz system). We summarize and conclude in Sec. IV.

## II. LEARNING TWO-DIMENSIONAL STOCHASTIC DYNAMICAL SYSTEMS

Zonal shear flows occur naturally in both oceans and the atmosphere. The Gulf Stream is a well-known example. Based on RAFOS floats observations of the Gulf Stream, Bower *et al.*^{4,5} viewed the fluid motion as a steady and eastward meander propagation in the moving frame and divided the velocity field into three regimes: a central jet, exterior retrograde motion, and intermediate closed circulations above meander troughs and below crests. There was no exchange occurring between the three regimes in this model. Fluid particles in the intermediate regime execute periodic motion but never escape. However, float observations reveal that several particle trajectories pass through meander crests and troughs and then leave jet, which indicate that exchange does occur across some part of the Gulf Stream. Samelson^{15} modified the basic model and considered three different types of variability: a time-dependent spatially uniform meridional velocity superimposed on the basic flow, a time-dependent meander amplitude, and a propagating plane wave superimposed on the basic flow. In particular, he only considered time-periodic variabilities. Subsequently, Samelson,^{16} del-Castillo-Negrete *et al.*,^{8} Pratt *et al.*,^{11} Beigie *et al.*,^{3} Duan *et al.*,^{6} and Brannan *et al.*^{7} have obtained a series of outcomes from the point of view of dynamical systems.

### A. Discovering a kinematical model for a two-dimensional meandering jet

As an approximate solution to the quasigeostrophic model for two-dimensional geophysical flows,^{13} the basic Bickley jet^{8} is

where $a=0.01$, $c=13(1+1\u221232\beta )$, and $k=6c$, $0\u2264\beta \u226423$. In this paper, we take $\beta =13$. Note $\beta =2\Omega rcos\u2061\theta $ is the meridional derivative of the Coriolis parameter, where $\Omega $ is the rotation rate of the earth, $r$ is the earth’s radius, and $\theta $ is the latitude. This stream function $\psi (x,y)$ defines the basic meandering jet system $x\u02d9=\u2212\psi y$, $y\u02d9=\psi x$. The phase portrait for this deterministic system is in Fig. 1 (left). As it shows no fluid exchange between the eddies and the jet, it is not an appropriate model for a meandering jet.

Thus, we incorporate random wind forcing and other fluctuations, either additive or multiplicative, in the model for the meandering jet, and consequently, we consider the following two stochastic dynamical systems:

and

where the noise intensity $\sigma $ satisfies $0<\sigma <1$, and $B1(t)$, $B2(t)$ are two independent Brownian motions. For both stochastic systems, a number of sample solution orbits are shown in Fig. 1 (middle, right). The central jet and two rows of recirculation eddies, which are called the northern and southern recirculation regions, are still visible. Outside the recirculation regimes are the exterior retrograde regimes. There is no exchange between regimes in unperturbed model, but exchange does occur in other two models with additive or multiplicative noises. These two stochastic models are more appropriate for modeling a meandering jet (as a simplified model for the Gulf Stream), as observations indicate fluid exchange.

A wide range of dynamical systems with intrinsic noise may be modeled with stochastic differential equations. However, identifying an appropriate SDE model from intermittent observations of the system is challenging, particularly if the dynamical process is nonlinear and the observations are noisy and indirect.^{10} Brunton *et al.*,^{20} Zhang *et al.*,^{18} and Dunker *et al.*^{21} presented some approaches to discover or learn governing physical laws from data, with underlying differential equations models.

We use our machine learning algorithm to data, from systems (10) and (11) or systems (12) and (13), to learn their governing laws in terms of the following model:

where $\Theta $ is a set of basis functions and $\Xi =[\xi 1,\xi 2]$ defines the coefficients (or weights). To achieve this, the idea is as follows: First, we collect a time series of the system state $(x(ti),y(ti)),i=1,\u2026,20000$, from systems (10) and (11) or systems (12) and (13) with noise intensity $\sigma =0.3$, numerically by the Euler method [as our observation data $(xTruei,yTruei)$]. Then, we construct a library $\Theta $ consisting of basis functions with polynomials ${1,x,y,\u2026,y5}$ and time derivative of Brownian motions ${dB1/dt,dB2/dt}$. We determine each column of coefficients $\xi k=[\xi 1k,\xi 2k,\u2026,\xi 23k],k=1,2$, by minimizing discrepancies $E\u2211i=120000\Vert xTruei\u2212xLearni\Vert 2$ and $E\u2211i=120000\Vert yTruei\u2212yLearni\Vert 2$. Here, $(xLearni,yLearni)$ is the numerical solution of the learned meandering jet model (14) and (15) by the Euler method (see Tables I and II).

Basis . | $x\u02d9$ . | $y\u02d9$ . |
---|---|---|

1 | 2.211 465 0 × 10^{+00} | 1.135 687 6 × 10^{+00} |

x | 9.167 092 9 × 10^{−01} | −3.105 340 6 × 10^{−01} |

y | −8.533 715 6 × 10^{+00} | −7.883 008 9 × 10^{+00} |

x^{2} | 2.158 193 3 × 10^{−01} | −5.888 226 7 × 10^{−02} |

xy | −3.540 679 3 × 10^{+00} | 1.208 688 5 × 10^{+00} |

y^{2} | 1.480 699 0 × 10^{+01} | 2.150 933 2 × 10^{+01} |

x^{3} | 2.316 303 0 × 10^{−02} | −2.068 737 6 × 10^{−03} |

x^{2} y | −6.440 459 8 × 10^{−01} | 1.655 640 0 × 10^{−02} |

xy^{2} | 5.216 326 3 × 10^{+00} | −2.304 683 1 × 10^{+00} |

y^{3} | −1.356 879 2 × 10^{+01} | −2.915 997 3 × 10^{+01} |

x^{4} | 6.686 967 2 × 10^{−04} | 1.118 142 7 × 10^{−02} |

x^{3} y | −7.380 868 5 × 10^{−02} | 6.441 070 0 × 10^{−02} |

x^{2} y^{2} | 5.297 775 4 × 10^{−01} | 2.256 811 7 × 10^{−01} |

xy^{3} | −3.578 395 4 × 10^{+00} | 2.294 521 4 × 10^{+00} |

y^{4} | 5.631 493 3 × 10^{+00} | 1.977 609 2 × 10^{+01} |

x^{5} | 2.851 040 0 × 10^{−04} | 1.238 256 6 × 10^{−03} |

x^{4} y | −8.195 872 2 × 10^{−04} | −8.028 851 6 × 10^{−04} |

x^{3} y^{2} | 4.227 503 3 × 10^{−02} | −3.688 687 4 × 10^{−02} |

x^{2} y^{3} | −1.161 596 8 × 10^{−01} | −1.736 499 4 × 10^{−01} |

xy^{4} | 9.944 926 4 × 10^{−01} | −9.095 955 7 × 10^{−01} |

y^{5} | −6.810 078 8 × 10^{−01} | −5.387 695 2 × 10^{+00} |

dB1/dt | 5.477 303 0 × 10^{−01} | −9.739 278 8 × 10^{−06} |

dB2/dt | 2.707 655 1 × 10^{−07} | 5.477 183 7 × 10^{−01} |

Basis . | $x\u02d9$ . | $y\u02d9$ . |
---|---|---|

1 | 2.211 465 0 × 10^{+00} | 1.135 687 6 × 10^{+00} |

x | 9.167 092 9 × 10^{−01} | −3.105 340 6 × 10^{−01} |

y | −8.533 715 6 × 10^{+00} | −7.883 008 9 × 10^{+00} |

x^{2} | 2.158 193 3 × 10^{−01} | −5.888 226 7 × 10^{−02} |

xy | −3.540 679 3 × 10^{+00} | 1.208 688 5 × 10^{+00} |

y^{2} | 1.480 699 0 × 10^{+01} | 2.150 933 2 × 10^{+01} |

x^{3} | 2.316 303 0 × 10^{−02} | −2.068 737 6 × 10^{−03} |

x^{2} y | −6.440 459 8 × 10^{−01} | 1.655 640 0 × 10^{−02} |

xy^{2} | 5.216 326 3 × 10^{+00} | −2.304 683 1 × 10^{+00} |

y^{3} | −1.356 879 2 × 10^{+01} | −2.915 997 3 × 10^{+01} |

x^{4} | 6.686 967 2 × 10^{−04} | 1.118 142 7 × 10^{−02} |

x^{3} y | −7.380 868 5 × 10^{−02} | 6.441 070 0 × 10^{−02} |

x^{2} y^{2} | 5.297 775 4 × 10^{−01} | 2.256 811 7 × 10^{−01} |

xy^{3} | −3.578 395 4 × 10^{+00} | 2.294 521 4 × 10^{+00} |

y^{4} | 5.631 493 3 × 10^{+00} | 1.977 609 2 × 10^{+01} |

x^{5} | 2.851 040 0 × 10^{−04} | 1.238 256 6 × 10^{−03} |

x^{4} y | −8.195 872 2 × 10^{−04} | −8.028 851 6 × 10^{−04} |

x^{3} y^{2} | 4.227 503 3 × 10^{−02} | −3.688 687 4 × 10^{−02} |

x^{2} y^{3} | −1.161 596 8 × 10^{−01} | −1.736 499 4 × 10^{−01} |

xy^{4} | 9.944 926 4 × 10^{−01} | −9.095 955 7 × 10^{−01} |

y^{5} | −6.810 078 8 × 10^{−01} | −5.387 695 2 × 10^{+00} |

dB1/dt | 5.477 303 0 × 10^{−01} | −9.739 278 8 × 10^{−06} |

dB2/dt | 2.707 655 1 × 10^{−07} | 5.477 183 7 × 10^{−01} |

Basis . | $x\u02d9$ . | $y\u02d9$ . |
---|---|---|

1 | 4.679 158 4 × 10^{+00} | −7.586 793 5 × 10^{+00} |

x | 6.300 678 4 × 10^{−01} | −2.108 732 6 × 10^{+00} |

y | −2.548 358 6 × 10^{+01} | 4.333 514 1 × 10^{+01} |

x^{2} | 2.610 926 3 × 10^{−01} | −1.595 241 1 × 10^{−01} |

xy | −1.742 074 1 × 10^{+00} | 1.056 439 5 × 10^{+01} |

y^{2} | 6.116 905 8 × 10^{+01} | −9.472 222 1 × 10^{+01} |

x^{3} | −3.182 419 0 × 10^{−03} | −1.005 654 0 × 10^{−01} |

x^{2} y | −9.661 535 4 × 10^{−01} | −2.467 371 8 × 10^{−01} |

xy^{2} | 8.853 546 5 × 10^{−01} | −2.153 239 1 × 10^{+01} |

y^{3} | −7.682 962 4 × 10^{+01} | 9.686 362 9 × 10^{+01} |

x^{4} | −1.734 162 0 × 10^{−04} | 4.713 976 6 × 10^{−03} |

x^{3} y | −6.630 334 9 × 10^{−03} | 2.637 329 7 × 10^{−01} |

x^{2} y^{2} | 1.143 603 0 × 10^{+00} | 1.300 516 0 × 10^{+00} |

xy^{3} | 1.074 806 0 × 10^{+00} | 2.069 413 0 × 10^{+01} |

y^{4} | 4.875 646 5 × 10^{+01} | −4.404 292 2 × 10^{+01} |

x^{5} | −1.691 796 4 × 10^{−04} | 2.235 049 5 × 10^{−04} |

x^{4} y | −3.953 219 7 × 10^{−03} | −5.468 421 9 × 10^{−03} |

x^{3} y^{2} | −1.339 436 5 × 10^{−02} | −1.884 372 4 × 10^{−01} |

x^{2} y^{3} | −4.833 223 9 × 10^{−01} | −1.007 315 1 × 10^{+00} |

xy^{4} | −8.735 959 0 × 10^{−01} | −7.798 103 5 × 10^{+00} |

y^{5} | −1.244 311 4 × 10^{+01} | 6.056 436 1 × 10^{+00} |

xdB1/dt | 5.477 237 4 × 10^{−01} | 4.326 505 7 × 10^{−05} |

ydB2/dt | 1.682 139 3 × 10^{−05} | 5.453 187 5 × 10^{−01} |

Basis . | $x\u02d9$ . | $y\u02d9$ . |
---|---|---|

1 | 4.679 158 4 × 10^{+00} | −7.586 793 5 × 10^{+00} |

x | 6.300 678 4 × 10^{−01} | −2.108 732 6 × 10^{+00} |

y | −2.548 358 6 × 10^{+01} | 4.333 514 1 × 10^{+01} |

x^{2} | 2.610 926 3 × 10^{−01} | −1.595 241 1 × 10^{−01} |

xy | −1.742 074 1 × 10^{+00} | 1.056 439 5 × 10^{+01} |

y^{2} | 6.116 905 8 × 10^{+01} | −9.472 222 1 × 10^{+01} |

x^{3} | −3.182 419 0 × 10^{−03} | −1.005 654 0 × 10^{−01} |

x^{2} y | −9.661 535 4 × 10^{−01} | −2.467 371 8 × 10^{−01} |

xy^{2} | 8.853 546 5 × 10^{−01} | −2.153 239 1 × 10^{+01} |

y^{3} | −7.682 962 4 × 10^{+01} | 9.686 362 9 × 10^{+01} |

x^{4} | −1.734 162 0 × 10^{−04} | 4.713 976 6 × 10^{−03} |

x^{3} y | −6.630 334 9 × 10^{−03} | 2.637 329 7 × 10^{−01} |

x^{2} y^{2} | 1.143 603 0 × 10^{+00} | 1.300 516 0 × 10^{+00} |

xy^{3} | 1.074 806 0 × 10^{+00} | 2.069 413 0 × 10^{+01} |

y^{4} | 4.875 646 5 × 10^{+01} | −4.404 292 2 × 10^{+01} |

x^{5} | −1.691 796 4 × 10^{−04} | 2.235 049 5 × 10^{−04} |

x^{4} y | −3.953 219 7 × 10^{−03} | −5.468 421 9 × 10^{−03} |

x^{3} y^{2} | −1.339 436 5 × 10^{−02} | −1.884 372 4 × 10^{−01} |

x^{2} y^{3} | −4.833 223 9 × 10^{−01} | −1.007 315 1 × 10^{+00} |

xy^{4} | −8.735 959 0 × 10^{−01} | −7.798 103 5 × 10^{+00} |

y^{5} | −1.244 311 4 × 10^{+01} | 6.056 436 1 × 10^{+00} |

xdB1/dt | 5.477 237 4 × 10^{−01} | 4.326 505 7 × 10^{−05} |

ydB2/dt | 1.682 139 3 × 10^{−05} | 5.453 187 5 × 10^{−01} |

### B. Mean residence time

In this subsection, we will consider stochastic quasigeostrophic meandering jet models, (10)–(11) and (12)–(13), to discover mean residence time by our approach. We take an eddy (Fig. 3) as the bounded domain $D$ (with boundary $\u2202D$ composed with meander trough and crest).

The mean residence time $u(x,y)$ of the stochastic system (10) and (11) with additive noise for a trajectory starting in the eddy satisfies the following elliptic partial differential equation as in (4) and (5):

The mean residence time, as solution of this elliptic partial differential equation, is denoted by $uTrue$. For the corresponding learned meandering jet model (14) and (15), with learned coefficients in Table I, we can also set up a similar elliptic partial differential equation for learned mean residence time $uLearn$.

Similarly, the mean residence time $u(x,y)$ of the stochastic system (12) and (13) with multiplicative noise, for a trajectory starting in the eddy, satisfies the following elliptic partial differential equation as in (4) and (5):

The mean residence time, as solution of this elliptic partial differential equation, is denoted by $uTrue$. For the corresponding learned meandering jet model (14) and (15), with learned coefficients in Table II, we can also set up a similar elliptic partial differential equation for learned mean residence time $uLearn$.

We use a finite element code to solve the preceding elliptic differential equations to get mean residence time $uTrue$ and $uLearn$, for both additive and multiplicative noise cases. Figure 4 shows the mean residence time $uLearn$, as $uTrue$ is barely distinguishable from $uLearn$ when they are plotted together. Thus, we compare the error between the maximal values of the mean residence times $uTrue$ and $uLearn$. In fact, the error of mean residence time between the learned system and the original system is ${error}(u)=max(x,y)\u2208D\Vert uLearn\u2212uTrue\Vert =2.2142\xd710\u22124$ for additive noise case, and ${error}(u)=max(x,y)\u2208D\Vert uLearn\u2212uTrue\Vert =1.3\xd710\u22123$ for multiplicative noise case.

### C. Escape probability

We again take the eddy (Fig. 3) to be the bounded domain $D$ [with boundary $\u2202D$ composed with meander trough (lower sub-boundary) and crest (upper sub-boundary)]. For the stochastic meandering jet system (10) and (11) with additive noise, the escape probability $p(x,y)$ of a fluid particle, starting at $(x,y)$ in an eddy and escaping through a boundary component $\Gamma $ (either the trough or crest), satisfies the elliptic partial differential equations (6)–(8),

For the corresponding learned meandering jet model (14) and (15), with learned coefficients in Table I, we can also set up a similar elliptic partial differential equation for the learned escape probability $pLearn$.

Similarly, for the stochastic meandering jet system (12) and (13) with multiplicative noise, the escape probability $p(x,y)$ of a fluid particle, starting at $(x,y)$ in an eddy and escaping through a boundary component $\Gamma $ (either the trough or crest), satisfies the elliptic partial differential Eqs. (6)–(8),

For the corresponding learned meandering jet model (14) and (15), with learned coefficients in Table II, we can also set up a similar elliptic partial differential equation for the learned escape probability $pLearn$.

Figures 5 and 6 show the finite element solution for the learned escape probability $pLearn$ alone, as $pTrue$ is barely distinguishable from $pLearn$ when they are plotted together. Thus, we compare the error between the average escape probability values.

For fluid particles initially uniformly distributed in $D$ (an eddy), we compute the average escape probability $P$ for particles leaving $D$ through the upper or lower sub-boundary $\Gamma $, given by $P=1D\u222b\u222bDp(x,y)dxdy$. The error of average escape probability, for fluid particles exiting from $D$ through the upper sub-boundary, between the learned system and the original system is ${error}(P)=\Vert PLearn\u2212PTrue\Vert =1.3949\xd710\u22125$ both for additive noise and multiplicative noise cases. The error of the average escape probability, for particles leaving $D$ through the lower sub-boundary, between the learned system and the original system is ${error}(P)=\Vert PLearn\u2212PTrue\Vert =7.4133\xd710\u22124$ both for additive noise and multiplicative noise cases.

## III. LEARNING THREE-DIMENSIONAL STOCHASTIC DYNAMICAL SYSTEMS

In this section, we will illustrate our method in two systems with additive noise. The first system is for a three-dimensional stochastic linear oscillator. The second is the more complex stochastic Lorenz system. In these two examples, data from direct numerical simulations (as “observation data”) are used to discover mean residence time and escape probability.

### A. A stochastic linear system

We consider a linear stochastic system

where $\u03f5$ is the noise intensity (taken to be $0.9$ here) and $B1(t)$, $B2(t)$, $B3(t)$ are three independent Brownian motions. With sample path data for this system, we try to discover the governing equation,

as in Brunton *et al.*,^{20} but we use a stochastic basis. Here, $\Theta $ is a set of basis functions and $\Xi =[\xi 1,\xi 2,\xi 3]$ are coefficients (or weights). As in Sec. II, we collect data $(x(ti),y(ti),z(ti)),i=1,\u2026,20000$, from samplewise simulations of the original system (26)–(28), and construct a library $\Theta $ consisting of polynomial ${1,x,y,z,\u2026,z4}$ and time derivatives of Brownian motions ${dB1/dt,dB2/dt,dB3/dt}$. We then solve a regression problem to determine the weights $\xi k=[\xi 1k,\u2026,\xi 38k],k=1,2,3$ (see Table III and Fig. 7).

Basis . | $x\u02d9$ . | $y\u02d9$ . | $z\u02d9$ . |
---|---|---|---|

1 | −3.901 872 6 × 10^{−06} | −2.527 055 3 × 10^{−06} | 9.962 128 3 × 10^{−06} |

x | −9.997 259 3 × 10^{−02} | 2.000 017 8 × 10^{+00} | −6.997 397 7 × 10^{−05} |

y | −1.999 946 3 × 10^{+00} | −9.996 519 8 × 10^{−02} | −1.371 972 5 × 10^{−04} |

z | 1.972 107 4 × 10^{−04} | 1.277 239 2 × 10^{−04} | −3.005 035 1 × 10^{−01} |

x^{2} | 4.195 535 5 × 10^{−05} | 2.717 246 7 × 10^{−05} | −1.071 189 8 × 10^{−04} |

xy | −2.942 892 2 × 10^{−06} | −1.905 969 8 × 10^{−06} | 7.513 692 0 × 10^{−06} |

xz | 1.225 505 7 × 10^{−04} | 7.937 011 3 × 10^{−05} | −3.128 919 4 × 10^{−04} |

y^{2} | 2.406 794 8 × 10^{−05} | 1.558 765 3 × 10^{−05} | −6.144 946 4 × 10^{−05} |

yz | 3.387 555 3 × 10^{−04} | 2.193 956 7 × 10^{−04} | −8.648 990 8 × 10^{−04} |

z^{2} | 4.021 428 3 × 10^{−03} | 2.604 485 8 × 10^{−03} | −1.026 737 4 × 10^{−02} |

x^{3} | −1.164 887 1 × 10^{−04} | −7.544 413 7 × 10^{−05} | 2.974 150 0 × 10^{−04} |

x^{2} y | −2.361 637 3 × 10^{−04} | −1.529 519 0 × 10^{−04} | 6.029 651 9 × 10^{−04} |

x^{2}z | −8.288 823 6 × 10^{−04} | −5.368 272 6 × 10^{−04} | 2.116 274 2 × 10^{−03} |

xy^{2} | −1.177 680 6 × 10^{−04} | −7.627 271 3 × 10^{−05} | 3.006 814 1 × 10^{−04} |

xyz | 3.517 902 6 × 10^{−05} | 2.278 376 4 × 10^{−05} | −8.981 788 9 × 10^{−05} |

xz^{2} | −1.828 761 2 × 10^{−03} | −1.184 400 7 × 10^{−03} | 4.669 130 7 × 10^{−03} |

y^{3} | −2.352 567 6 × 10^{−04} | −1.523 645 0 × 10^{−04} | 6.006 495 4 × 10^{−04} |

y^{2}z | −1.009 917 0 × 10^{−03} | −6.540 746 8 × 10^{−04} | 2.578 485 6 × 10^{−03} |

yz^{2} | −6.201 471 5 × 10^{−03} | −4.016 395 0 × 10^{−03} | 1.583 338 6 × 10^{−02} |

z^{3} | −1.242 848 7 × 10^{−02} | −8.049 333 6 × 10^{−03} | 3.173 198 9 × 10^{−02} |

x^{4} | −9.374 085 3 × 10^{−05} | −6.071 144 4 × 10^{−05} | 2.393 359 5 × 10^{−04} |

x^{3}y | −2.102 290 2 × 10^{−06} | −1.361 552 3 × 10^{−06} | 5.367 495 7 × 10^{−06} |

x^{3}z | 7.505 776 6 × 10^{−04} | 4.861 130 7 × 10^{−04} | −1.916 349 3 × 10^{−03} |

x^{2}y^{2} | −1.271 197 6 × 10^{−04} | −8.232 935 5 × 10^{−05} | 3.245 578 3 × 10^{−04} |

x^{2}yz | 1.801 066 1 × 10^{−03} | 1.166 463 9 × 10^{−03} | −4.598 420 6 × 10^{−03} |

x^{2}z^{2} | 2.218 039 3 × 10^{−03} | 1.436 517 4 × 10^{−03} | −5.663 022 3 × 10^{−03} |

xy^{3} | −8.741 210 4 × 10^{−06} | −5.661 261 8 × 10^{−06} | 2.231 776 1 × 10^{−05} |

xy^{2}z | 7.652 722 2 × 10^{−04} | 4.956 300 3 × 10^{−04} | −1.953 866 9 × 10^{−03} |

xyz^{2} | −3.053 939 8 × 10^{−05} | −1.977 890 0 × 10^{−05} | 7.797 215 1 × 10^{−05} |

xz^{3} | 4.838 899 2 × 10^{−05} | 3.133 922 4 × 10^{−05} | −1.235 451 2 × 10^{−04} |

y^{4} | −3.262 586 8 × 10^{−05} | −2.113 020 7 × 10^{−05} | 8.329 925 5 × 10^{−05} |

y^{3}z | 1.799 062 5 × 10^{−03} | 1.165 166 3 × 10^{−03} | −4.593 305 0 × 10^{−03} |

y^{2}z^{2} | 2.276 793 8 × 10^{−03} | 1.474 569 9 × 10^{−03} | −5.813 032 3 × 10^{−03} |

yz^{3} | 3.033 272 8 × 10^{−03} | 1.964 505 0 × 10^{−03} | −7.744 448 8 × 10^{−03} |

z^{4} | 5.742 310 1 × 10^{−03} | 3.719 018 2 × 10^{−03} | −1.466 107 0 × 10^{−02} |

dB1/dt | 9.486 847 6 × 10^{−01} | 9.444 878 4 × 10^{−07} | −3.723 349 0 × 10^{−06} |

dB2/dt | −1.808 827 7 × 10^{−05} | 9.486 715 8 × 10^{−01} | 4.618 237 2 × 10^{−05} |

dB2/dt | 4.246 885 2 × 10^{−05} | 2.750 503 4 × 10^{−05} | 9.485 748 7 × 10^{−01} |

Basis . | $x\u02d9$ . | $y\u02d9$ . | $z\u02d9$ . |
---|---|---|---|

1 | −3.901 872 6 × 10^{−06} | −2.527 055 3 × 10^{−06} | 9.962 128 3 × 10^{−06} |

x | −9.997 259 3 × 10^{−02} | 2.000 017 8 × 10^{+00} | −6.997 397 7 × 10^{−05} |

y | −1.999 946 3 × 10^{+00} | −9.996 519 8 × 10^{−02} | −1.371 972 5 × 10^{−04} |

z | 1.972 107 4 × 10^{−04} | 1.277 239 2 × 10^{−04} | −3.005 035 1 × 10^{−01} |

x^{2} | 4.195 535 5 × 10^{−05} | 2.717 246 7 × 10^{−05} | −1.071 189 8 × 10^{−04} |

xy | −2.942 892 2 × 10^{−06} | −1.905 969 8 × 10^{−06} | 7.513 692 0 × 10^{−06} |

xz | 1.225 505 7 × 10^{−04} | 7.937 011 3 × 10^{−05} | −3.128 919 4 × 10^{−04} |

y^{2} | 2.406 794 8 × 10^{−05} | 1.558 765 3 × 10^{−05} | −6.144 946 4 × 10^{−05} |

yz | 3.387 555 3 × 10^{−04} | 2.193 956 7 × 10^{−04} | −8.648 990 8 × 10^{−04} |

z^{2} | 4.021 428 3 × 10^{−03} | 2.604 485 8 × 10^{−03} | −1.026 737 4 × 10^{−02} |

x^{3} | −1.164 887 1 × 10^{−04} | −7.544 413 7 × 10^{−05} | 2.974 150 0 × 10^{−04} |

x^{2} y | −2.361 637 3 × 10^{−04} | −1.529 519 0 × 10^{−04} | 6.029 651 9 × 10^{−04} |

x^{2}z | −8.288 823 6 × 10^{−04} | −5.368 272 6 × 10^{−04} | 2.116 274 2 × 10^{−03} |

xy^{2} | −1.177 680 6 × 10^{−04} | −7.627 271 3 × 10^{−05} | 3.006 814 1 × 10^{−04} |

xyz | 3.517 902 6 × 10^{−05} | 2.278 376 4 × 10^{−05} | −8.981 788 9 × 10^{−05} |

xz^{2} | −1.828 761 2 × 10^{−03} | −1.184 400 7 × 10^{−03} | 4.669 130 7 × 10^{−03} |

y^{3} | −2.352 567 6 × 10^{−04} | −1.523 645 0 × 10^{−04} | 6.006 495 4 × 10^{−04} |

y^{2}z | −1.009 917 0 × 10^{−03} | −6.540 746 8 × 10^{−04} | 2.578 485 6 × 10^{−03} |

yz^{2} | −6.201 471 5 × 10^{−03} | −4.016 395 0 × 10^{−03} | 1.583 338 6 × 10^{−02} |

z^{3} | −1.242 848 7 × 10^{−02} | −8.049 333 6 × 10^{−03} | 3.173 198 9 × 10^{−02} |

x^{4} | −9.374 085 3 × 10^{−05} | −6.071 144 4 × 10^{−05} | 2.393 359 5 × 10^{−04} |

x^{3}y | −2.102 290 2 × 10^{−06} | −1.361 552 3 × 10^{−06} | 5.367 495 7 × 10^{−06} |

x^{3}z | 7.505 776 6 × 10^{−04} | 4.861 130 7 × 10^{−04} | −1.916 349 3 × 10^{−03} |

x^{2}y^{2} | −1.271 197 6 × 10^{−04} | −8.232 935 5 × 10^{−05} | 3.245 578 3 × 10^{−04} |

x^{2}yz | 1.801 066 1 × 10^{−03} | 1.166 463 9 × 10^{−03} | −4.598 420 6 × 10^{−03} |

x^{2}z^{2} | 2.218 039 3 × 10^{−03} | 1.436 517 4 × 10^{−03} | −5.663 022 3 × 10^{−03} |

xy^{3} | −8.741 210 4 × 10^{−06} | −5.661 261 8 × 10^{−06} | 2.231 776 1 × 10^{−05} |

xy^{2}z | 7.652 722 2 × 10^{−04} | 4.956 300 3 × 10^{−04} | −1.953 866 9 × 10^{−03} |

xyz^{2} | −3.053 939 8 × 10^{−05} | −1.977 890 0 × 10^{−05} | 7.797 215 1 × 10^{−05} |

xz^{3} | 4.838 899 2 × 10^{−05} | 3.133 922 4 × 10^{−05} | −1.235 451 2 × 10^{−04} |

y^{4} | −3.262 586 8 × 10^{−05} | −2.113 020 7 × 10^{−05} | 8.329 925 5 × 10^{−05} |

y^{3}z | 1.799 062 5 × 10^{−03} | 1.165 166 3 × 10^{−03} | −4.593 305 0 × 10^{−03} |

y^{2}z^{2} | 2.276 793 8 × 10^{−03} | 1.474 569 9 × 10^{−03} | −5.813 032 3 × 10^{−03} |

yz^{3} | 3.033 272 8 × 10^{−03} | 1.964 505 0 × 10^{−03} | −7.744 448 8 × 10^{−03} |

z^{4} | 5.742 310 1 × 10^{−03} | 3.719 018 2 × 10^{−03} | −1.466 107 0 × 10^{−02} |

dB1/dt | 9.486 847 6 × 10^{−01} | 9.444 878 4 × 10^{−07} | −3.723 349 0 × 10^{−06} |

dB2/dt | −1.808 827 7 × 10^{−05} | 9.486 715 8 × 10^{−01} | 4.618 237 2 × 10^{−05} |

dB2/dt | 4.246 885 2 × 10^{−05} | 2.750 503 4 × 10^{−05} | 9.485 748 7 × 10^{−01} |

We take cuboid $D=[\u22122,2]\xd7[\u22122,2]\xd7[0,1]$, with boundary $\u2202D$, containing a sub-boundary $\Gamma $ to be the surface $z=1$ (top boundary) or the surface $z=0$ (bottom boundary).

For the original system (26)–(28), the mean residence time $u$ (for solutions with initial points in $D$) and escape probability $p$ (for solutions exiting through a sub-boundary $\Gamma $) satisfy the elliptic partial differential Eqs. (4)–(5) and (6)–(8), respectively. Hence,

These equations can also be solved by a finite element method to get the mean residence time $uTrue$, the escape probability $pTrue$, and the average escape probability $PTrue$.

For the learned model (29)–(31), we also have similar elliptic partial differential equations for the mean residence time $uLearn$, the escape probability $pLearn$ and the average escape probability $PLearn$.

The mean residence time $u(x,y,z)$ of the learned model is shown in Fig. 8. It is barely distinguishable from the mean residence time for the original system, and we thus do not show the latter. Instead, we compute the error (using 20 000 uniformly distributed points in $D$) between the maximal values of mean residence time for the learned and original systems: ${error}(u)=max(x,y,z)\u2208D\Vert uLearn\u2212uTrue\Vert =3.7495\xd710\u22124$.

The escape probability $p(x,y,z)$ of the learned system is shown in Fig. 9. It is barely distinguishable from the escape probability for the original system, and we thus do not show the latter. Instead, we calculate the error between the average escape probability values of the learned system and original system to be ${error}(P)=\Vert PLearn\u2212PTrue\Vert =9.9478\xd710\u22127$ (escaping from the top boundary) and ${error}(P)=\Vert PLearn\u2212PTrue\Vert =1.0783\xd710\u22124$ (escaping from the bottom boundary).

### B. Lorenz system with random noises

We consider a stochastic Lorenz system,^{20}

where $\sigma =10$, $\beta =83$, $\rho =28$ are standard parameters and $B1(t)$, $B2(t)$, $B3(t)$ are independent scalar Brownian motions. We take $\u03f5=0.9$ in the following computations.

As before, a learned model is identified in the following form:

where $\Theta $ is a set of basis function consisting of polynomials in $(x,y,z)$ up to fourth order and time derivatives of three independent Brownian motions. Then, we solve a regression problem to determine the weights $\Xi =[\xi 1,\xi 2,\xi 3]$ (see Table IV).

Basis . | $x\u02d9$ . | $y\u02d9$ . | $z\u02d9$ . |
---|---|---|---|

1 | 7.519 769 2 × 10^{−03} | 3.352 837 9 × 10^{−03} | −1.369 622 8 × 10^{−02} |

x | −1.000 889 5 × 10^{+01} | 2.799 603 4 × 10^{+01} | 1.620 035 0 × 10^{−02} |

y | 1.000 469 0 × 10^{+01} | −9.979 090 7 × 10^{−01} | −8.541 396 7 × 10^{−03} |

z | −2.704 581 5 × 10^{−03} | −1.205 891 2 × 10^{−03} | −2.661 740 6 × 10^{+00} |

x^{2} | −4.563 336 2 × 10^{−04} | −2.034 653 7 × 10^{−04} | 8.311 490 8 × 10^{−04} |

xy | −6.632 591 2 × 10^{−04} | −2.957 272 1 × 10^{−04} | 1.001 208 0 × 10^{+00} |

xz | 1.013 984 2 × 10^{−03} | −9.995 479 0 × 10^{−01} | −1.846 833 1 × 10^{−03} |

y^{2} | 4.180 401 0 × 10^{−04} | 1.863 914 5 × 10^{−04} | −7.614 026 9 × 10^{−04} |

yz | −5.534 438 9 × 10^{−04} | −2.467 639 1 × 10^{−04} | 1.008 022 1 × 10^{−03} |

z^{2} | 3.086 039 3 × 10^{−04} | 1.375 971 7 × 10^{−04} | −5.620 797 2 × 10^{−04} |

x^{3} | 3.801 123 8 × 10^{−05} | 1.694 806 2 × 10^{−05} | −6.923 225 5 × 10^{−05} |

x^{2} y | −5.625 070 2 × 10^{−05} | −2.508 048 9 × 10^{−05} | 1.024 529 4 × 10^{−04} |

x^{2}z | 4.395 476 0 × 10^{−05} | 1.959 810 0 × 10^{−05} | −8.005 756 6 × 10^{−05} |

xy^{2} | 1.107 879 7 × 10^{−05} | 4.939 700 9 × 10^{−06} | −2.017 850 9 × 10^{−05} |

xyz | 1.607 076 0 × 10^{−05} | 7.165 466 2 × 10^{−06} | −2.927 068 4 × 10^{−05} |

xz^{2} | −3.438 403 0 × 10^{−05} | −1.533 080 0 × 10^{−05} | 6.262 579 5 × 10^{−05} |

y^{3} | 4.308 422 0 × 10^{−06} | 1.920 995 2 × 10^{−06} | −7.847 199 7 × 10^{−06} |

y^{2}z | −1.191 378 3 × 10^{−05} | −5.311 996 1 × 10^{−06} | 2.169 932 2 × 10^{−05} |

yz^{2} | 2.043 225 2 × 10^{−05} | 9.110 124 0 × 10^{−06} | −3.721 454 4 × 10^{−05} |

z^{3} | −1.364 567 7 × 10^{−05} | −6.084 195 1 × 10^{−06} | 2.485 372 9 × 10^{−05} |

x^{4} | 2.471 056 3 × 10^{−07} | 1.101 769 3 × 10^{−07} | −4.500 690 0 × 10^{−07} |

x^{3}y | −1.531 861 5 × 10^{−06} | −6.830 107 7 × 10^{−07} | 2.790 075 6 × 10^{−06} |

x^{3}z | −6.061 635 5 × 10^{−07} | −2.702 700 1 × 10^{−07} | 1.104 043 8 × 10^{−06} |

x^{2}y^{2} | 1.996 378 6 × 10^{−06} | 8.901 248 8 × 10^{−07} | −3.636 129 7 × 10^{−06} |

x^{2}yz | 1.158 815 1 × 10^{−06} | 5.166 806 5 × 10^{−07} | −2.110 622 8 × 10^{−06} |

x^{2}z^{2} | −1.034 200 2 × 10^{−06} | −4.611 186 3 × 10^{−07} | 1.883 653 8 × 10^{−06} |

xy^{3} | −7.478 363 0 × 10^{−07} | −3.334 376 1 × 10^{−07} | 1.362 081 2 × 10^{−06} |

xy^{2}z | −2.273 687 1 × 10^{−07} | −1.013 768 4 × 10^{−07} | 4.141 209 1 × 10^{−07} |

xyz^{2} | 2.575 782 3 × 10^{−07} | 1.148 463 5 × 10^{−07} | −4.691 434 1 × 10^{−07} |

xz^{3} | 3.373 253 2 × 10^{−07} | 1.504 031 7 × 10^{−07} | −6.143 917 9 × 10^{−07} |

y^{4} | −2.952 088 7 × 10^{−09} | −1.316 247 2 × 10^{−09} | 5.376 824 5 × 10^{−09} |

y^{3}z | −1.101 172 7 × 10^{−07} | −4.909 796 3 × 10^{−08} | 2.005 635 0 × 10^{−07} |

y^{2}z^{2} | −6.271 069 7 × 10^{−08} | −2.796 080 5 × 10^{−08} | 1.142 189 3 × 10^{−07} |

yz^{3} | −2.384 655 6 × 10^{−07} | −1.063 245 9 × 10^{−07} | 4.343 323 0 × 10^{−07} |

z^{4} | 2.034 981 4 × 10^{−07} | 9.073 367 2 × 10^{−08} | −3.706 439 4 × 10^{−07} |

dB1/dt | 9.458 149 5 × 10^{−01} | −1.278 908 2 × 10^{−03} | 5.224 296 3 × 10^{−03} |

dB2/dt | −3.127 910 5 × 10^{−03} | 9.472 886 6 × 10^{−01} | 5.697 059 9 × 10^{−03} |

dB3/dt | 3.494 294 4 × 10^{−03} | 1.558 000 3 × 10^{−03} | 9.423 189 2 × 10^{−01} |

Basis . | $x\u02d9$ . | $y\u02d9$ . | $z\u02d9$ . |
---|---|---|---|

1 | 7.519 769 2 × 10^{−03} | 3.352 837 9 × 10^{−03} | −1.369 622 8 × 10^{−02} |

x | −1.000 889 5 × 10^{+01} | 2.799 603 4 × 10^{+01} | 1.620 035 0 × 10^{−02} |

y | 1.000 469 0 × 10^{+01} | −9.979 090 7 × 10^{−01} | −8.541 396 7 × 10^{−03} |

z | −2.704 581 5 × 10^{−03} | −1.205 891 2 × 10^{−03} | −2.661 740 6 × 10^{+00} |

x^{2} | −4.563 336 2 × 10^{−04} | −2.034 653 7 × 10^{−04} | 8.311 490 8 × 10^{−04} |

xy | −6.632 591 2 × 10^{−04} | −2.957 272 1 × 10^{−04} | 1.001 208 0 × 10^{+00} |

xz | 1.013 984 2 × 10^{−03} | −9.995 479 0 × 10^{−01} | −1.846 833 1 × 10^{−03} |

y^{2} | 4.180 401 0 × 10^{−04} | 1.863 914 5 × 10^{−04} | −7.614 026 9 × 10^{−04} |

yz | −5.534 438 9 × 10^{−04} | −2.467 639 1 × 10^{−04} | 1.008 022 1 × 10^{−03} |

z^{2} | 3.086 039 3 × 10^{−04} | 1.375 971 7 × 10^{−04} | −5.620 797 2 × 10^{−04} |

x^{3} | 3.801 123 8 × 10^{−05} | 1.694 806 2 × 10^{−05} | −6.923 225 5 × 10^{−05} |

x^{2} y | −5.625 070 2 × 10^{−05} | −2.508 048 9 × 10^{−05} | 1.024 529 4 × 10^{−04} |

x^{2}z | 4.395 476 0 × 10^{−05} | 1.959 810 0 × 10^{−05} | −8.005 756 6 × 10^{−05} |

xy^{2} | 1.107 879 7 × 10^{−05} | 4.939 700 9 × 10^{−06} | −2.017 850 9 × 10^{−05} |

xyz | 1.607 076 0 × 10^{−05} | 7.165 466 2 × 10^{−06} | −2.927 068 4 × 10^{−05} |

xz^{2} | −3.438 403 0 × 10^{−05} | −1.533 080 0 × 10^{−05} | 6.262 579 5 × 10^{−05} |

y^{3} | 4.308 422 0 × 10^{−06} | 1.920 995 2 × 10^{−06} | −7.847 199 7 × 10^{−06} |

y^{2}z | −1.191 378 3 × 10^{−05} | −5.311 996 1 × 10^{−06} | 2.169 932 2 × 10^{−05} |

yz^{2} | 2.043 225 2 × 10^{−05} | 9.110 124 0 × 10^{−06} | −3.721 454 4 × 10^{−05} |

z^{3} | −1.364 567 7 × 10^{−05} | −6.084 195 1 × 10^{−06} | 2.485 372 9 × 10^{−05} |

x^{4} | 2.471 056 3 × 10^{−07} | 1.101 769 3 × 10^{−07} | −4.500 690 0 × 10^{−07} |

x^{3}y | −1.531 861 5 × 10^{−06} | −6.830 107 7 × 10^{−07} | 2.790 075 6 × 10^{−06} |

x^{3}z | −6.061 635 5 × 10^{−07} | −2.702 700 1 × 10^{−07} | 1.104 043 8 × 10^{−06} |

x^{2}y^{2} | 1.996 378 6 × 10^{−06} | 8.901 248 8 × 10^{−07} | −3.636 129 7 × 10^{−06} |

x^{2}yz | 1.158 815 1 × 10^{−06} | 5.166 806 5 × 10^{−07} | −2.110 622 8 × 10^{−06} |

x^{2}z^{2} | −1.034 200 2 × 10^{−06} | −4.611 186 3 × 10^{−07} | 1.883 653 8 × 10^{−06} |

xy^{3} | −7.478 363 0 × 10^{−07} | −3.334 376 1 × 10^{−07} | 1.362 081 2 × 10^{−06} |

xy^{2}z | −2.273 687 1 × 10^{−07} | −1.013 768 4 × 10^{−07} | 4.141 209 1 × 10^{−07} |

xyz^{2} | 2.575 782 3 × 10^{−07} | 1.148 463 5 × 10^{−07} | −4.691 434 1 × 10^{−07} |

xz^{3} | 3.373 253 2 × 10^{−07} | 1.504 031 7 × 10^{−07} | −6.143 917 9 × 10^{−07} |

y^{4} | −2.952 088 7 × 10^{−09} | −1.316 247 2 × 10^{−09} | 5.376 824 5 × 10^{−09} |

y^{3}z | −1.101 172 7 × 10^{−07} | −4.909 796 3 × 10^{−08} | 2.005 635 0 × 10^{−07} |

y^{2}z^{2} | −6.271 069 7 × 10^{−08} | −2.796 080 5 × 10^{−08} | 1.142 189 3 × 10^{−07} |

yz^{3} | −2.384 655 6 × 10^{−07} | −1.063 245 9 × 10^{−07} | 4.343 323 0 × 10^{−07} |

z^{4} | 2.034 981 4 × 10^{−07} | 9.073 367 2 × 10^{−08} | −3.706 439 4 × 10^{−07} |

dB1/dt | 9.458 149 5 × 10^{−01} | −1.278 908 2 × 10^{−03} | 5.224 296 3 × 10^{−03} |

dB2/dt | −3.127 910 5 × 10^{−03} | 9.472 886 6 × 10^{−01} | 5.697 059 9 × 10^{−03} |

dB3/dt | 3.494 294 4 × 10^{−03} | 1.558 000 3 × 10^{−03} | 9.423 189 2 × 10^{−01} |

A single trajectory of this stochastic system with initial condition $(x,y,z)T=(\u22128,7,27)T$ is shown in Fig. 10 (left), together with the same trajectory captured by the learned model in Fig. 10 (right).

For the stochastic Lorenz system, we take a cuboid $D$ to be either $D1=[\u22129,\u22128]\xd7[\u22129,\u22128]\xd7[27,28]$ containing the left saddle as the left residence region or $D2=[8,9]\xd7[8,9]\xd7[27,28]$ containing the right saddle as the right residence region. Then, mean residence time and escape probability satisfy the following elliptic partial differential equations, respectively:

The finite element numerical solutions to these partial differential equations are denoted by $uTrue$ and $pTrue$.

For the learned Lorenz system (35)–(37), with basis and coefficients in Table IV, we can also set up the partial differential equations for the learned mean residence time $uLearn$, and the learned escape probability $pLearn$ together with the average escape probability $PLearn$.

The learned mean residence time $u$ of the learned Lorenz system (35)–(37) is shown in Fig. 11. The maximal mean residence time in the left region $D1$ is $uLearn(D1)=0.1238$, while the maximal mean residence time in the right region $D2$ is $uLearn(D2)=0.1236$. Then, we compute the error of maximal mean residence time between the learned and the original systems, ${error}(u(D1))=max(x,y,z)\u2208D1\Vert uLearn\u2212uTrue\Vert =1.0023\xd710\u22124$ and ${error}(u(D2))=max(x,y,z)\u2208D2\Vert uLearn\u2212uTrue\Vert =1.2284\xd710\u22125$.

We take $\Gamma $ to be each surface of the region. The learned escape probability $pLearn$ of particles exiting through each surface $\Gamma $ of left and right regions are shown in Figs. 12 and 13, respectively, and the average escape probability of left region is $PLearn(L1)=0.1065$ (from $x=\u22129$), $PLearn(L2)=0.1167$ (from $x=\u22128$), $PLearn(L3)=0.0955$ (from $y=\u22129$), $PLearn(L4)=0.3485$ (from $y=\u22128$), $PLearn(L5)=0.2343$ (from $z=27$), $PLearn(L6)=0.1244$ (from $z=28$). Similarly, the average escape probability of the right region is $PLearn(R1)=0.1161$ (from $x=8$), $PLearn(R2)=0.1061$ (from $x=9$), $PLearn(R3)=0.3469$ (from $y=8$), $PLearn(R4)=0.0947$ (from $y=9$), $PLearn(R5)=0.2367$ (from $z=27$), $PLearn(R6)=0.1264$ (from $z=28$).

Moreover, for example, we compute the error of the average escape probability between the learned and the original systems only for two cases: ${error}(P)(L4)=\Vert PLearn(L4)\u2212PTrue(L4)\Vert =6.5838\xd710\u22123$ and ${error}(P)(R3)=\Vert PLearn(R3)\u2212PTrue(R3)\Vert =6.431\xd710\u22123$ (L4).

## IV. CONCLUSION

In summary, we have demonstrated a new data-driven method to determine dynamical quantities, mean residence time and escape probability, based on sample path data of complex systems, as long as these systems are modeled by stochastic differential equations.^{2,12,14,19,22,23}

A novelty of this method is that it contains noise terms (in terms of Brownian motions) in the basis, in order to discover the governing stochastic differential equation. With the governing stochastic differential equation as the model for the sample path data, we can then compute the mean residence time and escape probability, which are deterministic quantities that carry significant dynamical information.

This method combines machine learning tools (i.e., extracting governing stochastic differential equations from data) and stochastic dynamical systems techniques (i.e., understanding random phenomena with deterministic quantities) to provide an efficient approach in examining stochastic dynamics from data.

To demonstrate that our method is effective, we have illustrated it on several stochastic systems with Brownian motions. It is expected to be applicable to data sets from other stochastic dynamical systems.

## ACKNOWLEDGMENTS

We would like to thank Min Dai and Jian Ren for helpful comments. The first author was supported by the China Scholarship Council under No. 201808220008, and the second author was supported by the China Scholarship Council under No. 201808220009. The second author was also supported by the Jilin Province Industrial Technology Research and Development Special Project under No. 2016C079. The third author was partly supported by the National Science Foundation (NSF) Grant No. 1620449.