Investigating the interplay between mechanisms of anomalous diffusion via fractional Brownian walks on a comb-like structure

The comb model is a simplified description for anomalous diffusion under geometric constraints. It represents particles spreading out in a two-dimensional space where the motions in the x-direction are allowed only when the y coordinate of the particle is zero. Here, we propose an extension for the comb model via Langevin-like equations driven by fractional Gaussian noises (long-range correlated). By carrying out computer simulations, we show that the correlations in the y-direction affect the diffusive behavior in the x-direction in a non-trivial fashion, resulting in a quite rich diffusive scenario characterized by usual, superdiffusive or subdiffusive scaling of second moment in the x-direction. We further show that the long-range correlations affect the probability distribution of the particle positions in the x-direction, making their tails longer when noise in the y-direction is persistent and shorter for anti-persistent noise. Our model thus combines and allows the study/analysis of the interplay between different mechanisms of anomalous diffusion (geometric constraints and long-range correlations) and may find direct applications for describing diffusion in complex systems such as living cells.


Introduction
Understanding diffusive motions is a long-standing problem in the physicists' agenda. Since the works of Einstein, Smoluchowski and Langevin, we know that one of the most striking patterns of classical (usual) free diffusion is the linear growth with time of the second moment of the particle positions. While this feature is ubiquitous in nature, many other patterns have been also observed for the second moment, both from analytical calculations and experimental data [1]. These deviants behaviors are often called anomalous diffusion.
Nowadays, efforts have been mainly focused on understanding the underlying physical mechanisms that lead to such deviations, since only the classification of the diffusive behaviors of a given system is not enough to fully characterize it. In this context, two common pathways for anomalous diffusion are long-range correlations in the particle positions and geometric constraints related with the structural complexity of the environment where the particles are moving. The former may be related to memory affects and an example is the fractional Brownian motion [2,3], where the second moment has a power-law dependence with time, that is, ⟨x 2 (t)⟩ ∝ t 2h . For this case, the Hurst exponent h classifies the anomalous diffusion in subdiffusion (0 < h < 1 2) or superdiffusion (h > 1 2) [4,5] and also recovers usual (h = 1 2) and ballistic (h = 1) diffusions as limiting cases. Anomalous diffusion due to geometric constraints is well exemplified by the comb model. In this model, diffusive particles walk on a two-dimensional space; however, motions in the x-direction are only allowed when the y coordinate of the particle positions is zero. The result of this simple constraint is a subdiffusive motion (precisely, ⟨x 2 (t)⟩ ∼ t α with α = 1 2) and the appearance of a backbone (at y = 0) and teeth (along the y−axis direction), which were originally proposed to mimic the quasilinear structure and dangling ends of percolation clusters [6,7,8].
Mainly due to its simplicity and ability of mimicking diffusive aspects of highly disordered systems, the comb model has been extensively studied and extended by means of analytical calculations [9,10,11,12,13,14,15] and also employed as simplified description of natural phenomena such as cancer proliferation [16,17], transport of spiny dendrites [18] and diffusion of ultra cold atoms [19,20]. Moreover, as crowded environments tend to slow down the diffusive motion of particles by obstruction and/or trapping [21], the comb model could also be used as a toy model for describing diffusion in intracellular processes, where several researchers have reported the existence of subdiffusion and ergodicity breaking in the random walk of different biochemical compounds [22,23,24,25,27,26,28,21,29,30,31,32] by using single particle tracking [33]. In the case of intracellular processes, the power-law exponents α describing the second moment ranges from α ≈ 0.1 to α ≈ 0.9 (evaluated both from ensemble and time average) and, in the case of ultra cold atoms, α appears as a function of lattice depth and assumes values larger than one (α > 1, superdiffusion).
In addition to the anomalous scaling of the second moment, the work of Weiss [34] has also shown that the diffusive processes of intracellular fluids may present long-range corrections as in the fractional Brownian motion, an empirical finding that have not been explored within the comb model framework. Also, recently, Yamamoto et al. [35] have provided evidence (from molecular dynamics simulation) that the anomalous diffusion of water molecules on cell membrane surface arises from both divergent mean trapping time and long-range correlated noises. Here, we extend the comb model by considering a Langevin-like equation where the noises are long-range correlated. Our model thus combines two of the main mechanisms of anomalous diffusion: geometric constraints and long-range correlations, allowing a more complete description of complex diffusive motions and also the study of the interplay between these two mechanisms. By carrying out computer simulations, we show that the correlations in the y−direction affect the diffusive behavior in the x−direction in a non-trivial fashion and that a quite rich diffusive scenario emerges from these long-range correlations. We report that diffusive process in the x−direction of this generalized version of the comb model can be usual, superdiffusive or subdiffusive, depending on both corrections in the x and y directions. We further show that the long-range correlations affect the probability distribution of the particle positions in the x−direction: the distribution becomes more leptokurtic when the noise in the y−direction is persistent whereas, for anti-persistent noise, it becomes less leptokurtic and approaches the Gaussian distribution as a limiting case.

Generalized comb and long-range correlations
We start by writing the generalized diffusion equation describing the motion under the comb structure [36]: here ρ(x, y; t) is the joint probability of the particle positions (x and y coordinates) as a function of time t, D x and D y are the diffusive coefficients in the x and y directions. Note the presence of the Dirac delta δ(y) multiplying the spatial derivative with respect to x and, consequently, limiting the diffusion in the x−direction only for y = 0. In order to investigate the role of long-range correlations in this diffusive process, we propose to numerically solve the following coupled Langevin-like equations where β x and β y are constants related to D x and D y , η x (t) and η y (t) are fractional Gaussian noises [2,3], both with zero mean (⟨η x (t)⟩ = 0 and ⟨η y (t)⟩ = 0), unitary variances (⟨η 2 x (t)⟩ = 1 and ⟨η 2 y (t)⟩ = 1) and correlations functions (for t 1 ≠ t 2 ) given by Notice that the noises η x (t) and η y (t) are power-law correlated with scaling exponents (2h x −1) and (2h y −1), where h x and h y are the so-called Hurst exponents. For h x,y < 1 2, the correlation functions have negative signs and the noises are anti-persistent, meaning that positive values are followed by negative values (or vice-versa) more frequently than by chance; while for h x,y > 1 2 the correlations are positive, meaning that positive values are followed by positive values and negative values are followed by negative values more frequently than by chance. Fully persistent noises occur as a limiting case when h x,y = 1, whereas normal Brownian (white) noises correspond to the limit of h x,y = 1 2. When running out the simulations of the Langevin equations, we have considered a narrow band of thickness ε along the x−axis, inside which the diffusion in the x−direction occurs. This strip mimics the effect of the Dirac delta δ(y) (that appears multiplying the noise η x (t) in Eq. 2) and we have verified that the value ε has no influence in the diffusive process, as long as ε and noises amplitudes β x and β y are of the same order of magnitude. This restriction ensures that the variable y presents no dynamics inside the strip. For our proposes, we have fixed ε = β x = β y = 0.1. The fractional Gaussian noises η x (t) and η y (t) were generated by following the Hosking method [37]. Figure 1 shows examples of simulated trajectories for three different values of h y and keeping h x = 0.5. By visual inspection, we note that depending on whether the motion in y is anti-persistent [Fig 1(a)] or persistent [Fig 1(c)], the diffusive behavior in the x−direction drastically changes. The anti-persistence in y(t) causes the particle to return more often to the strip where the diffusive motion in the x−direction is allowed; consequently, x(t) covers a larger interval along the x−axis when compared with the usual comb [h x = h y = 0.5, Fig 1(b)]. On the other hand, persistence in y(t) causes the particle to return less often to the strip and, consequently, x(t) covers a small interval along the x−axis.
Before quantifing the effects of the long-range correlations on the variances of x(t) and y(t), we have first verified whether the Langevin equations [Eq. (2)] are actually equivalent to the generalized diffusion equation [Eq. (1)] when the noises η x (t) and η y (t) are uncorrelated (h x = h y = 0.5). We thus simulate an ensemble of 10 5 particle positions considering dt = 1 and h x = h y = 0.5. By using these data, we evaluate the temporal evolution of the variances where ⟨. . . ⟩ denotes ensemble average. Figures 2(a) and (b) show the behavior of σ 2 x (t) and σ 2 y (t) versus time for several values of the maximum integration time t max in the Langevin equations, ranging from 2 7 to 2 14 as indicated in these plots. In both cases, we note a finite-size-like effect characterized by a crossover time t × where σ 2 x (t) and σ 2 y (t) change their behaviors. Despite that and also discarding an initial transient regime (t ≲ 10), we observe that the variance profiles of σ 2 x (t) and σ 2 y (t) are in well agreement with the power laws predicted by the generalized diffusion equation of the usual comb [Eq. (1)], that is, σ 2 x (t) ∼ t αx and σ 2 y (t) ∼ t αy with α x = 0.5 and α y = 1 [36,38,39]. We have also analyzed the behavior of the crossover time t × in function of the maximum integration time t max . Figures 2(c) and (d) show that t × grows linearly with t max for σ 2 x (t) and σ 2 y (t), confirming that the simulated power-law regimes prevail when t max tends to infinity.
In order to further strengthen the equivalence between our Langevin equations (2) and the diffusion equation (1), we have also evaluated the marginal probability functions P (x; t) = ∫ ρ(x, y; t)dy and P (y; t) = ∫ ρ(x, y; t)dx from the simulated data.  Figures 3(c) and 3(d) reveal a good collapse of the empirical distributions and also a good agreement with the analytical distributions (dashed lines) predicted by Eq. (1), which is a standard Gaussian for P (ξ y ) and can be expressed in terms of the Fox H function for P (x; t) [36,38,39] or by the following summation [5,30]  x (t) ∼ t αx with α x = 0.5 and σ 2 y (t) ∼ t αy with α y = 1 [36,38,39]. We note a finite-size-like effect in both plots. In (c) and (d), we show the crossover time t × , where the variances σ 2 x (t) and σ 2 y (t) start to change their behaviors, as a function of the maximum integration time t max . We observe that t × grows linearly with t max in both cases; therefore, in the limit as t max approaches infinity, the power-law regimes must prevail.
Thus, for the case where no correlations are present in the noises η x (t) and η y (t), the results shown in Fig. 3 corroborate with the equivalence between our simulations of the Langevin equations and the diffusion equation for the comb model.
Convinced of the equivalence in the usual case, we now address the role of the long-range correlations on the diffusive properties of a random walker described by our Langevin equations (2). We first investigate how the memory effects in the y−direction modify the diffusion in the x and y directions. In order to do so, we have built an ensemble of 10 5 particle positions considering dt = 1 and h x = 0.5 for each one of nineteen values of h y homogeneously distributed between 0.05 and 0.95. For each value of h y , we evaluate the temporal evolution of the variances σ 2 x (t) and σ 2 y (t). Figures 4(a) and 4(b) show the variance profiles in the x and y directions, respectively. In these plots, the dashed lines represent the behavior for the usual comb, that is, when no correlations are present. In the case of σ 2 x (t), we note that the more the noise in y is anti-persistent (h y < 0.5), the larger is the power-law exponent α x describing the behavior of σ 2 x (t); in contrast, persistent noise in y slows down the diffusion and decreases the power-law . In addition to the good (visual) agreement between simulations and analytical expressions, we have found that the p-values of the Kolmogorov-Smirnov test (testing the equality of the empirical and analytical distributions) are all larger than 0.05, indicating that the hypothesis of the analytical distributions describing our data cannot be rejected at a confidence level of 95%.
exponent α x . The extreme cases where h y = 0.05 and h y = 0.95 are characterized by α x ≈ 0.95 and α x ≈ 0.19, respectively. Thus, as we have first noticed in Fig. 1, antipersistent noise in y forces the walker to return more frequently to the strip where it can move in the x−direction and, consequently, the diffusion is enhanced when compared with the case usual comb (h y = 0.5). On the other hand, persistent noise in y drives the walker away from strip and, consequently, the diffusion slows down. In the case of σ 2 y (t), the behavior is straightforward: anti-persistent noise in y slows down the diffusion while persistent noise enhances it [see Fig. 4(b)]. In order to quantify the role of long-range correlations in y on σ 2 x (t) and σ 2 y (t), we have evaluated the dependence of the power-law exponents α x and α y on h y . In Fig. 4(c), we observe a monotonic but nonlinear decay of α x as h y increases, while in Fig. 4(d), we confirm the straightforward behavior of α y  x (t) and σ 2 y (t). Profile of (a) σ 2 x (t) and (b) σ 2 y (t) for nineteen values of h y homogeneously distributed between 0.05 and 0.95. In both cases, the dashes lines represent the behavior of the usual comb: σ 2 x (t) ∼ t αx with α x = 0.5 and σ 2 y (t) ∼ t αy with α y = 1. Here, we have omitted the finite-size-like effect on the variances for better visualization of the power-law behavior. For σ 2 x (t), we observe that the diffusion is enhanced (when compared with the usual comb) for an anti-persistent noise in y (h y < 0.5) and that is slows down when the noise in y is persistent. In the case of σ 2 y (t), anti-persistent noise slows down the diffusion while persistent noise enhances it. In (c) and (d), we show the dependence of α x and α y on h y . The values of α x and α y were obtained by least squares fitting a linear model to the log-log relationships between the variances and the time t. The error bars (which are very tiny) are the standard errors of α x and α y . When adjusting the data, we have selected only the values σ 2 x (t) and σ 2 y (t) within the interval 30 ≤ t ≤ 1600 in order to avoid the initial transient regime and the finite-size-like effect. We note the monotonic but nonlinear decay of α x in function of h y and the straightforward behavior of α y versus h y described by α y = 2h y (dashed line).
as function of h y (α y = 2h y ).
Another interesting question is whether the combinations of long-range correlated noises in the x and y directions affect the diffusive behaviors in a nontrivial way. To answer this question, we have built an ensemble of 10 5 particle positions for each possible pair of h x and h y homogeneously distributed between 0.05 and 0.95 with step size of 0.05. For each combination of h x and h y , we calculate the profile of the variances σ 2 x (t) and σ 2 y (t) and also the power-law exponents α x and α y describing the main tendency of these curves (by using the same procedure of Fig. 4). Figure 5(a) shows a contour plot of α x as a function of h x and h y and Fig. 5(b) represents the same for α y . For the Fig. 5(a), we note a quite rich diffusive scenery where, depending on the values of h x and h y , we may have subdiffusion slower (α x < 0.5) or faster (0.5 < α x < 1) than the usual comb, usual diffusion (α x = 1) and even superdiffusion (α x > 1). We further observe from this figure that α x cannot be written as a linear combination of h x and h y (which would be represented by straight lines in this plot) and thus, the correlations in x and y present a nontrivial interplay in the behavior of α x . On the other hand, the case of Fig 5(b) is rather simple because the correlations in x do not affect the diffusive behavior in the y−direction, this becomes evident by noting that the level curves of α y are horizontal straight lines. (α = 1) x Figure 5. (color online) The interplay between the long-range correlated noises in the x and y directions on the variances σ 2 x (t) and σ 2 y (t). Contour plots of (a) α x and (b) α y as a function of h x and h y . We note the rich diffusive scenery in the x−direction where we may have subdiffusion slower (α x < 0.5) or faster (0.5 < α x < 1) than the usual comb, usual diffusion (α x = 1) and even superdiffusion (α x > 1). We further observe that the power-law exponent α x cannot be expressed as a linear combination of h x and h y , suggesting a more complex coupling between the Hurst exponents and diffusion in the x−direction.
It is also important to investigate whether the long-range correlations affect the distribution profiles of the particle positions. In order to do so, we first calculate the marginal distributions P (ξ x ) and P (ξ y ) for there different values of h y (0.05, 0.50 and 0.95) and keeping h x = 0.50. Figures 6(a) and 6(b) show these distributions. We note that the distribution P (ξ x ) displays tails longer than those of usual comb (dark shaded region) when the noise in the y−direction is persistent (h y = 0.95) whereas, for antipersistent noise (h y = 0.05), the tails are shorter than those of usual comb and close to the Gaussian distribution (light shaded region). The distributions P (ξ y ) of the normalized positions ξ y are not affected by long-range correlated noises in the y−direction. We thus confirm that the long-range correlations in y have influence on the profiles of P (ξ x ). For a better and quantitative characterization of the role of the correlated noise η y (t) on the shape of P (ξ x ), we evaluate the coefficient of kurtosis κ x = ⟨ξ 4 x ⟩ ⟨ξ 2 x ⟩ 2 in function of h x and h y . For κ x > 3 the distribution is leptokurtic (peaked) while for κ x < 3 the distribution is platykurtic (flat) also, κ x = 3 for a Gaussian distribution and κ x ≈ 4.71 for the usual comb (distribution of the Eq. 5). Figure 6(c) shows a contour plot of κ x as a function of h x and h y , where we observe that the kurtosis ranges from very close to 3 (κ x = 3.04 ± 0.04 for h x = h y = 0.05) up to almost 13 (κ x = 12.55 ± 0.87 for h x = h y = 0.95). It is worth noting that, similar to the variance case [ Fig. 5(a)], the kurtosis κ x is not a trivial combination of h x and h y .  5)]. For P (ξ x ), we note that the distribution tails decay slower than those of usual comb when the noise in the y−direction is persistent (h y = 0.95); whereas for anti-persistent noise (h y = 0.05), the tails decay faster and approach the Gaussian distribution. In the case P (ξ y ), the long-range correlations in y have practically no effect on the profile of these distributions. (c) Contour plot of the coefficient of kurtosis κ x = ⟨ξ 4 x ⟩ ⟨ξ 2 x ⟩ 2 for the normalized position ξ x as a function of h x and h y . The value of κ x is an average over t within the interval 30 ≤ t ≤ 1600. We note that, depending on the values of h x and h y , the distribution P (ξ x ) exhibits different profiles, which can be mesokurtic (κ x = 3, as in the Gaussian case), leptokurtic (κ x > 3) or platykurtic (κ x < 3). We also observe that, similarly to variance case, κ x cannot be expressed in terms of linear combination of h x and h y .

Summary and Conclusions
We have proposed an extension for the comb model via Langevin equations driven by long-range correlated noises (fractional Gaussian noises). We initially showed that our Langevin equations are equivalent to the generalized diffusion equation (describing the usual comb) through the comparison of the numerical-obtained variances and marginal distributions with the analytical results predicted by the diffusion equation. Next, we have presented an extensive characterization of the diffusive properties of the particle positions in function of the Hurst exponents h x and h y defining the correlated noises. Our results show that the noise in the y−direction affects the variances and marginal distributions along the x−direction. Specifically, long-range persistence in y slows down the diffusion in x by reducing the power-law exponent α x and makes the tails of the marginal distribution P (ξ x ) longer as h y > 0.5 increases; whereas long-range antipersistence in y enhances the diffusion in x by increasing the value of α x and makes the tails of P (ξ x ) shorter as h y < 0.5 decreases. We have also investigated the interplay between the long-range correlations in the x and y directions, where we found that both α x and the kurtosis coefficient κ x of P (ξ x ) are not a trivial (such as a linear combination) function of h x and h y , but instead our results suggest a more complex coupling between the Hurst exponents and the diffusion in the x−direction. In summary, a quite rich diffusive scenery emerges from our model which may be directly applied for describing (at least in a first approximation) some of the recent empirical findings related to diffusion in complex systems such as living cells. Furthermore, our model combines and allows the study/analysis of the interplay between different mechanisms of anomalous diffusion, geometric constraints and long-range correlations, an important fact that had not been explored within the comb model framework and that may provide a better understanding of these processes as a combination of different mechanisms.