First passage time distribution in heterogeneity controlled kinetics: going beyond the mean first passage time

The first passage is a generic concept for quantifying when a random quantity such as the position of a diffusing molecule or the value of a stock crosses a preset threshold (target) for the first time. The last decade saw an enlightening series of new results focusing mostly on the so-called mean and global first passage time (MFPT and GFPT, respectively) of such processes. Here we push the understanding of first passage processes one step further. For a simple heterogeneous system we derive rigorously the complete distribution of first passage times (FPTs). Our results demonstrate that the typical FPT significantly differs from the MFPT, which corresponds to the long time behaviour of the FPT distribution. Conversely, the short time behaviour is shown to correspond to trajectories connecting directly from the initial value to the target. Remarkably, we reveal a previously overlooked third characteristic time scale of the first passage dynamics mirroring brief excursion away from the target.

Scientific RepoRts | 6:20349 | DOI: 10.1038/srep20349 to describe the FPT statistic [17][18][19][20] . Quite generally, the first passage statistic of Brownian motion has the generic asymptotic behaviour 18 where μ is the so-called persistence exponent 13 and a and b are dimension and geometry specific parameters. Eq. (1a) encodes the fact that it takes a finite minimum time to reach the target followed by a power law decay of FPTs on a time scale on which the searcher does not yet feel the presence of the boundary. In addition, Eq. (1b) states that the searcher will eventually find the target in a finite system of linear dimension R on a time scale up to ( ) b R f r a d a 2 , where ( ) f x d is a dimension dependent function, which diverges as the target radius r a goes to zero in two and three dimensions 18 .
The above results hold for translation invariant systems. Yet, numerous real systems such as biological cells are spatially heterogeneous and therefore display fundamentally different dynamics [21][22][23][24] . Various aspects of diffusion in heterogeneous media have already been addressed [25][26][27][28][29][30][31] but the rôle of spatial heterogeneity in the FPT statistics beyond the MFPT 32 remains elusive. Moreover, the results in ref. 32 suggest that the MFPT in (hyper)spherically symmetric domains is apparently independent of indirect trajectories, these are those that interact with the confining boundary, in contrast to direct trajectories, that head swiftly to the target.
Here we present exact results for the full FPT statistics in a simple heterogeneous model system. Based on a rigorous asymptotic analysis of the FPT statistics of Brownian motion in a confined spherically symmetric domain with a piece wise constant diffusion coefficient, see Fig. 1 and ref. 32, we here demonstrate the emergence of a new time scale in the FPT dynamics, which is controlled by the spatial heterogeneity. More precisely, we prove that the intermediate time power law asymptotics in Eq. (1a) breaks down in a sufficiently heterogeneous medium. For such heterogeneity-controlled kinetics we derive analytical asymptotic results for the short, intermediate, and long time FPT statistics for an arbitrary degree of heterogeneity. We also quantify the most likely (typical) FPT and the width of the FPT distribution. We demonstrate that the MFPT is dominated by long and unlikely indirect trajectories, while the overall relative contribution to the MFPT of the latter remains coupled to the most likely, direct trajectories. Finally, we discuss the implications of our results for more general systems.

Results
System setup and general result. We consider a spherically symmetric and potential free system with a perfectly absorbing central target of finite radius r a and a perfectly reflecting boundary at radius > R r a 32 . The system is in contact with a heat bath at constant and uniform temperature T. The particle experiences a space dependent friction Γ( ) r originating from spatial variations in the long range hydrodynamic coupling to the motion of the medium 33 . We focus on the high friction limit corresponding to overdamped motion and assume that the particle diffuses with the isotropic position dependent diffusion coefficient ( . Simultaneously, the diffusing particle experiences the fluctuation induced thermal drift ( ) ∇Γ( ) F r k T r B ensuring thermodynamic consistency in the sense that ( ) D r has a purely stochastic origin and does not reflect any heterogeneity in the entropic potential of mean force 34 . More precisely, in the absence of the target the system relaxes to the correct Boltzmann-Gibbs equilibrium-a spatially uniform probability. The thermodynamically consistent theory of diffusion in inhomogeneous media 34 corresponds to the so-called kinetic interpretation of the underlying multiplicative-noise Langevin equation; see, for instance, refs 35,36. We are interested in the evolution of the probability density function ( , ) P r t r 0 in dependence on the particle radius r at time t after starting from the initial radius r 0 at = t 0. Due to the symmetry of the system the angular Figure 1. Schematic of the model system with absorbing, finite central target of radius r a . The concentric shell of radius r i separates regions of different diffusivities D 1 and D 2 . The initial radius of the particle (dashed line) is r 0 , and the outer shell at radius R is reflective. The black and red lines denote direct trajectories and the green line denotes an indirect trajectory, see text for details.
Scientific RepoRts | 6:20349 | DOI: 10.1038/srep20349 co-ordinate is not of interest, and we average over the space angle. The diffusion equation governing the radial probability density function ( , ) P r t r 0 is then given by In our analysis we consider the particular case of a piece wise constant diffusion coefficient of magnitude ( ) = D r D 1 for ≤ r r i and ( ) = D r D 2 otherwise-see also Fig. 1. Physically, this form should be viewed as an ideal limit of a two phase system with a sharp interface. A similar limit of an infinitely sharp interface can be taken for the Langevin equation of the process, as well.
The exact solution for the Laplace transform ( , )  P r s r 0 was derived in ref. 32. From this result the Laplace transform of the FPT density is obtained from the corresponding probability flux into the target, π ℘( ) = × ∂ ( , We introduce dimensionless variables = / x r R for the particle position, for the target radius, = / where ν I and ν K denote the modified Bessel functions of the first and second kind, respectively 37 . In order to allow for a meaningful comparison of the FPT kinetics at various degrees ϕ of heterogeneity we introduce a constraint on the conservation of the spatially averaged diffusion coefficient, . Note, however, that the constraint does not introduce any additional information which would affect the qualitative picture of our results. The general case for arbitrary D 1 and D 2 is recovered trivially by treating Q as an independent parameter or by replacing / → Q is the starting point of our asymptotic analysis. In addition, in order to validate the analytical results we numerically invert ℘( )  s using the fixed Talbot method 38 .

Short time asymptotic.
Starting from the central result (3a) the short time behaviour of the FPT distribution ℘( ) t is obtained from the asymptotic behaviour of the respective modified Bessel functions for large argument (large Laplace variable) 37 . First it can be shown for ,  z 1 , respectively. Combined with Eq. (3a) we obtain a limiting expression for ℘( ) t which can be inverted exactly-by performing the contour integral within the domain of validity of Eq. (4)-leading us to the Lévy-Smirnov density and θ ϕ ϕ ( ) = + / 1 1 otherwise. It is easy to see that all x i -dependent terms vanish for a homogeneous systems with ϕ = 1. Note that Eq. (5a) obeys the generic behaviour given in Eq. (1a) with a persistence exponent μ = / 1 2 and . Moreover, Eqs (5a) and (5b) have an intuitive physical meaning: they show that the earliest FPTs will be observed on a time scale on which the particle diffuses over a distance corresponding to the initial separation to the target, with the respective diffusion coefficients. Moreover, we note that the FPT density to a finite radius x a for three dimensional radial Brownian motion in our heterogeneous system obeys the − / t 3 2 scaling Sparre Andersen theorem which needs to generally hold for one dimensional Markov processes with symmetric jumps 39,40 In other words, the first passage behaviour of particles released inside the interface radius x i is dominated by trajectories, which head straight for the target and do not venture into the outer part of the system. These are the direct trajectories introduced in ref. 32, see also below.
The FPT densities for various degrees of heterogeneity ϕ are shown in Figs 2 and 3, corresponding to initial positions x 0 in the inner and outer regions, respectively. Note the excellent agreement between the exact numerical result for ℘( ) t and the short-time asymptotics in Eq. (5a) as seen from comparison of the symbols with the black lines.
Long time asymptotic. The long time asymptotic behaviour of the FPT distribution ℘( ) t is obtained from Eq. (3a) by using the expansions of the respective modified Bessel functions for small argument (small Laplace variable s) 37 . Here we strictly note that corrections to leading order terms in ν I and ν K must be retained in order to obtain the correct behaviour of  ν and ν  for small s. We find that where t x x 0 0 a denotes the MFPT from x 0 to x a in a homogeneous sphere with unit radius and with unit diffusion coefficient, x a a 0 2 2 a As before, all x i -dependent terms vanish for a homogeneous system. Eq. (7) has the generic form of the exponential long time tail (1b) of the FPT distribution in a finite system. From Eq. (7) we identify the inverse of the characteristic time as = ( ) Moreover, Eq. (7) has an intuitive meaning: it demonstrates that first FPTs are exponentially unlikely beyond a time scale corresponding to the MFPT to arrive at the target from the external boundary. In turn, the long time exponential region evidently corresponds to trajectories, which are reflected from the external boundary. Therefore, for = x 1 0 we can no longer distinguish between direct and indirect trajectories. Note the excellent agreement between the exact numerical result for ℘( ) t and the long time asymptotics (7) in Figs 2 and 3 as shown by the blue lines. Moreover, we emphasise another observation. For a homogeneous system with ϕ = 1 the short and long time asymptotics together fully describe the FPT density. Put differently the overlap region between the regimes is extremely narrow. This holds true in general up to some critical heterogeneity ϕ ⁎ , which will be specified in the following section. Beyond this value ϕ ⁎ a new time scale emerges, as seen in Figs 2b-d and 3b-d) which is not captured by the the short and long time asymptotics and does not correspond to an overlap regime, as we now explain.  Emergence of a new time scale. Here we focus on the regime ϕ > 1 when the inner region has the higher diffusivity. This is the scenario which we would naively expect to enhance the FPT kinetics. The opposite case ϕ < 1 is physically less interesting but can be obtained analogously to the steps presented below. Note that in this regime a resolved third time scale does not exist because short excursions into the outer region are always faster. However, as ϕ → 0 the correct short time asymptotics is obtained in the limit ϕ  s Qx 1 k and  sQ x 1 i , and these are in general not of the Lévy-Smirnov type (5a). Considering the two different types of argument in Eq. (3a), sQ x k and ϕ s Qx l it becomes obvious that an additional separation of time scales occurs in the limit −   Qx s Qx k l 2 1 2 , where k, l stand for the different indices used in our model and depend on the starting position relative to the interface position. In other words, there exists a time scale separation between direct trajectories corresponding to the short time asymptotic (5a) and reflected trajectories accounted for by the long time asymptotic (7).
This new time scale corresponds to trajectories, which are much longer than the direct ones yet much shorter than the reflected ones. In such trajectories the particle ventures into the outward direction away from the target with respect to its initial position. However, this excursion is much shorter than the average time needed to reach the interface. The result are terms of mixed order in s having the form   are given in Eqs (20) and (21) in the Appendix. Details of the calculation will be reserved for a separate longer publication. Both forms, Eq. (11a) and Eq. (11b) hold for ϕ ϕ independently of x 0 . The transition as ϕ crosses the critical value ϕ ⁎ is discontinuous, the factor involved changing , and the functional dependence on t changes concurrently as well. Below ϕ ⁎ the FPT distribution ℘( ) t is completely specified in terms of Eq. (5a) and Eq. (7). This is what we may call a threshold heterogeneity. Moreover, while Eq. (11a) does not reduce to a simple form as ϕ → ∞, Eq. (11b) reduces to a Lévy-Smirnov density, in Eq. (5a). Other than that, Eq. (11b) continuously interpolates between the two different Lévy-Smirnov densities Eqs (5a) and (13) as ϕ increases. Intuitively, this latter limit demonstrates the fact, that reaching the interface from > x x i 0 becomes rate limiting for large ϕ. We also note that the infinite series in Eqs (5a) and (13)  To gain more intuition on how exactly the time scale separation arises we plot the exact numerical result for the FPT distribution ℘( ) t at ϕ = 30 for different initial positions, as depicted in Fig. 4. Starting in the inner region the scale separation emerges as we continuously move the initial position towards the target. Conversely, if starting in the outer region a scale separation emerges as we continuously move the initial position away from the reflecting surface. The necessary and sufficient requirement for heterogeneity controlled kinetics is therefore a Scientific RepoRts | 6:20349 | DOI: 10.1038/srep20349 large heterogeneity and a corresponding existence of two length scales inwards and outwards from the initial position towards the closest boundary or interface.
Mean first passage times are not typical. We now quantify the most likely or typical FPT times, i.e., those which occur most frequently. Already from Figs 2 and 3 it is apparent that there is a large discrepancy in the likelihood of typical and mean FPTs, compare the maximum of ℘( ) t with the dashed vertical line denoting the MFPT. In many cases the most likely FPT is the more relevant quantity. For instance, consider a certain species of bacteria, in which genetic regulation can be viewed as an FPT problem 5 . When we compare the fitness of an individual bacteria in a colony, those who are among the first to respond to an external challenge will be of advantage. Similarly, those predators that first discover a prey have a larger chance of survival.
The MFPT was defined in the previous section, while the typical FPT is corresponds to the extremum whereas in the case > x x i 0 we obtain closed form expressions only in the limits ϕ  1 and ϕ  1, these being 15a

Discussion
We analysed a simple prototype model for Brownian motion in heterogeneous environments in terms of a spherically symmetric geometry with two concentric regions of different diffusivity. We obtained rigorous, asymptotically exact results for the FPT distribution and identified a short time Lévy-Smirnov as well as a long time exponential behaviour. Moreover, we demonstrated the existence of a hitherto overlooked intermediate time scale. All three time scales were interpreted in terms of direct and indirect trajectories. The distinction between the latter also effects the major discrepancy between the mean and the most likely FPTs. Cognisance of this difference is important in many systems. We first focus on the implications of our findings for Brownian motion in spherically symmetric homogeneous media. Our results suggest that the prevailing paradigm of first passage kinetics for Brownian motion 1,14,15 becomes even richer. First, as highlighted already in previous works [17][18][19] the MFPT is often a rough measure of the first passage kinetics: the FPT distribution is typically positively skewed and therefore asymmetric. According to the results in refs 17-19 any two arbitrary trajectories are often more likely to be very different than similar. From our results reported here this result is substantiated in the sense that according to our interpretation an arbitrary trajectory will typically be direct and will not interact with the reflecting boundary. For homogeneous Brownian motion this holds Scientific RepoRts | 6:20349 | DOI: 10.1038/srep20349 strictly for all starting positions satisfying ℘( ) 〈 ( )〉 > 〈 ( )〉 t Q t t x 1 . This statement is particularly important for single molecule observations, where a finite number of trajectories will more likely reveal the typical and not the mean behaviour. We have shown that there is a large discrepancy between the typical FPT and the MFPT. Moreover, an upper bound for the FPT-if starting from an arbitrary position-is set by the MFPT from the confining surface. In the latter case, the distinction between direct and indirect trajectories obviously ceases to exist, giving rise to dominantly exponential statistics. Conversely, the MFPT, while indeed dominated by the rare long time behaviour, is remarkably fully specified by the typical direct trajectories. Hence, the most likely but less significant trajectories turn out to determine the unlikely but dominant trajectories. This adds another surprising feature to the first passage behaviour of Brownian motion, starting with the Lévy arcsine laws for one dimensional free Brownian motion 41,42 .
We demonstrated that a sufficient heterogeneity in the diffusion coefficient gives rise to an additional, intermediate time scale, on which trajectories contain a short excursion towards the external surface. This excursion, however, is much shorter than the typical time needed to diffuse across the entire domain. In this heterogeneity-controlled kinetic domain there thus exist three distinct classes of trajectories: the direct ones, the indirect ones and those which initially make a short indirect excursion and then go directly to the target, as shown in Fig. 1. This latter class of trajectories in fact represents the typical trajectories for starting positions in the outer region of our model system.
The results for our idealised two component model studied in the present work have important consequences for an arbitrary spherically symmetric modulation of the diffusivity. Namely, in such a case there might exist several distinct time scales in the heterogeneity controlled kinetic regime. In fact, heterogeneity controlled kinetics will generically be observed in the presence of sharp changes in the local diffusivity ( ) D r . Such sharp modulations are indeed present in cellular signalling processes when a particle starts in the cytoplasm and searches for its target in the nucleus, where the diffusivity is different. In particular, signalling particles synthesised inside the cell as part of a particular signalling cascade will inherently start away from the cellular membrane, and a separation of scales in the FPT is therefore expected to exist. Moreover, as a large discrepancy is expected to exist in the likelihood for observing typical FPTs with respect to the MFPT, the MFPT is a fairly poor measure for the kinetics at low copy numbers of signalling molecules and for a finite number of possible realisations. Conversely, smooth modulations of ( ) D r are not expected to change the qualitative two time scale picture of FPT kinetics but will of course alter the coefficients in the short and long time asymptotic regimes. We can extend the discussion also to systems with off-centre targets under the condition that the searching particle does not start too close to the external boundary. In this case the long and intermediate time scale asymptotics would remain unchanged-as long as there exist a separation of scales, of course-but the long time asymptotics would be altered. Thus, our rigorous results for our idealised model system are relevant for FPT kinetics in generic heterogeneous media.
We also briefly comment on alternative forms of the diffusion equation (2), the so called Ito and Stratonovich forms 35 . Here we focus on physical stochastic dynamics satisfying the fluctuation-dissipation theorem. In this view, the Ito and Stratonovich equations are nothing but the corresponding Fokker-Planck equations with diffusion coefficient ( ) D r in an external potential ( ) = − ( ) U r k T D r ln , respectively. The MFPT problem in such scenarios including a spurious drift was studied recently in 43 . The FPT kinetics in the presence of such effective external fields will therefore be fundamentally different. We finally note that it will be of further interest to include intermittent active motion in the analysis 44