Noisy continuous time random walks
Experimental studies of the diffusion of biomolecules in the environment of biological cells are routinely confronted with multiple sources of stochasticity, whose identification renders the detailed data analysis of single molecule trajectories quite intricate. Here we consider subdiffusive continuous time random walks, that represent a seminal model for the anomalous diffusion of tracer particles in complex environments. This motion is characterized by multiple trapping events with infinite mean sojourn time. In real physical situations, however, instead of the full immobilization predicted by the continuous time random walk model, the motion of the tracer particle shows additional jiggling, for instance, due to thermal agitation of the environment. We here present and analyze in detail an extension of the continuous time random walk model. Superimposing the multiple trapping behavior with additive Gaussian noise of variable strength, we demonstrate that the resulting process exhibits a rich variety of apparent dynamic regimes. In particular, such noisy continuous time random walks may appear ergodic while the naked continuous time random walk exhibits weak ergodicity breaking. Detailed knowledge of this behavior will be useful for the truthful physical analysis of experimentally observed subdiffusion.
In a series of experiments in the Weitz lab at Harvard, Wong et al. wong () observed that the motion of plastic tracer microbeads in a reconstituted mesh of cross-linked actin filaments is characterized by anomalous diffusion of the form
of the ensemble averaged mean squared displacement (MSD). Here, is the anomalous diffusion exponent of physical dimension , and the anomalous diffusion exponent is in the subdiffusive range . is the probability density function to find the test particle at position at time . Wong et al. demonstrated that the motion of the microbeads is represented by a random walk with subsequent immobilization events of the beads in ‘cages’ within the network wong (). The durations of these immobilization periods were shown to follow the distribution
with the scaling exponent . The exponent turns out to be a function of the ratio between the bead size and the typical mesh size. Thus, when the bead size is larger than the mesh size, the bead becomes fully immobilized, . Conversely, when the bead is much smaller than the typical mesh size, it moves like a Brownian particle [, in Eq. (1)], almost undisturbed by the actin mesh. However, when the bead size is comparable to the mesh size, the motion of the bead is impeded by the mesh, giving rise to a saltatory bead motion in between cages. While the measured distribution (2) of immobilization times indeed captures the statistics of the long time movement of the beads, the measured trajectories also exhibit an additional noise, superimposed to the jump motion according to the law (2) wong (). This additional noise is visible as jiggling around a typical position during sojourn periods of the bead, and appears like regular white noise.
Above example is representative for the task of physical analysis of single particle tracking experiments. Indeed, following single particles such as large labeled molecules, viruses, or artificial tracers in the environment of living cells has become a standard technique in many laboratories, since it reveals the individual trajectories without the problem of ensemble averaging. However, the large amount of potentially noisy data also poses a practical challenge. In physics and mathematics ideal stochastic processes have been investigated for many years, including Brownian motion, fractional Brownian motion, continuous time random walks, Lévy flights, etc. While these can be reasonable approximations for the physical reality, they rarely represent the whole story. In many cases, a superposition of at least two types of stochastic motion are found, such as the ‘contamination’ of the pure stop-and-go motion described by the waiting time distribution (2) with additional noise in the above example of the tracer beads in the actin mesh. Similarly, Weigel et al. in their study of protein channel motion in the walls of living human kidney cells observed a superposition of stochastic motion governed by the law (2) and motion patterns corresponding to diffusion on a fractal weigel (). Tabei et al. show that the motion of insulin granules in living MIN6 cells is best explained by a hybrid model in which fractional Brownian motion is subordinated to continuous time random walk subdiffusion tabei (). Finally, Jeon et al. find that the motion of lipid granules in living yeast cells is governed by continuous time random walk motion at short times, turning over to fractional Brownian motion at longer times lene (). Even if the trajectories are ideal, in the lab we always have to deal with additional sources of noise, for instance, the motion within the sample may become disturbed by a slow, random drift of the object stage in the microscope setup, or simply by the diffusive motion on the cover slip of a living cell, inside which the actual motion occurs that we want to follow.
In recent years a tool box for data analysis was developed by several groups, in particular, in the context of diffusion of tracers in the microscopic cellular environment. These also include fundamental aspects of statistical physics such as the ergodic properties of the underlying process and the irreproducible nature of the time averages of the data pt (). However, to the best of our knowledge these analysis methods and the fundamental properties of statistical mechanics (that is, ergodicity) have not been tested in the presence of additional noise. For instance, if we encounter a non-ergodic process superimposed with an ergodic one, as defined precisely below, what will we find as the result of such an experiment? Here we develop and study in detail a new framework for the motion of a test particle, the noisy continuous time random walk (nCTRW) subject simultaneously to a power-law distribution (2) of sojourn times and additional Gaussian noise. Interestingly, the noise turns out to have a dramatic effect on the time averaged MSD, typically evaluated in single particle tracking experiments, but a trivial effect on the ensemble averaged MSD. We hope that this contribution will be a step towards more realistic modeling of trajectories of single molecules, but also of other diffusive motions of particles in complex environments.
Our analysis is based on two different scenarios for the additional Gaussian noise. Thus, we will consider (i) a regular diffusive motion on top of the power-law sojourn times (2). This scenario mirrors effects such as the diffusion on the cover slip of the biological cell, in which the actual motion occurs that we want to monitor. (ii) We study an Ornstein-Uhlenbeck noise with a typical amplitude which may reflect an intrinsically noisy environment, such as in the case of the tracer beads in the actin network. We will study the MSD of the resulting motion both in the sense of the conventional ensemble average and, for its relevance in the analysis of single particle tracking measurement, the time average. Moreover we demonstrate how a varying strength of the additional noise may blur the result of stochastic diagnosis methods introduced recently, such as the scatter of amplitudes of time averaged MSDs around the average over an ensemble of trajectories, or the -variation method.
The paper is organized as follows. In Sec. II.1, we briefly discuss the concepts of anomalous diffusion and single particle tracking along with the role of (non-)ergodicity in the context of ensemble versus time averaged MSDs. In Sec. II.2 we review some methods of single particle tracking analysis. Section III then introduces our nCTRW model consisting of a random walk with long sojourn times superimposed with Gaussian noise, and we explain the simulation scheme. The main results and discussions are presented in Secs. IV and V: in Sec. IV, we study the statistical behavior of nCTRW motion in the presence of superimposed free Brownian motion, while Sec. V is devoted to the motion disturbed by Ornstein-Uhlenbeck noise. Finally, the conclusions and an outlook are presented in Sec. VI. In the Appendix, we provide a brief outline of the stochastic analysis tools used to quantify the nCTRW processes studied in Secs. IV and V. Note that in the following, for simplicity we concentrate on the one-dimensional case; all results can easily be generalized to higher dimensions.
Ii Anomalous diffusion and single particle tracking
In this Section we briefly review the continuous time random walk (CTRW) model and the behavior of CTRW time series. Moreover, we present some common methods to analyze traces obtained from single particle tracking experiments or simulations.
ii.1 Continuous time random walks and time averaged mean squared displacement
Anomalous diffusion of the subdiffusive form (1) with is quite commonly observed in a large variety of systems and on many different scales havlin (); bouchaud (); report (); report1 (). Examples range from microscopic systems such as the motion of small tracer particles in living biological cells and similar crowded systems golding (); hoefling (); saxton1 (), over the motion of charge carriers in amorphous semiconductors scher () to the dispersion of chemical tracers in groundwater systems grl ().
A prominent model to describe the subdiffusion law (1) is given by the Montroll-Weiss-Scher continuous time random walk (CTRW) montroll (); scher (). Originally applied to describe extensive data from the stochastic motion of charge carriers in amorphous semiconductors scher (), the CTRW model finds applications in many areas. Inter alia, CTRW subdiffusion was shown to underlie the motion of microbeads in reconstituted actin networks wong (), the subdiffusion of lipid granules in living yeast cells lene (); natascha () of protein channels in human kidney cells weigel (), and insulin granules in MIN6 cells tabei (), as well as the temporal spreading of tracer chemicals in groundwater aquifers grl (); brian (). While the sojourn times in a subdiffusive CTRW follow the law (2), the length of individual jumps is governed by a distribution of jump lengths, with finite variance . In the simplest case of a jump process on a lattice, each jump is of the length of the lattice constant. Physically, the power-law form (2) of sojourn times may arise due to multiple trapping events in a quenched energy landscape with exponentially distributed depths of traps monthus (); burov (). Subdiffusive CTRWs macroscopically exhibit long-tailed memory effects, as characterized by the dynamic equation for the probability density function containing fractional differential operators, see below mebakla ().
Usually, we think in terms of the ensemble average (1) when we talk about the MSD of a diffusive process. However, starting with Ivar Nordlund’s seminal study of the Brownian motion of a slowly sedimenting mercury droplet nordlund () and now routinely performed even on the level of single molecules wong (); golding (); lene (); natascha (); weigel (); tabei (); weber (); sm (); bronstein (); szymansky (); natascha1 (); caspi (); banks (), the time series obtained from measuring the trajectory of an individual particle is evaluated in terms of the time averaged MSD
where is the so-called lag time, and is the overall measurement time. For Brownian motion, the time averaged MSD (3) is equivalent to the ensemble averaged MSD (1), , if only the measurement time is sufficiently long to ensure self-averaging of the process along the trajectory. This fact is a direct consequence of the ergodic nature of Brownian motion pt (); stas2 (). To obtain reliable behaviors for the time averaged MSD for trajectories with finite one often introduces the average
over several individual trajectories .
What happens when the process is anomalous, of the form (1)? There exist anomalous diffusion processes, that are ergodic in the above sense that for sufficiently long we observe the equivalence . Prominent examples are given by diffusion on random fractal geometries yasmin (), as well as fractional Brownian motion and the related fractional Langevin equation motion for which ergodicity is approached algebraically deng (); igorg (); jae (). Such behavior was in fact corroborated in particle tracking experiments in crowded environments szymansky (); bronstein (); natascha1 () and simulations of lipid bilayer systems kneller (); matti (). In superdiffusive Lévy walks time and ensemble averaged MSDs differ merely by a factor lws (); zukla (); froemberg (). However, the diverging sojourn time in the subdiffusive CTRW processes naturally causes weak ergodicity breaking in the sense that even for extremely long measurement times ensemble and time averages no longer coincide bouchaud1 (); barkai (); zaid (). For the time averaged MSD the consequences are far-reaching. Thus, for subdiffusive CTRW we find the somewhat counterintuitive result that the time averaged MSD for free motion scales linearly with the lag time, for , and thus He (); lubelski (). In particular, the result decreases with the measurement time , an effect of ageing. Concurrently, the MSD of subdiffusive CTRWs depends on the time difference between the preparation of the system and the start of the measurement johannes (). Under confinement, while the ensemble averaged MSD saturates towards the thermal value, the time averaged MSD continues to grow as as long as stas1 (); neusius (). Another important consequence of the weakly non-ergodic behavior is that the time averaged MSD (3) remains irreproducible even in the limit of : the amplitude of individual time averaged MSD curves scatters significantly between different trajectories, albeit with a well-defined distribution He (); igor (). In other words, this scatter of amplitudes corresponds to a distribution of diffusion constants saxton2 (). Active transport of molecules in the cell exhibits superdiffusion with non-reproducible results for the time averages caspi2 (), a case treated theoretically only recently lws (); froemberg (); akimoto ().
To pin down a given stochastic mechanism underlying some measured trajectories, several diagnosis methods have been developed. Thus, one may analyze the first passage behavior olivier (), moment ratios and the statistics of mean maximal excursions vincent (), the velocity autocorrelation stas2 (); matti (), the statistics of the apparent diffusivities radons (), or the -variation of the data marcin (); kepten (). In the following we analyze the sensitivity of the time averaged MSD and its amplitude scatter as well as the -variation method to noise, that is superimposed to naked subdiffusive CTRWs. In particular, we show that at larger amplitudes of the additional noise the CTRW-inherent weak non-ergodicity may become completely masked.
ii.2 A primer on single particle trajectory analysis
Before proceeding to define the nCTRW process, we present a brief review of several techniques developed recently for the analysis of single particle trajectories.
ii.2.1 Mean squared displacement
Already from the time averaged MSD of individual particle traces important information on the nature of the process may be extracted, provided that the measurement time is sufficiently long. Thus, for ergodic processes the time and ensemble averaged MSD are equivalent, . A significant difference between both quantities points at non-ergodic behavior. In particular, the scatter of amplitudes (distribution of diffusion constants) between individual traces turns out to be a quite reliable measure for the (non-)ergodicity of a process. Thus, in terms of the dimensionless parameter the distribution of amplitudes has a Gaussian profile centered on the ergodic value for finite-time ergodic processes jeon (). Its width narrows with increasing , eventually approaching the sharp distribution at , i.e., each trajectory gives exactly the same result, equivalent to the ensemble averaged MSD. For subdiffusive CTRW processes is quite broad and has a finite value at for completely stalled trajectories stas2 (); He (); jeon (); johannes (); heinsalu (). For the maximum of the distribution is at , while for the maximum is located at . The non-Gaussian distribution for subdiffusive CTRWs is almost independent of , and its shape is already well established even for relatively short trajectories jeon ().
ii.2.2 -variation test
The -variation test was recently promoted as a tool to distinguish CTRW and FBM-type subdiffusion marcin (). It is defined in terms of the sum of increments of the trajectory on the interval as
where . The quantity has distinct properties for certain stochastic processes. Thus, for both free and confined Brownian motion and for any . For fractional Brownian motion, , while . Finally, for CTRW subdiffusion, features a step-like, monotonic increase as function of time , and .
Fig. 1 shows the -variation results for Brownian noise , Ornstein-Uhlenbeck noise , and the naked subdiffusive CTRW . We plot the sum (5) for various, finite values of . For the Brownian noise , monotonically decreases with growing , a signature of the predicted convergence to . The sum appears independent of and proportional to time , as expected. For the Ornstein-Uhlenbeck noise we observe that scales linearly with and the slope decreases with growing , indicating a convergence to zero. Also is linear in . For smaller , the slope increases with and saturates at large . Finally, for CTRW subdiffusion appears to converge to zero for increasing , as predicted. In contrast, has the distinct, monotonic step-like increase expected for CTRW subdiffusion. A more detailed description of the -variation for Ornstein-Uhlenbeck noise and CTRW subdiffusion is found in Appendix A.
Iii Noisy continuous time random walk
In this Section we define the nCTRW model and describe our simulations scheme used in the following Sections.
iii.1 Two nCTRW models
We consider an nCTRW process in which ordinary CTRW subdiffusion with anomalous diffusion exponent is superimposed with the Gaussian noise ,
That means that the Gaussian process is additive and thus independent of . The relative strength of the additional Gaussian noise is controlled by the amplitude parameter . In this study we consider the following two Gaussian processes: (1) in the first case is a simple Brownian diffusive process with zero mean and variance . As mentioned, physically this could represent the (slow) diffusion of a living bacteria or endothelial cell on the cover slip while we want to record the motion of a tracer inside the cell, or the random drifting of the experimental stage. (2) In the second case represents Ornstein-Uhlenbeck process , the confined Brownian motion in an harmonic potential (see Eq. (18) for the definition of the process). Stochastic processes of this second kind, for example, likely mimic the motion of a bead confined in a polymer network or immersed in a macromolecularly crowded environment, where the thermal agitation of the confining environment gives rise to the random fluctuation of the position of the trapped tracer particle around some average value.
According to the definition (6), nCTRW possesses have the following generic properties. Its probability density function (PDF) is given by the convolution of the individual PDFs of a CTRW subdiffusion process and of the Gaussian process,
This chain rule states that a given position of the combined process is given by the product of the probability that the CTRW process has reached the position and the Gaussian process contributes the distance , or vice versa. Here the PDF satisfies the fractional Fokker-Planck equation report ()
where the Riemann-Liouville fractional derivative of is
Physically, this fractional operator thus represents a memory integral with a slowly decaying kernel. From definition (6) it follows that the ensemble averaged MSD is given by
The characteristic function of according to Eq. (7) is given by the product
of the characteristic functions of the individual processes, and . We here use the simplified notation that the Fourier transform of a function is expressed by its explicit dependence on the Fourier variable . With the Mittag-Leffler function we find that report ()
assuming that the CTRW process starts at with initial conditions . The characteristic function initially decays like a stretched exponential and has the asymptotic power-law decay .
The Brownian noise with initial condition has the characteristic function
In the nCTRW process this Brownian noise always dominates the dynamics of the process at long times, since the exponential relaxation (13) dominates the characteristic function . The Ornstein-Uhlenbeck noise with , defined in Eq. (18), has the characteristic function
At short times when confinement by the harmonic Ornstein-Uhlenbeck potential is negligible, Eq. (14) reduces to the result (13) for Brownian noise, and thus the characteristic functions of the two nCTRW processes are identical. At longer times the characteristic function (14) saturates to . This means that the long-time behavior of the nCTRW with superimposed Ornstein-Uhlenbeck noise, largely reflects the naked CTRW process , if the noise level is not too high. Keeping these general features in mind, we further study the statistical quantities of the two nCTRW processes numerically.
iii.2 Simulation of the nCTRW process
To simulate the nCTRW process we independently obtain time traces of the subdiffusive CTRW motion and the additional Gaussian noise. The motion with is generated on a lattice of spacing from the normalized waiting time distribution
with the power-law asymptotic scaling . Here is a scaling constant of dimension . The jump lengths are determined by the -distribution lene (); jae (). Within this construction the CTRW process is associated with the fractional Fokker-Planck equation (8) with the anomalous diffusion exponent eli ()
In the simulations we choose , and consequently in the following times are given in units of , which is also chosen equal to the time increments , at which the system is updated. The lattice spacing is . We simulate the cases and .
The added Gaussian noise is obtained as follows. Brownian noise is obtained at discrete times (with ), in terms of the Brownian walk
where represents white Gaussian noise of zero mean and unit variance with our choice . In the simulations we take , such that . The Ornstein-Uhlenbeck noise is obtained by integration of the Langevin equation
with the initial condition . Here, the coefficient of the restoring force is chosen as . For both Brownian and Ornstein-Uhlenbeck noise, the following values for the noise strength are used: , , , and .
Iv nCTRW with Brownian noise
In Fig. 2 we show typical examples of simulated trajectories of the nCTRW process with added Brownian noise. The trajectories for different noise strengths is constructed for the same CTRW process , i.e., only the noise strength varies in between the panels. This way it is easier to appreciate the influence of the added noise. We generally observe that the trajectories preserve the profile of the underlying CTRW process at small , and they become quite distorted from the original CTRW trajectory for larger values of . In particular, for the largest noise strength the character of the pure CTRW with its pronounced stalling events is completely lost, and visually one might judge these traces to be pure Brownian motion, at least when the length of the time series is not too large. Moreover, the effect of the noise is stronger for smaller . This is because long stalling events occur more frequently as decreases, and thus the actual displacement of the process is also smaller and the influence of the Brownian motion becomes relatively more pronounced.
iv.1 Ensemble-averaged mean squared displacements
From the simulated trajectories we evaluate the ensemble-averaged MSD for the nCTRW and study how its scaling behavior is affected by the Brownian noise . Fig. 3 summarizes the results for the nCTRW process with and . In both cases, a common feature is that the ensemble averaged MSD exhibits a continuous transition from subdiffusion with anomalous diffusion exponent to normal diffusion. This occurs either as the noise strength is increased at fixed time , or as time is increased at a fixed . Due to the additivity of the two contributions we obtain
which is valid as long as is considerably larger than the time increment . Eq. (19) demonstrates that for the subdiffusive CTRW processes with the ensemble averaged MSD of the nCTRW has a crossover in its scaling from at short times to at long times, with the crossover time scale . That is, the effect of the Brownian noise emerges only at long times . Conversely, below the process appears to behave as the bare CTRW process. Note that the crossover time rapidly decreases with increasing as . This explains why we only observe normal diffusion behavior without crossover for the largest noise strength . The crossover time also rapidly increases as the exponent approaches to one [as in our choice of , see Eq. (16)]. Accordingly, in Fig. 3 the ensemble averaged MSDs for do not fully reach the linear regime within the time window of our simulation, in contrast to the case for with .
iv.2 Time-averaged mean squared displacement
We now consider the individual time averaged MSDs of the nCTRW process from single trajectories according to our definition in Eq. (3). Fig. 3 presents the trajectory-to-trajectory average (4). In all cases we find that the time averaged MSDs grow linearly with lag time in the entire range of , showing a clear disparity from the scaling behavior of the ensemble averaged MSDs above. In the time averaged MSD the Brownian noise simply affects the apparent diffusion constant, that is, the effective amplitude of the linear curves. To understand this phenomenon quantitatively, we obtain the analytical form of the trajectory-averaged time averaged MSD,
valid for . Eq. (20) shows that both contributions the naked CTRW and the Brownian process are linearly proportional to the lag time . Note that the exponent and the noise strength only enter into the apparent diffusion constant
where . This result indicates that the effect of the Brownian noise cannot be noticed when we exclusively consider the scaling behavior of the time averaged MSD. Also in terms of the apparent diffusion constant, one hardly notices the presence of the Brownian noise component as long as is small. This agrees with our observations for the time averaged MSD curves in Fig. 3 for noise strengths and .
We also check the fluctuations between individual time averaged MSD curves . Each panel in Fig. 4 plots ten individual time averaged MSDs for the nCTRW process. In all cases, the individual time averaged MSDs display linear scaling with lag time , namely, the scaling behavior of . The individual amplitudes scatter, that is, the apparent diffusion constant fluctuates. With increasing strength of the Brownian component the relative scatter between individual trajectories dramatically diminishes, leading to an apparently ergodic behavior. Thus, for the largest noise strength , the ten trajectories almost fully collapse onto a single curve for .
iv.3 Scatter distribution
We quantify the amplitude scatter of the individual time averaged MSDs in terms of the normalized scatter distributions , where the dimensionless variable stands for the ratio of individual traces versus the trajectory average (compare the derivations in Refs. pt (); He (); stas2 ()). Fig. 5 shows for several values of the lag time . When the noise is negligible (), the observed broad distribution is nearly that of the pure CTRW process. For the distribution has the expected Gaussian profile centered at He (), indicating that long stalling events of the order of the entire measurement time occur with appreciable probability. In the opposite case , the Brownian noise results in a relatively sharply peaked, bell-shaped distribution typical for ergodic processes, at all lag times. An interesting effect of the Brownian noise is that at intermediate strengths it only tends to suppress the contribution at around zero while it does not significantly change the overall profile of the distribution compared to the noise-free case. This implicates that the trajectories share non-ergodic and ergodic elements. For instance, the trajectory of the nCTRW process itself for in Fig. 2 shows a substantially blurred profile of the underlying CTRW process due to the relatively strong Brownian noise; however, the time averaged MSD and its distribution do exhibit non-ergodic behavior, as shown in Fig. 5. The case shown in Fig. 5 features similar albeit less pronounced effects. In particular, for the naked CTRW has its maximum at the ergodic value . Again, the effect of the Brownian noise is to suppress the contribution at .
iv.4 -variation test
We now turn to the -variation test and investigate its sensitivity to the additional noise in the nCTRW process, in comparison to the established results for the naked CTRW . From the trajectory , the partial sum is calculated for finite according to definition (5), where we choose and , compare Section II.2. In Fig. 6, we plot the results of the -variation at increasing for the nCTRW process. We observe that for both cases and , the -sums behave analogously to the predictions for the naked CTRW process, as long as the noise strength remains sufficiently small, according to Fig. 6 this holds for and 0.005 (not shown). In this case, monotonically decreases with increasing , indicating the limiting behavior for large . Meanwhile, approaches the monotonic step-like behavior typical for the CTRW process , as increases. Note that the -sums have plateaus in their increments in analogy to the noise-free case due to the long stalling events in the trajectory.
As the magnitude of the Brownian noise grows larger ( and 0.1), however, the behavior of the -sums changes significantly. While decreases with increasing as for the weaker noise case, it increases linearly with nearly without any sign of plateaus for both nCTRW processes of and 0.8 when the noise strength is increased to . This new feature is the expected behavior of with for a Brownian diffusive process (see Sec. II.2 and Fig. 1). Indeed, it can be shown that for the -sum of the Brownian noise behaves as . We note that for the Brownian noise the -sum with always decays out to zero as . Therefore, the nCTRW process will always have the same variation result of (in the limit of ) as the pure CTRW , even in case that its profile is dominated by large noise.
In contrast, the -sum exhibits a more distinguished effect of the Brownian noise due to the fact that the Brownian noise has . Especially, we find that the noise effect is pronounced for the nCTRW with , when the underlying CTRW process features only few jumps. In this case, the Brownian noise of moderate strength () causes an incline with almost identical slope to the step-like profiles of . At the strongest noise , the step-like behavior typical for the naked CTRW process is nearly masked and the overall tendency follows that of Brownian motion shown in Fig. 1. Accordingly, the -variation test does not properly pin down the underlying CTRW process and thus potentially produces inconsistent conclusion for the nCTRW process. A qualitatively identical behavior is obtained for the nCTRW process when the underlying CTRW process performs relatively frequent jumps (corresponding to the case of ). However, here the effect of the Brownian noise appears weak, because the contribution of relative to the magnitude of becomes smaller at larger values, see the trajectories in Fig. 2.
V nCTRW with Ornstein-Uhlenbeck noise
We now turn to nCTRW processes , in which the superimposed noise is of Ornstein-Uhlenbeck form (18). In this case the influence of the added noise is expected to diminish as the process develops, according to our discussion in Sec. III.
Fig. 7 shows simulated trajectories for the nCTRW process with two different anomalous diffusion exponents, (a) and (b) , for different noise strengths . Indeed we observe that, compared to case of Brownian noise depicted in Fig. 2 the profiles of the jumps and rests of the naked CTRW motion are relatively well preserved despite the mixing with the Ornstein-Uhlenbeck noise . Note that the simulated trajectories with moderate noise strength appear quite similar to the experimental traces of the microbeads in the reconstituted actin network reported by Wong et al. wong (). The noise interference appears considerably lesser for the case , when the magnitude of the net displacements is large relative to the noise contribution due to frequent jumps in the trajectory when is closer to unity. Only for the highest noise strength the pure CTRW behavior with its distinct sojourns becomes blurred by the Ornstein-Uhlenbeck noise.
v.1 Ensemble-averaged mean squared displacement
In Fig. 8 we plot the ensemble averaged MSD of the nCTRW process with Ornstein-Uhlenbeck noise for different noise strengths . We note that, regardless of the intensity , the ensemble averaged MSDs follow the scaling law of the noise-free CTRW process , in particular, at long times. Moreover, all MSD curves at different almost collapse onto each other, although small differences are discernible at short times. These results suggest that, in contrast to the Brownian noise case discussed in Sec. IV, the Ornstein-Uhlenbeck noise does not critically interfere with the diffusive behavior of the noise-free CTRW motion, as expected. To obtain a quantitative understanding of these results we derive the analytic form of the ensemble averaged MSD,
Here the last term stems from the contribution of the Ornstein-Uhlenbeck noise. For it saturates to the constant , which is typically small relative to . Hence, at long times , the noise term is negligible, and the ensemble averaged MSD grows as , consistent with the observations in Fig. 8. In the opposite case for , the contribution of the noise is expanded to obtain . As discussed for the characteristic function (14), on these time scales the Ornstein-Uhlenbeck noise has the same form as the Brownian noise , leading to the same scaling form (19) for the ensemble averaged MSD. However, as studied in the previous case, the effect of the linear scaling is irrelevant at short times.
v.2 Time-averaged mean squared displacement
We now turn to the time averaged MSD curves from individual trajectories of the nCTRW process. Fig. 8 shows the trajectory-averaged time averaged MSD for different noise strengths. For both anomalous diffusion exponents and , we observe qualitatively the same behavior. On the one hand, the CTRW motion superimposed with moderate noise ( and ) leads to a linear scaling of the time averaged MSD with lag time , with almost identical amplitude. On the other hand, when the noise amplitude becomes large, the scaling of the time averaged MSD is significantly affected. We find that these results are consistent with the analytical form of the trajectory-averaged time averaged MSD,
valid at lag times . In this expression it is worthwhile to point out that the contribution of the naked CTRW process is decreased as the length of the trajectory becomes longer, due to the aging effect of the decreasing effective diffusion constant , while the noise is independent of . Due to this effect the time averaged MSD (23) has three distinct scaling regimes: (i) at lag times , the time averaged MSD is linearly proportional to the lag time with apparent diffusion constant . (ii) At lag times , the time averaged MSD is again proportional to . On this timescale, however, the noise part cannot be ignored, and the apparent diffusion constant is given by . Note that the diffusion constant at short lag times is larger than the one at long times. (iii) For lag times , the time averaged MSD is that of confined Brownian diffusion where . These three regimes are expected to occur only in the presence of large noise strengths when (see the case in Fig. 8). In the opposite case, when is negligible compared to , the three regimes are indistinguishable, and only one scaling law is observed at all lag times. This behavior is seen for the two cases of and in Fig. 8.
Interestingly, the time averaged MSD for the nCTRW with anomalous diffusion exponent and noise strength is reminiscent of the MSD curves observed for micron-sized tracer particles immersed in wormlike micellar solutions, which are known to behave as a viscoelastic polymer network when the micelles concentration is above a critical value dreiss (). Experiments using diffusing wave spectroscopy galvan (); bellour () and single-particle tracking natascha1 () revealed that the immersed particles exhibit three distinct diffusive behaviors in different time windows. It was shown that particles surrounded by the micellar network undergo a Brownian diffusion at short (sub-milliseconds) times until they engage with the caging effects of the micellar network, while at later times (milliseconds to sub-seconds) one observes a seemingly confined Brownian motion. It turns out that this confined diffusion is in fact a pronounced subdiffusive motion characterized by anti-persistent spatial correlation induced by the polymer network, governed by the fractional Langevin equation natascha1 (). At macroscopic times, when the wormlike micelle solution behaves as a viscous fluid, the particle again shows a Brownian diffusion, albeit with a significantly reduced diffusion constant. The results obtained here suggest that due to the presence of the Ornstein-Uhlenbeck noise the CTRW process could be mistakenly interpreted to conform to a physically different system, i.e., Brownian motion, and, thus, one needs to be careful in analyzing the data with several possible models, and to use several complementary diagnosis tools.
The expression (23) for the time averaged MSD suggests that stops aging and reflects almost entirely the character of the Ornstein-Uhlenbeck noise process if . This is indeed shown in Fig. 9 where the time averaged MSD is plotted as a function of the overall measurement time at a fixed lag time for nCTRW with , together with the theoretical prediction Eq. (23). For the weakest noise whose crossover time is beyond , the time averaged MSD only displays aging of the naked CTRW process with scaling . When the noise strength is increased to , the time averaged MSD starts to show ergodic behavior as the measurement time gets larger than the crossover time . In the extreme case when the nCTRW process is dominated by the noise (here ), effectively no aging is observed in and the process appears ergodic.
In Fig. 10 we plot ten individual time averaged MSD curves . For the case of more pronounced subdiffusion (), the trajectory-to-trajectory variations are significant, due to the combined effect of long-time stalling events and the Ornstein-Uhlenbeck noise. Thus, for the smallest noise strength some time averaged MSDs exhibit a large deviation from the linear scaling expected from the trajectory-averaged time averaged MSD. In this case, the long stalling events, that are of the order of the measurement time , lead to the plateaus in . Intriguingly, such plateaus also appear in the presence of the largest noise strength, . Here, however, they represent the confined diffusion of the Ornstein-Uhlenbeck noise in which one or few long-stalling events occur in the the CTRW process. For the CTRW process with more frequent jumps () the individual time averaged MSDs follow the expected scaling behavior with smaller amplitude fluctuations in , as now the noise is relatively stronger.
v.3 Scatter distribution
Fig. 11 illustrates the normalized scatter distribution of the individual time averaged MSDs obtained from trajectories, as a function of the dimensionless variable . The overall behavior is qualitatively consistent with those of the scatter distributions for added Brownian noise, as shown in Fig. 5. The distribution becomes increasingly ergodic as the Gaussian noise increases, especially for the case of the more pronounced subdiffusive process (). Here, the finite contribution at in is gradually suppressed and the peak of the distributions is approaching as the strength of the noise is increased.
From the time averaged MSD (23) and its aging property in Fig. 9 it is also expected that the distribution attains apparent ergodic features as the observation time is increased. We show this in Fig. 12 where the scatter distributions for the same nCTRW process simulated up to are compared to the previous ones with . The general trend is that becomes sharper as increases. Importantly, this ergodic effect becomes relevant provided that is increased to at least be comparable with the crossover time defined above. As seen for the case for , the distribution is almost unaffected, except at where . In contrast to this, the distribution at is noticeably narrower around when is increased ().
v.4 -variation test
In Fig. 13 we show the variation of the -sums and at increasing for the simulated trajectories of Fig. 7. For the smallest noise strength shown in Fig. 13 the -variation test produces the result expected for the naked CTRW process. In the case of nCTRW with this behavior persists in the presence of large noise strengths up to . This is due to the fact that the naked CTRW process features large displacements relative to those of the Ornstein-Uhlenbeck process. In contrast, for the profiles of are mildly affected by the Ornstein-Uhlenbeck noise. Namely, the plateaus are tilted, with a somewhat larger slope at higher , while their overall profiles preserve features of the monotonic, step-like behavior of the CTRW process .
Concurrently, we note that at larger values of neither nor show any indication of CTRW for both and . At these noise strengths, the Ornstein-Uhlenbeck noise dominates the -variation results. As explained in Section II.2 and shown in Fig. 1, in the limit the -sums for the Ornstein-Uhlenbeck noise converge to the results for Brownian noise. Thus, while rapidly decreases with increasing in the case of Brownian noise, such a tendency is also present for the Ornstein-Uhlenbeck noise. Moreover, the spike-like profiles in at for reflects the property of the Ornstein-Uhlenbeck noise (see Fig. 1). For , we find that the -sum of the Ornstein-Uhlenbeck noise behaves as (see App. A). This result explains why monotonically increases with up to , and for larger it grows as , as shown in Fig. 13. As in the case of the above nCTRW in the presence of Brownian noise , we find that the -variation result may be substantially affected by the added Ornstein-Uhlenbeck noise, and the identification of the underlying CTRW process become impossible.
We introduced the noisy CTRW process, in which an ordinary CTRW process is superimposed with Gaussian noise, representing physically relevant cases when the pure CTRW motion becomes distorted by a noisy environment. We investigated how the additional ergodic noise interferes with the non-ergodic behaviors of the underlying subdiffusive CTRW motion. Considering the two types of Gaussian noise, Brownian noise and Ornstein-Uhlenbeck noise, we simulated the resulting nCTRW motion and studied physical quantities such as the ensemble and time averaged MSDs, the amplitude scatter distribution, and the behavior of the -variation.
The analysis demonstrates that the influence of the Gaussian noise on these statistical quantities is highly specific to the quantity of concern. Moreover, it depends not only on the type of the noise and its strength but also on the length of the trajectory and the time scale. Depending on those specific conditions a quantitative analysis of the nCTRW process may reveal or mask the underlying non-ergodicity of the naked CTRW process. Thus care is needed when we want to diagnose the stochastic nature of a physical process based on experimental data. One way to avoid wrong conclusions is to apply complementary analysis techniques, such as the quantities used herein, or moment ratios, mean maximal excursion methods, first passage dynamics, or others. The other necessary ingredient is a good physical intuition for the observed process. It would be also interesting to find analytical expressions for the -variations, of mixed processes, as those we have considered here. Our simulations results show that on finite measurement time the noise is crucial, and could easily destroy our basic understanding of the underlying process.
From the present study an experimentally relevant inverse problem can be posed. Can one filter out the Gaussian noise from the experimentally obtained nCTRW process? Although obtaining the noise-cleansed profile from a given nCTRW trajectory may appear infeasible, one could in principle obtain noise-free contribution in some ensemble- or time-averaged physical quantities of the nCTRW process provided one is able to attain a sufficiently long trajectory. In the ensemble averaged MSD the noise survives if the noise is Brownian motion while the CTRW process wins if the noise is of Ornstein-Uhlenbeck type [see Eq. (22)]. This is of course what we expect, since the MSD of the Brownian motion increases linearly with time, for CTRW like , and is a constant for Ornstein-Uhlenbeck—hence it is not surprising to see this behavior. In contrast for the time averaged MSD the dominant contribution comes always from the noise [Eqs. (20) and (23)] in the sense that when the measurement time is very large even the bounded Ornstein-Uhlenbeck noise wins over. This is due to an aging effect, the time averaged MSD of the CTRW process decreases with measurement time . We thus see that the influence of the noise on time averages is fundamentally different from that on ensemble averages. As an example, one can extract almost solely the noise contribution in the time-averaged MSD (and not the bare CTRW itself) from a very long nCTRW trajectory. For an ergodic process this corresponds to the almost identical noise contribution in the ensemble-averaged MSD, namely, . How to subtract this noise from real data is left for future work.
The nCTRW process developed herein is a physical extension of pure CTRW dynamics. We believe that it represents an important advance in the truthful description of anomalous diffusion data in thermal microscopic systems, where the environment is noisy by definition. Mathematically, the nCTRW process is quite intuitive, due to the additivity of the Gaussian noise.
Our present study can naturally be extended to more complicated noise sources, such as fractional Brownian motion (FBM), in order to obtain insight into intracellular anomalous diffusion that show both CTRW and FBM behaviors. It is expected that although FBM-like noise should lead to similar effects on the statistical behavior of the nCTRW process, the quantitative results will be profoundly different due to the scaling law of the FBM-like noise with exponent .
Acknowledgements.We acknowledge funding from the Academy of Finland within the Finland Distinguished Professor (FiDiPro) scheme and from the Israel Science Foundation.
Appendix A -variation of CTRW subdiffusion and Ornstein-Uhlenbeck noise
We here discuss the -variation properties of CTRW subdiffusion and the Ornstein-Uhlenbeck noise in some more detail.
CTRW subdiffusion. The subdiffusive CTRW process with its PDF governed by the fractional Fokker-Planck equation (8) can be described through the subordinated Brownian motion , where defined with the internal time is an ordinary Brownian motion satisfying and . Here, is the so called inverse subordinator matching the laboratory time to the internal time . For the CTRW process , the -sum (t) as grows to infinity satisfies marcin ()
As shown in Fig. 1, has a step-like incremental profile and its jump times represent those for a given realization of CTRW process . On the other hand, the -sum (with ) decreases with increasing , finally at marcin (). This is also shown in the simulations in Fig. 1.
Ornstein-Uhlenbeck noise. The ensemble average of the -variation sum at finite becomes
neglecting an additional term which becomes negligible for large . For small , the above sum simplifies to
In this case the linear slope of increased with up to values (for the given parameter values used in our simulation). This is shown in Fig. 1. In the other case, when , the sum converges to the result of the Brownian noise . Hence, for large , is proportional to with an -independent slope. From the fact that on short-time scales (as ) the Ornstein-Uhlenbeck process behaves like a free Brownian process, it can be inferred that its -variation results are identical with those for simple Brownian motion. Therefore, is expected to converge to zero as increases. This is indeed observed in the simulations result in Fig. 1. We find that when is small, the increment of exhibits a spike-like profile. This behavior presumably occurs due to the fact that the quartic moment of the displacement becomes very small when the lag time is larger than the relaxation time of the process for small .
- (1) I. Y. Wong, M. L. Gardel, D. R. Reichman, E. R. Weeks, M. T. Valentine, A. R. Bausch, and D. A. Weitz, Phys. Rev. Lett. 92, 178101 (2004).
- (2) A. V. Weigel, B. Simon, M. M. Tamkun, and D. Krapf, Proc. Nat. Acad. Sci. USA 108, 6438 (2011).
- (3) S. M. A. Tabei, S. Burov, H. Y. Kim, A. Kuznetsov, T. Huynh, J. Jureller, L. H. Philipson, A. R. Dinner, and N. F. Scherer, Proc. Natl. Acad. Sci. USA 110, 4911 (2013).
- (4) J.-H. Jeon, V. Tejedor, S. Burov, E. Barkai, C. Selhuber-Unkel, K. Berg-Sørensen, L. Oddershede, and R. Metzler, Phys. Rev. Lett. 106, 048103 (2011).
- (5) E. Barkai, Y. Garini, and R. Metzler, Physics Today 65(8), 29 (2012).
- (6) S. Havlin and D. Ben-Avraham, Adv. Phys. 36, 695 (1987).
- (7) J.-P. Bouchaud and A. Georges, Phys. Rep. 195, 127 (1990).
- (8) R. Metzler and J. Klafter, Phys. Rep. 339, 1 (2000).
- (9) R. Metzler and J. Klafter, J. Phys. A 37, R161 (2004).
- (10) I. Golding and E. C. Cox, Phys. Rev. Lett. 96 098102 (2006).
- (11) F. Höfling and T. Franosch, Rep. Prog. Phys. 76, 046602 (2013).
- (12) M. J. Saxton and K. Jacobson, Annu. Rev. Biophys. Biomol. Struct. 26, 373 (1997).
- (13) H. Scher and E. W. Montroll, Phys. Rev. B 12, 2455 (1975).
- (14) H. Scher, G. Margolin, R. Metzler, J. Klafter, and B. Berkowitz, Geophys. Res. Lett. 29, 1061 (2002).
- (15) E. W. Montroll and G. H. Weiss, J. Math. Phys. 6, 167 (1965).
- (16) N. Leijnse, J.-H. Jeon, S. Loft, R. Metzler, and L. B. Oddershede, Eur. Phys. J. Special Topics 204, 75 (2012).
- (17) B. Berkowitz, A. Cortis, M. Dentz, and H. Scher, Rev. Geophys. 44, RG2003 (2006).
- (18) C. Monthus and J.-P. Bouchaud, J. Phys. A 29, 3847 (1996).
- (19) S. Burov and E. Barkai, Phys. Rev. Lett. 98, 250601 (2007).
- (20) R. Metzler, E. Barkai, and J. Klafter, Phys. Rev. Lett. 82, 3563 (1999); Europhys. Lett. 46, 431 (1999); E. Barkai, R. Metzler, and J. Klafter, Phys. Rev. E 61, 132 (2000).
- (21) I. Nordlund, Z. Phys. Chem. 87, 40 (1914).
- (22) C. Bräuchle, D. C. Lamb, and J. Michaelis, Single Particle Tracking and Single Molecule Energy Transfer (Wiley-VCH, Weinheim, Germany, 2012); X. S. Xie, P. J. Choi, G.-W. Li, N. K. Lee, and G. Lia, Annu. Rev. Biophys. 37, 417 (2008).
- (23) I. Bronstein, Y. Israel, E. Kepten, S. Mai, Y. Shav-Tal, E. Barkai, and Y. Garini, Phys. Rev. Lett. 103, 018102 (2009).
- (24) J. Szymanski and M. Weiss, Phys. Rev. Lett. 103, 038102 (2009).
- (25) J.-H. Jeon, N. Leijnse, L. B. Oddershede, and R. Metzler, New J. Phys. 15, 045011 (2013).
- (26) A. Caspi, R. Granek and M. Elbaum, Phys. Rev. Lett. 85, 5655 (2000).
- (27) D. Banks and C. Fradin, Biophys. J. 89, 2960 (2005).
- (28) S. C. Weber, A. J. Spakowitz, and J. A. Theriot, Phys. Rev. Lett. 104, 238102 (2010).
- (29) S. Burov, J.-H. Jeon, R. Metzler, and E. Barkai, Phys. Chem. Chem. Phys. 13, 1800 (2011).
- (30) Y. Meroz, I. Eliazar, and J. Klafter, J. Phys. A 42, 434012 (2009); Y. Meroz, I. M. Sokolov, and J. Klafter, Phys. Rev. Lett. 110, 090601 (2013).
- (31) W. Deng and E. Barkai, Phys. Rev. E 79, 011112 (2009).
- (32) J.-H. Jeon and R. Metzler, Phys. Rev. E 81, 021103 (2010). Note, however, the occurrence of transient non-ergodic relaxation to equilibrium under confinement, as shown in J.-H. Jeon and R. Metzler, Phys. Rev. E 85, 021147 (2012).
- (33) I. Goychuk, Phys. Rev. E 80, 046125 (2009); Adv. Chem. Phys. 150, 187 (2012).
- (34) G. R. Kneller, K. Baczynski, and M. Pasenkiewicz-Gierula, J. Chem. Phys. 135, 141105 (2011).
- (35) J.-H. Jeon, H. Martinez-Seara Monne, M. Javanainen, and R. Metzler, Phys. Rev. Lett. 109, 188103 (2012).
- (36) G. Zumofen and J. Klafter, Physica D 69, 436 (1993).
- (37) A. Godec and R. Metzler, Phys. Rev. Lett. 110, 020603 (2013).
- (38) D. Froemberg and E. Barkai, Phys. Rev. E. Rapid Com. 87, 031104(R) (2013).
- (39) J.-P. Bouchaud, J. Phys. I (Paris) 2, 1705 (1992).
- (40) G. Bel and E. Barkai, Phys. Rev. Lett. 94, 240602; A. Rebenshtok and E. Barkai, ibid. 99, 210601 (2007); J. Stat. Phys. 133, 565 (2008).
- (41) M. A. Lomholt, I. M. Zaid, and R. Metzler, Phys. Rev. Lett. 98, 200603 (2007); I. M. Zaid, M. A. Lomholt, and R. Metzler, Biophys. J. 97, 710 (2009).
- (42) Y. He, S. Burov, R. Metzler, and E. Barkai, Phys. Rev. Lett. 101, 058101 (2008).
- (43) A. Lubelski, I. M. Sokolov, and J. Klafter, Phys. Rev. Lett. 100, 250602 (2008).
- (44) J. H. P. Schulz, E. Barkai, and R. Metzler, Phys. Rev. Lett. 110, 020602 (2013); E. Barkai, ibid. 90, 104101 (2003); E. Barkai and Y-C. Cheng, J. Chem. Phys. 118, 6167 (2003).
- (45) S. Burov, R. Metzler, and E. Barkai, Proc. Natl. Acad. Sci. USA 107, 13228 (2010).
- (46) T. Neusius, I. M. Sokolov, and J. C. Smith, Phys. Rev. E 80, 011109 (2009).
- (47) I. M. Sokolov, E. Heinsalu, P. Hänggi and I. Goychuk, Europhys. Lett. 86, 30009 (2009).
- (48) M. J. Saxton, Biophys. J. 72, 1744 (1997).
- (49) A. Caspi, R. Granek, and M. Elbaum, Phys. Rev. E 66m 011916 (2002).
- (50) T. Akimoto, Phys. Rev. Lett. 108, 164101 (2012).
- (51) S. Condamin, V. Tejedor, R. Voituriez, O. Bénichou, and J. Klafter, Proc. Natl. Acad. Sci. USA 105, 5675 (2008).
- (52) V. Tejedor, O. Bénichou, R. Voituriez, R. Jungmann, F. Simmel, C. Selhuber-Unkel, L. Oddershede, and R. Metzler, Biophys. J. 98, 1364 (2010).
- (53) M. Bauer, R. Valiullin, G. Radons, and J. Kaerger, J. Chem. Phys. 135, 144118 (2011).
- (54) M. Magdziarz, A. Weron, K. Burnecki, and J. Klafter, Phys. Rev. Lett. 103, 180602 (2009); M. Magdziarz and J. Klafter, Phys. Rev. E 82, 011129 (2010).
- (55) K. Burnecki, E. Kepten, J. Janczura, I. Bronshtein, Y. Garini, and A. Weron, Biophys. J. 103, 1839 (2012).
- (56) I.M. Sokolov, E. Heinsalu, P. Hänggi, and I. Goychuk, EPL 86, 30009 (2009).
- (57) J.-H. Jeon and R. Metzler, J. Phys. A: Math. Theor. 43, 252001 (2010).
- (58) E. Barkai, R. Metzler, and J. Klafter, Phys. Rev. E. 61, 132 (2000). Note the occurrence of the factor here due to the different choice of in Eq. (15).
- (59) C. A. Dreiss, Soft Matter 3, 956 (2007).
- (60) J. Galvan-Miyoshi, J. Delgado, R. Castillo, Eur. Phys. J. E 26, 369 (2008).
- (61) M. Bellour, M. Skouri, J.-P. Munch, and P. Hébraud, Eur. Phys. J. E 8, 431 (2002).