Mean escape time for randomly switching narrow gates in a cellular flow

TL;DR

This study analyzes how the mean escape time of particles through narrow gates in a cellular flow is affected by gate switching, flow dynamics, and diffusion, revealing optimal conditions for escape efficiency.

Contribution

We derive coupled PDEs for the mean escape time in a switching gate system and quantify the effects of system parameters through numerical simulations.

Findings

Mean escape time decreases with switching rate, flow velocity, and diffusion coefficient.

Escape time is minimized when gates are evenly spaced.

For small gate size, escape time is nearly independent of gate size.

Abstract

The escape of particles through a narrow absorbing gate in confined domains is a abundant phenomenon in various systems in physics, chemistry and molecular biophysics. We consider the narrow escape problem in a cellular flow when the two gates randomly switch between different states with a switching rate k between the two gates. After briefly deriving the coupled partial differential equations for the escape time through two gates, we compute the mean escape time for particles escaping from the gates with different initial states. By numerical simulation under nonuniform boundary conditions, we quantify how narrow escape time is affected by the switching rate k between the two gates, arc length s between two gates, angular velocity w of the cellular flow and diffusion coefficient D. We reveal that the mean escape time decreases with the switching rate k between the two gates, angular velocity w and diffusion coefficient D for fixed arc length, but takes the minimum when the two gates are evenly separated on the boundary for any given switching rate k between the two gates. In particular, we find that when the arc length size {\epsilon} for the gates is sufficiently small, the average narrow escape time is approximately independent of the gate arc length size. We further indicate combinations of system parameters (regions located in the parameter space) such that the mean escape time is the longest or shortest. Our findings provide mathematical understanding for phenomena such as how ions select ion channels and how chemicals leak in annulus ring containers, when drift vector fields are present.

Equations

$$A=\frac{-w{r}_1^2}{r_2^2-r_1^2},$$

$$B=\frac{w{r}_1^2{r}_2^2}{r_2^2-r_1^2}.$$

$$b_1(x,y)=-y\frac{-w{r}_1^2(x^2+y^2)+w{r}_1^2{r}_2^2}{\sqrt{x^2+y^2}(r_2^2-r_1^2)},$$

$$b_2(x,y)=x\frac{-w{r}_1^2(x^2+y^2)+w{r}_1^2{r}_2^2}{\sqrt{x^2+y^2}(r_2^2-r_1^2)}.$$

$$\{\begin{array}{ll}d{x}_t& =b_1(x_t,y_t)dt+\sqrt{2D}d{B}_t^1-n_1(x_t,y_t)d{L}_t,\\ d{y}_t& =b_2(x_t,y_t)dt+\sqrt{2D}d{B}_t^2-n_2(x_t,y_t)d{L}_t,\end{array}$$

$$L_t={\int }_0^t{\mathbb{I}}_{(x,y)\in \partial \Omega }d{L}_s.$$

$$N_{t}=\left\{\begin{array}[]{ll}1,&\hbox{if $T_{2j}\leq t<T_{2j+1}$ for some $j$,}\\ 2,&\hbox{if $T_{2j+1}\leq t<T_{2j+2}$ for some $j$,}\end{array}\right.$$

$$-\left(\begin{array}{c}1\\ 1\end{array}\right)=k\left(\begin{array}{cc}-1& 1\\ 1& -1\end{array}\right)\left(\begin{array}{c}{u}_1(x,y)\\ u_2(x,y)\end{array}\right)+D\Delta \left(\begin{array}{c}{u}_1(x,y)\\ u_2(x,y)\end{array}\right)+\left(\begin{array}{c}b\cdot \nabla u_1(x,y)\\ b\cdot \nabla u_2(x,y)\end{array}\right)$$

$${\partial }_n{u}_j{|}_{\partial {\Omega }_c}=0,j=1,2,$$

$$u_1{|}_{\partial {\Omega }_1}=0,u_2{|}_{\partial {\Omega }_2}=0,$$

$${\partial }_n{u}_1{|}_{\partial {\Omega }_2}=0,{\partial }_n{u}_2{|}_{\partial {\Omega }_1}=0,$$

$$\mathbb{E}[u({\mathbf{X}}_{t\wedge \mathcal{T}},N_{t\wedge \mathcal{T}})-u({\mathbf{X}}_{s\wedge \mathcal{T}},N_{s\wedge \mathcal{T}})-{\int }_{s\wedge \mathcal{T}}^{t\wedge \mathcal{T}}\mathcal{A}u({\mathbf{X}}_r,N_r)dr|{\mathcal{F}}_{s\wedge \mathcal{T}}]$$

$$\mathbb{E}[u({\mathbf{X}}_{t\wedge \mathcal{T}},N_{t\wedge \mathcal{T}})-u({\mathbf{X}}_{s\wedge \mathcal{T}},N_{s\wedge \mathcal{T}})+(t\wedge \mathcal{T}-s\wedge \mathcal{T})|{\mathcal{F}}_{s\wedge \mathcal{T}}]=0,$$

$$\mathbb{E}[u({\mathbf{X}}_{\mathcal{T}},N_{\mathcal{T}})-u({\mathbf{X}}_0,N_0)+\mathcal{T}|{\mathcal{F}}_0]=0.$$

$$\mathbb{E}[\mathcal{T}|{\mathcal{F}}_0]=\mathbb{E}[u({\mathbf{X}}_0,N_0)|{\mathcal{F}}_0].$$

$${\mathbb{E}}_{k,\mathbf{x}}\mathcal{T}=\mathbb{E}u({\mathbf{x}}_0,N_0)=\mathbb{E}u(\mathbf{x},k)=u_k(\mathbf{x}).$$

$${\partial }_{\nu }u={\partial }_{x}u\cos \theta +{\partial }_{y}u\sin \theta =0,$$

$${\partial }_{x}u={\partial }_{r}u{\partial }_{x}r+{\partial }_{\theta }u{\partial }_x\theta =\cos \theta {\partial }_{r}u-\frac{\sin \theta }{r}{\partial }_{\theta }u,$$

$${\partial }_{y}u={\partial }_{r}u{\partial }_{y}r+{\partial }_{\theta }u{\partial }_y\theta =\sin \theta {\partial }_{r}u+\frac{\cos \theta }{r}{\partial }_{\theta }u,$$

$$\Delta u={\partial }_{xx}u+{\partial }_{yy}u={\partial }_{rr}u+\frac{1}{r^2}{\partial }_{\theta \theta }u+\frac{1}{r}{\partial }_{r}u.$$

$$\frac{1}{\pi (r_2^2-r_1^2)}\int {\int }_{{\Omega }^{\prime }}u(r,\theta )drd\theta ,$$

Full text

Key words: Narrow escape; Brownian motion; small gate; mixed boundary conditions; biophysical modeling.

1 Introduction

The narrow escape problems have attracted a lot of attention in recent years, due to its significant relevance in cellular and molecular biology, and chemical physics [1, 2, 3, 4, 5, 6]. The heart of a narrow escape problem is to compute the narrow escape time (NET), namely the mean time for a particle starting in the domain before exiting through a narrow gate on the boundary.

In the biophysics context, when an ion searches for an open ion channel within the membrane, neurotransmitter receptors transmit signal within a synapse of a neuron, a reactive particle to activate a given protein, or newly transcribed mRNA transports from the nucleus to the cytoplasm through nuclear pores [7, 8, 9]. These examples have something in common: ‘Particles’ are all confined to a bounded domain with small exits (or holes, or gates) on the boundary of the domain. Meanwhile, in order to fulfill their biological function, particles must exit from the small gates on the boundary. For instance, viruses must enter the cell nucleus through small nanopores in order to replicate [10].

The narrow escape problems have been studied via various stochastic models recently. Most of the scenarios considered are two-dimensional [11, 12, 13, 14, 15], and some are one-dimensional [16, 17] or three-dimensional [17, 18, 19]. In addition, some authors consider the boundary with one gate [7, 20, 21, 22], with two gates[16, 20], while others multiple gates [10, 20, 23].

Several papers [14, 17, 18, 24, 25] have been devoted to the asymptotic approximate expressions for the narrow escape time. Since the mixed boundary conditions, there is no exact solution for a narrow escape problem. The narrow escape problems are mostly studied analytically [12, 14, 15, 16, 18, 20, 21, 23, 26, 27, 28, 29, 30, 31, 32]. For example, a general exact formula for the mean first passage time from a fixed point inside a planar domain to an escape region on its boundary is attained in [12]; high-order asymptotic formulas for the mean first passage time are derived in [14]; approximations for the average mean first passage time are found in[18]; the upper and lower bounds of the mean first passage time for mortal walkers are derived in[26]. Moreover, a spectral approach to derive an exact formula for the mean exit time of a particle through a hole on the boundary is developed in[27]; partial differential equation and probabilistic methods are applied to solve escape problems in [28]; a generalized Kramers formula for the mean escape time through a narrow window is obtained in[30]; and boundary homogenization is used when the boundary contains non-overlapping identical absorbing arcs in [31].

A few authors numerically studied the narrow escape problems, including Brownian simulations [17] and Monte Carlo simulations [33]. Think of the narrow escape problem the other way around, Agranov and Meerson [22] developed a formalism to evaluate the nonescape probability of a gas of diffusing particles.

Most works approximately estimate the narrow escape time of particles under Brownian diffusion (i.e., no drift or no external vector field). It turns out that the escape time is described by a mixed Dirichlet-Neumann boundary value problem for the Laplace operator in a bounded domain [1, 34]. Very recently, Lagache and Holcman [10] considered the narrow escape problem for particles in a constant drift, in order to analyze viral entry into the cell nucleus.

In this paper, we consider a two dimensional model of narrow escape problem for particles in a unit annulus, with two randomly switching gates on the outer boundary, under a vector field (i.e., a cellular flow). The circular annulus may represent the cell cytoplasm, with outer radius $r_{2}$ and inner radius $r_{1}$. The two dimensional model may represent flat culture cells. Since the adhesion to the substrate, they stay flat. The thickness can be neglected in this idealized model.

This paper is organized as follows. In Section 2, we first describe a stochastic narrow escape model, then briefly drive the coupled partial differential equations to be satisfied by the mean escape time through two gates. In Section 3, we present numerical experiments and analysis. We end the paper with a discussion in Section 4.

2 Materials and methods

2.1 Stochastic narrow escape model

We consider particles moving in an annular domain with inner radius $r_{1}$ and outer radius $r_{2}$, in a cellular flow. The steady cellular flow [35] is induced by the rotation of the inner boundary with angular velocity $w$, and is expressed in polar coordinates as $(0,Ar+\frac{B}{r})$ with

[TABLE]

By the transformation of polar coordinates, $x=r\cos\theta$, $y=r\sin\theta$ with $r=\sqrt{x^{2}+y^{2}}$, the cellular flow in the Cartesian coordinates is

[TABLE]

With the driving cellular flow $b(x,y):=(b_{1}(x,y),b_{2}(x,y))$, we consider the two dimensional domain $\Omega$, whose boundary is decomposed as $\partial\Omega=\partial\Omega_{c}\bigcup\partial\Omega_{1}\bigcup\partial\Omega_{2}$. The boundary has two small gates at the arcs $\partial\Omega_{1}$ and $\partial\Omega_{2}$, together with the remaining part $\Omega_{c}$. A particle’s trajectory is reflected on $\Omega_{c}$, but absorbed at gates $\Omega_{1}$ or $\Omega_{2}$. Each gate is a small absorbing arc. We consider the following stochastic narrow escape model, i.e., the equations of motion for a particle

[TABLE]

where $(b_{1}(x,y),b_{2}(x,y))$ is a drift vector field (a cellular flow), $(B_{t})_{t\geq 0}=(B_{t}^{1},B_{t}^{2})$ is a standard Brownian motion, $n=(n_{1},n_{2})$ is the outward unit normal vector on the boundary, $(L_{t})_{t\geq 0}$ is a continuous nondecreasing process (with $L_{0}=0$) which increases only when $(x,y)$ is on the boundary $\partial\Omega$:

[TABLE]

Assume that the two gates at the arcs $\partial\Omega_{1}$ and $\partial\Omega_{2}$ open alternatively [11, 15] in terms of a telegraph process $(N_{t})_{t\geq 0}$ with parameter $k$, where $k$ is the rate of switching between the two gates. When the particle in the cytoplasm hits the open gate it escapes, if not it will be reflected. The telegraph process $(N_{t})_{t\geq 0}$ with parameter $k$ is a Markov process which takes on only two values, $1$ and $2$. The different values of Markov process corresponding to different states, different function of cell [36]. The process can be described as follows

[TABLE]

where $T_{0}=0$, $T_{j}=\sum_{l=1}^{j}\tau_{l}$, and $\tau_{l}$ is a sequence of independent and identically distributed exponential random variables with parameter $k>0$. That is, $\mathbb{E}[\tau_{l}]=\frac{1}{k}$. Denote $\textbf{X}_{t}=(x_{t},y_{t})$ and $\textbf{x}=(x,y)$. We denote $\mathcal{T}$ as the stopping time corresponding to the escape of the particle: $\mathcal{T}=\mathcal{T}_{1}\wedge\mathcal{T}_{2}$, $\mathcal{T}_{j}=\inf\{t\geq 0,N_{t}=j\ and\ \mathbf{X}_{t}\in\partial\Omega_{j}\}$, j=1, 2.

2.2 Narrow escape time with drift

Our goal is to compute the expected quantities $u_{1}(\textbf{x})=\mathbb{E}_{1,\textbf{x}}[\mathcal{T}]$ and $u_{2}(\textbf{x})=\mathbb{E}_{2,\textbf{x}}[\mathcal{T}]$ of the stopping time given the initial state is $N_{0}=j$ and $\mathbf{X}_{0}=\mathbf{x}$, j=1, 2.

More specifically, $u_{1}(x,y)$ is the mean time that the particle starting at $(x,y)$ in the domain escape from any gate with the first gate open initially. Similarly, $u_{2}(x,y)$ is the mean time that the particle starting at $(x,y)$ in the domain escape from any gate with the second gate open initially. We now prove the following theorem.

Theorem 1**.**

*The mean narrow escape times $u_{1},u_{2}$ satisfy the following coupled partial differential equations

[TABLE]

with nonuniform boundary conditions:

[TABLE]

[TABLE]

[TABLE]

where $\partial_{n}$ represents the normal derivative and $n$ is outward unit normal vector on the outer boundary.

Proof.

Note that the right-hand side of the (7) is the infinitesimal generator $\mathcal{A}$ of the process $(\mathbf{X}_{t},N_{t})_{t\geq 0}$. Denoting $u(\mathbf{x},j)=u_{j}(\mathbf{x})$, by Itô formula [37], we have the following relation, for $s<t$,

[TABLE]

where $(\mathcal{F}_{t})_{t\geq 0}$ is the filtration of the process $(\mathbf{X}_{t},N_{t})_{t\geq 0}$, and $L_{t}$ is defined as above. Neumann boundary conditions (8) and (10) show us that $\partial_{\nu}u(\mathbf{X}_{r},N_{r})=0$ if $r<\mathcal{T}$. By the given condition, we know $\mathcal{A}u(\mathbf{X}_{r},N_{r})=-1$. This leads to

[TABLE]

which means the process $(u(\mathbf{X}_{t\wedge\mathcal{T}},N_{t\wedge\mathcal{T}})+t\wedge\mathcal{T})_{t\geq 0}$ is a martingale. Then applying this equality with $s=0$ and letting $t\rightarrow\infty$, we get

[TABLE]

From the Dirichlet boundary conditions in (9), we have $u(\mathbf{X}_{\mathcal{T}},N_{\mathcal{T}})=0$. Thus we obtain

[TABLE]

Finally, we take the expectation of both sides of the equation (12),

[TABLE]

This gives the desired result when the initial distribution is such that $N_{0}=k$ and $\mathbf{X}_{0}=\mathbf{x}$. ∎

Note that Ammari et al. [11] showed this theorem in the case of no drift.

In the polar coordinates, the NET $u(r,\theta)$ is the solution of the preceding boundary value problem given the initial point $(r,\theta)$. It is a reflected diffusion process to the absorbing boundary $\partial\Omega_{j},(j=1,2)$. In the context of polar coordinates: $x=r\cos\theta$, $y=r\sin\theta$, we obtain that

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

We use a finite difference method to solve the Dirichlet-Neumann boundary value problem (7)-(10) in the polar coordinates.

3 Results and Discussion

We now describe our numerical simulation results of the narrow escape time (NET). The Dirichlet-Neumann boundary value problem in Theorem 1 is solved by a finite difference method. The goal here is to present the results for the averaged NET from the numerical simulation, then conduct comparisons for different parameters.

The average NET is the average escape time that all the particles in the annular domain take to escape the annular domain from any gate. It is also called average mean first passage time [38, 39]. In the context of polar coordinates, the average NET is given by

[TABLE]

where $\Omega^{\prime}=[r_{1},r_{2}]\times[0,2\pi]$.

First, we show the relationship of average NET with switching rate $k$ between the two gates and arc length (‘distance’) $s$ between two gates in Fig 2 and Fig 3.

In Fig 2, the average NET is plotted as a function of switching rate between the two gates $k$. We compare average NET for different arc length $s$ between two gates. when the arc length between two gates $0<s\leq\pi$, the average NET decreases with the increase of arc length $s$ between two gates , while for $\pi<s<2\pi$, the average NET is longer than that of $s=\pi$. To understand the reason for this at $s=\pi$, note that the two gates are evenly distributed on the cell boundary, so the particle may choose the nearest side to escape quickly; but for $0<s<\pi$ or $\pi<s<2\pi$, the two gates are on the same side of the semicircle, thus the particle in the other semicircle has to cross over to the other side then escape. Fig 3 plots the effects of arc length $s$ between two gates on average NET for various switching rate $k$ between the two gates. As we explained above, the average NET is smallest when arc length $s$ between two gates is around $\pi$. Meanwhile, Fig 2 and Fig 3 both suggest that the average NET decrease as the switching rate $k$ between the two gates increases, this is consistent with known theoretical result [11].

Each curve in Fig 4 predicts an decrease in the average NET with angular velocity $w$. We plot here for switching rate $k=0.1$ between the two gates, for $0<k<1$, it has similar trends.

Fig 5 demonstrates that the effects of diffusion coefficient $D$ on the average NET. We see that the average NET has a monotonic behavior with the diffusion coefficient $D$. The average NET declines when $D$ is small, then when average NET crosses the inflection point, it changes more modestly. In this case, small diffusion coefficient $D$ is not good for particles to escape quickly.

Fig 6 presents the combined effects of angular velocity $w$ and switching rate $k$ between the two gates on average NET. In both Fig 6(a) and Fig 6(b), the red region presents the larger average NET domain, while the blue region represents the smaller parts. The smaller average NET values occur when $(w,k)$ is in the blue domains at the top. Different initial states have very little influence on average NET.

In Fig 7, we display average NET as a function of arc length $s$ between two gates and switching rate $k$ between the two gates. Interestingly, we find that the average NET $u_{1}$ is larger than $u_{2}$ in many domains. If the initial state with the first gate open, escape from the first gate takes longer time than from the second gate. Similarly, we show the effects of arc length $s$ between two gates and angular velocity $w$ on average NET in Fig 8. We find a counterintuitive phenomenon, when the initial state is that the second gate is open, escape from the second gate is faster than from the first gate.

Fig 9 shows the effects of the gate arc length size $\varepsilon$ on average NET. As expected, the average NET decreases with the switching rate $k$ between the two gates, but when $\varepsilon$ is sufficiently small, the average NET is independent of gate arc length size $\varepsilon$.

4 Conclusion

In this work , we have analysed the average narrow escape time through two randomly switching gates in an annular domain, under a drift or a vector field, while most existing works are for particles following Brownian diffusion (i.e., no drift). We have conducted numerical experiments to reveal escape time’s dependence on angular velocity $w$ for the background cellular flow, switching rate $k$ between the two gates, diffusion coefficient $D$, arc length distance $s$ between two gates, and gate arc length size $\varepsilon$.

We have found that the mean escape time decreases with the switching rate $k$ between the two gates, angular velocity $w$ and diffusion coefficient $D$ for fixed arc length, but takes the minimum when the two gates are evenly separated on the boundary for given switching rate $k$ between the two gates. In particular, we have verified that when the arc length size $\varepsilon$ for the gates is sufficiently small, the average narrow escape time is independent of the gate arc length size. Our approach selects combinations of system parameters (located in the parameter space) such that the mean escape time is the longest or shortest.

The narrow escape problems are relevant in certain crucial processes in chemical physics and cellular biology. Our research sheds some lights on the narrow escape time in a vector field, when system parameters vary.

Acknowledgements

We would like to thank Xi Chen, Xiaoli Chen and Xiaofan Li for helpful discussions for numerical schemes. This work was partly supported by the NSF, USA grant 1620449, and the NSFC, China grants 11531006 and 11771449.