Analyze PFG Anomalous Diffusion via Real Space 2 and Phase Space Approaches 3

Pulsed-field gradient (PFG) diffusion experiments can be used to measure anomalous 7 diffusion in many polymer or biological systems. However, it is still complicated to analyze PFG 8 anomalous diffusion, particularly the finite gradient pulse width (FGPW) effect. In practical 9 applications, the FGPW effect may not be neglected such as in clinical diffusion magnetic resonance 10 imaging (MRI). Here, two significantly different methods are proposed to analyze PFG anomalous 11 diffusion: the effective phase shift diffusion equation (EPSDE) method and an observing the signal 12 intensity at the origin method. The EPSDE method describes the phase evolution in virtual phase 13 space, while the method to observe the signal intensity at the origin describes the magnetization 14 evolution in real space. However, these two approaches give the same general PFG signal 15 attenuation including FGPW effect, which can be numerically evaluated by a direct integration 16 method. The direct integration method is fast and without overflow. It is a convenient numerical 17 evaluation method for Mittag-Leffler function type PFG signal attenuation. The methods here 18 provide a clear view of spin evolution under field gradient, and their results will help the analysis 19 of PFG anomalous diffusion. 20


Introduction
Anomalous dynamic behavior exists in many polymer or biological systems [1,2,3,4,5].These systems often consist of various molecules such as macromolecules and small penetrant molecules like water.These molecules often demonstrate anomalous dynamics behavior: their rotational motion time constant usually obey a stretched exponential distribution, namely the Kohlrausch-Williams-Watts (KWW) function distribution [6,7], and their translational diffusion jump time and jump length could follow power law distributions [2].Nuclear magnetic resonance (NMR) is one of the important techniques to detect these anomalous dynamic behaviors.For instance, anomalous diffusion can be detected by pulsed-field gradient (PFG) NMR experiments.
The history of using field gradient to measure diffusion can be tracked back to Hahn, who first observed the influence of the molecule diffusion under gradient magnetic field upon echo amplitudes in 1950 [16].PFG diffusion measurements have broad applications in NMR and magnetic resonance imaging (MRI) [12][13][14][15].For instance, the water diffusion difference in different part of the brain can be used as an important contrast factor to build imaging of acute stroke [15].The PFG technique, being an ultra-valuable tool for normal disunion, could play an increasingly important role in monitoring anomalous diffusion occurring in many polymer and biological systems.
PFG anomalous diffusion [16,17,18,19,20,21,22] is different from PFG normal diffusion.Unlike normal diffusion, anomalous diffusion has a non-Gaussian probability distribution [2,23], and its mean square displacement is not linearly proportional to diffusion time [2].These non-Gaussian characteristics make it complicated to analyze PFG anomalous diffusion.Although PFG anomalous diffusion can be approximately analyzed by traditional methods such as the apparent diffusion coefficient method, PFG theories based on the fractional derivative could not only improve the analyzing accuracy on diffusion domain size, diffusion constant and other variables [ [16,179,24,25,26,27], but also yield additional information such as the time derivative order α and space derivative order β that are related with diffusion jump time and jump length distributions determined by material properties.
Many efforts have been devoted to studying PFG anomalous diffusion based on fractional calculus, which includes the propagator method [28], the modified-Bloch equation method [16,24,29], the effective phase shift diffusion equation method [17], the instantaneous signal attenuation method [25], the modified-Gaussian or non-Gaussian distribution method [26], etc.Additionally, PFG anomalous diffusions in restricted geometries such as plate, sphere, and cylinder have been investigated [27].These theoretical methods analyze PFG anomalous diffusion from different angles.
Therefore, each of them can have its own advantages in handling certain types of PFG anomalous diffusion.To better apply the PFG technique to study anomalous diffusion, it is still valuable to develop new theoretical treatments for PFG anomalous diffusion.
In this paper, two methods based on the fractional derivative [1,2,5,30,31,32] are proposed to give general analytical PFG signal attenuation expressions for anomalous diffusion.The first method is the recently developed effective phase shift diffusion equation (EPSDE) method [17].This method describes the spin phase evolution by an effective phase diffusion process in virtual phase space [17].
Solving the effective phase shift diffusion equation gives valuable information about the phase evolution process such as the phase probability distribution function and the moment of mean phase displacement.While, other conventional methods are difficult to get these phase information, and usually assume an approximate phase distribution such as Gaussian phase distribution [12].In this paper, it will be shown that a solvable PFG signal attenuation equation can be derived by applying Fourier transform to the effective phase shift equation.The second method is to observe the signal intensity at the origin method, which is an ultra-simple new method.For a homogeneous diffusion spin system, although the magnetization amplitude attenuates because of the gradient magnetic field effect, the phase of magnetization keeps constant at the origin of the gradient field.Such a specific phase property is employed to derive a PFG signal attenuation equation in this paper.The above two methods give the same signal attenuation equation, from which the general PFG signal attenuation expression can be derived by the Adomian decomposition method [33,34,35,36,37].Besides the Adomian decomposition method, a direct integration method was proposed for the numerical evaluation of the PFG signal attenuation, which is a fast and simple method.The results include finite gradient pulse width (FGPW) effect [12][13][14], namely the signal attenuation during each gradient pulse applying period.Theoretically, during a short gradient pulse the PFG signal attenuation can be neglected; nevertheless, the gradient pulse used in a clinical MRI is usually long [15,38 ].
Additionally, a longer gradient pulse allows the measuring of slower diffusion under the same gradient maximum intensity, which matters in studying polymer and biological systems where the molecule diffusion is often slow.Therefore, it needs to consider the FGPW effect in many real applications.The two methods here give the same results on the FGPW effect.These results agree with reported results from some other methods [17,25,26,] and the continuous time random walk simulation [25,28].Furthermore, PFG anomalous diffusion of intramolecular multiple quantum coherence (MQC) is also discussed [26].The MQC effect has the benefit to enhance the gradient effect on PFG signal attenuation [39].The two methods provide a compliment view of PFG anomalous diffusion from both the real space and phase space.

The phase space method-the effective phase shift diffusion method
For the sake of simplicity, only one-dimensional anomalous diffusion is considered in this paper.
In PFG experiments, the field gradient pulse results in a time and space dependent magnetic field , where 0 B is the exterior magnetic field, z is the position, and ) (t g is the timedependent gradient [12][13][14].The magnetic field exerts a torque on each spin moment.The torque changes the spin angular momentum direction, which leads the spin to precess about the magnetic field with Larmor frequency . In a rotating frame with angular frequency , a diffusing spin at position , and its phase accumulated along the diffusion path is [12][13][14]18]  where is the net accumulating phase.The range of )) ( cos( t φ is the projection factor of the spin magnetization to the observing coordinate axis.The PFG signal comes from the ensemble contribution from all spins by averaging over all possible phases [12][13][14] Because the NMR signal comes from numerous spins, the percentage of spins with potential infinite  3), but the phase evolution process based on the phase path integral equation ( 1) is complicated.Nevertheless, the recently developed effective phase shift diffusion equation method provides a simple way to describe the phase evolution based on an effective phase diffusion process [17].In the following, the effective phase diffusion equation method [17] will be briefly reintroduced.Additionally, a general solution including the finite gradient pulse width effect will be given that has not been reported in [17].
For simplicity, only the diffusion with symmetric probability distribution is studied here.The self-diffusion process can be described by a random walk, which consists of a sequence of independent random jumps with waiting times , and corresponding displacement lengths Based on the random walk, Equation (1) can be rewritten as [17] [ ] [ ] where is the wavenumber [12][13][14].In Equation (4), The term in Equation ( 6) is an effective phase random walk with a jump length and a jump waiting time m t Δ .The phase random walk can be treated as an effective phase diffusion, which can be obtained from the corresponding spin diffusion in the real space by scaling the jump length , the effective phase diffusion constant is ) [17], where where is the Caputo fractional derivative defined as [5,[30][31][32]] and is the space fractional derivative [5,[30][31][32] defined in Appendix A. By replacing coordinate z with φ , and respectively from real space diffusion equation ( 9), the effective phase diffusion equations can be obtained [17] as As ( ) [30][31][32] where q is wavenumber, applying Fourier transform on both sides of Equation ( 9) gives which is a q space phase diffusion equation.Here, q is the wavenumber dedicated to the above Fourier transform, which differs from the field gradient pulse induced wavenumber ) (t K defined by Equation ( 5).When q = 1 , ) , ( t q P in Equation ( 10) equals the attenuated signal amplitude because in PFG experiments the signal attenuation can be calculated based on Equation (2) as [17] ( ) The normalized value of ) 0 ( S equaling 1 will be used and dropped out throughout this paper. Substituting Equation (11) and q = 1 into Equation ( 10), we get which is the PFG signal attenuation equation.
The solution of Equation ( 12) will be given in section 2.3 as the same equation will be obtained by the observing the signal intensity at the origin method in the subsequent section.
and the effective phase diffusion constant is also demonstrated below PGSE pulse sequence.

Observing the signal intensity at the origin in real space
For a non-diffusing spin, in a rotating frame with angular frequency , its precession phase can be described as [12][13][14] where ) , ( t z ψ is the precession phase, which is different from the accumulated phase ) (t φ described by (1).The non-diffusing spin system has a time and space modulated magnetization, which can be described as ( ) where is the magnetization amplitude which can be referred to as signal intensity, and the ( ) is the phase term.For a pulsed gradient spin echo (PGSE) or the pulsed gradient stimulated-echo (PGSTE) experiment with a constant gradient as shown in Figure 1, the wavenumber is where δ is the gradient pulse length, and Δ is the diffusion delay starting from the beginning of the first dephasing gradient pulse to the end of the last rephrasing gradient pulse.Because of the two counteracting gradient pulses-a dephasing pulse and a rephasing pulse, ) (t K returns to 0 at the end of the rephasing pulse, and the phase of non-diffusing spins will be refocused.Therefore, the signal of the non-diffusing spin system does not attenuate, when the T2 relaxation is neglected. For the spin diffusion in a homogeneous sample, its magnetization still can be described by Equation ( 14), ( ) . At a random position z , the possibilities of spins diffusing to opposite directions z z Δ ± are equal [25].The opposite movements yield , which has no effect on the phase, ( ) , but does make the signal intensity decay by a factor ( ) for these spins.By substituting ( ) into diffusion Equation ( 7), and applying At the origin, , and ( ) , we thus have Equation ( 17) is identical to Equation ( 12) obtained by the effective phase shift diffusion equation method in the previous section 2.1.
This observing signal intensity at the origin method could be understood in another way: if only the signal intensity or amplitude is observed from a selected position at different diffusion time in a PFG anomalous diffusion experiment, the observed signal intensity or amplitude shall still obey Equation ( 17) regardless whether the phase is observed or not.It is reasonable that both ways can give the same signal attenuation expression, because the signal amplitude is homogeneous throughout the sample at each instant, and selecting the origin at a random position in the sample does not affect the result.

Analytical solution by the Adomian decomposition method
The same equation as Equation ( 12) or (17) has been obtained by the fractional integral modified-Bloch equation method.The three different methods get the same PFG signal attenuation equation.
Reference [29] has employed the Adomian decomposition method to give a general analytical solution to Equation (12).In the following, the general solution from the Adomian decomposition will be given in more detail for the cases of normal diffusion and general anomalous diffusion.

normal diffusion.
First, the signal attenuation equation, Equation ( 12) or ( 17) can be used to obtain the PFG signal attenuation expression for normal diffusion, which is a specific case of anomalous diffusion with where D is the diffusion constant of normal diffusion.The solution of Equation ( 18) is which is the PFG signal attenuation expression for normal diffusion.The same result as expression (19) has been obtained in [40].For PGSE or PGSTE experiments, the PFG signal attenuation calculated based on Equation ( 19) is , which is a routinely used expression for PFG normal diffusion [12][13][14].

general anomalous diffusion.
Second, the PFG signal attenuation for general anomalous diffusion can be obtained.The same equation as Equation ( 17) and its solution have been obtained by modified-Bloch equation method in [29].According to [29], Equations ( 12) or ( 17) can be written equivalently as ( ) where Based on the Adomian decomposition method [33][34][35][36][37], the solution of Equation ( 20) is [29]  where In PGSE or PGSTE experiments, the time t can be separated into three periods: , is used [30][31][32] for a free diffusion in a homogeneous sample, we can get the following: (a) PFG signal attenuation under SGP approximation: δ is short enough and the diffusion inside each gradient pulse can be neglected.We get where ) . Equation (23) replicates the SGP approximation result obtained in references [17,25,26].
(b) Single pulse attenuation: this is an ideal situation, the first gradient pulse is regular, but the second gradient pulse is infinitely narrow whose purpose is to counteract the first gradient pulse.We get where 1 ) ( is a Mittag-Leffler type function [ 41 ].Equation ( 24) is consistent with the results obtained by the modified Bloch equation proposed in [24] (c) General PFG signal attenuation: the PGSE or PGSTE experiment includes three periods: The integration in Equation (22c) during these three periods is tedious, which can be calculated with computer assistance.Nevertheless, we can get the first and second terms of Equation (22a) as the following: where B(x,y) and B(a;x,y) are the Beta function and incomplete Beta function.When Equation ( 26) agrees with the signal attenuation expression obtained by the instantaneous signal attenuation method [25].At small signal attenuation, can be approximately expanded as Equation ( 27) is the same as that given by the combination of Equations (25a) and (26) (Note ( ) ( ) ), however, at large signal attenuation, may deviate from the combination of Equations (25a) and (26) as the higher order terms in the expansion cannot be omitted.This agreement can be explained by the following: in typical PGSE or PGSTE experiments, both the gradient intensity and pulse length of the dephasing pulse and the rephrasing pulse are identical, , the general anomalous diffusion reduces to time-fractional diffusion and Equation (25b) can be further calculated as , the general anomalous diffusion reduces to space-fractional diffusion.

Numerical evaluation of Mittag-Leffler function based PFG signal attenuation by the direct integration method
Although the Adomian decomposition method can be used to numerically evaluate the PFG signal attenuation, its calculation speed is slow and it could cause overflow at large signal attenuation.References [29,42 ] proposed an alternative method, a direct integration method that can give the same numerical results, but with a much fast speed and without overflow.From Equation ( 20), if we set where 30) can be further written as , which is a Mittag-Leffler function; when 0 , . Hence, we have ( ) ( ) Equation ( 30) can be rewritten in a discrete form as where [ ] This method is simple, and it has no overflow issue in the calculation.

CTRW in a lattice model
Random walk simulation is a powerful numerical method that employs a stochastic jump process to model normal and anomalous diffusion in physics, chemistry, biology and many other disciplines.It can be used to simulate the PFG signal attenuation in mathematically tractable or intractable systems.In this paper, the CTRW simulation method developed in [25] is used to verify the theoretical results, which is based on two models: the CTRW model [43] and the Lattice model [44 ,45 ].In the simulation, a sequence of independent random waiting time and jump length combinations ( ) is produced by computer program.
The individual waiting time t Δ and jump length ξ Δ , are given according to [25] by the following expressions: and where t η and z η are scale constants, are two independent random numbers.The CTRW simulation based on the above waiting time and jump length satisfies the timespace fractional diffusion equation under the diffusive limit, which provides a simple way to simulate anomalous diffusion in various research areas, such as physics, and economics [43].
Although continuous waiting time and jump length are used in the simulation, the accumulating spin phase associated with the diffusion path is recorded in a discrete manner.Such a discrete phase recording manner is convenient and reasonable.The simulation can be viewed as a numerical pulsedfield gradient experiment, whose observables such as phase can be observed in a discrete time selected by an experiment observer; the observing or recording manner will not affect the fundamental numerical experiment process, namely the producing of the continuous random walk sequence: ).The spin phase was recorded by the lattice model developed in [27,28], which has been applied to simulate PFG diffusion in polymer system [46].The spin phase in the simulation is ( ) ( ) where , and j t ′ ′ is the discrete record time, which takes place between the l th and l+1 th steps of the random walk.The second term ( ) on the right-hand side of Equation ( 37) corresponds to the partial phase evolution of the 1 + l th jump step that needs to be recorded at the time j t ′ ′ .The PFG signal attenuation in the simulation can be obtained by averaging over all the walkers in the simulation [25,45,46] [ ] [ ] where walks N is the total number of the walks.The total walks in each simulation are 1,000,000.
Because the CTRW model proposed in [25] is only for subdiffusion, the simulation here is limited to the subdiffusion.Interested readers are referred to references [25,29,[43][44][45][46] for more detailed information.

Results and discussion
Table 1.Comparison PFG anomalous diffusion results by the effective phase shift diffusion method, the observing the signal intensity at the origin method and other methods.

Virtual phase space diffusion Real space spin diffusion Effective phase shift diffusion equation
Observing the signal intensity at the origin PFG signal attenuation expression by the Adomian decomposition method [29]*: Under SGP approximation: At small attenuation: Numerical evaluation by the direct integration method [42] [17,25,26] *The PFG signal attenuation equation and expression obtained in this paper are the same as that obtained by the integral modified-Bloch equation method [29].

[ ] [ ]
This paper uses two different methods to describe PFG anomalous diffusion: the effective phase shift diffusion equation method [17] and the observing the signal intensity at the origin method.The major results in this paper are summarized in Table 1.The EPSDE method [17] uses an effective phase diffusion process to describe the phase evolution, while the observing the signal intensity at the origin method only monitors the real space signal amplitude change at the origin where the phase of the magnetization keeps constant.Each of these two approaches has its unique advantages.The effective phase shift diffusion equation method can directly obtain the phase distribution, which greatly simplifies the PFG diffusion analyzing [17]; in contrast, the traditional methods often need to get the real space spin particle probability distribution first, which is then used to obtain the phase distribution.The ultra-simple observing signal intensity at the origin method obtains the signal attenuation equation by substituting the magnetization ( ) directly into the diffusion equation and using the phase property  These two approaches view the spin evolution process from two different spaces−the real space and the phase space, which provides a clear picture of spin dynamics in PFG diffusion experiments.
These two methods agree with each other.They get the same PFG signal attenuation equation, (11) or (16), which can be solved by the Adomian decomposition method.The solution gives PFG signal attenuation expressions (21a)-(21c), which include the finite gradient pulse width effect.
Understanding the FGPW effect is important as the clinical MRI applications often use long gradient pulse.Additionally, using long gradient pulses allows researchers to monitor slower diffusion in polymer or biological systems where the molecules or ions often diffuse slowly.The consistent results from two different methods help us better understand the FGPW effect.
Additionally, the two methods agree with other methods.First, the PFG signal attenuation equation ( 17) obtained here is the same as that obtained by the integral modified-Bloch equation method [29].Second, the results agree with the CTRW simulation as shown in Figure 2. In Figure 2a, the fractional diffusion constant f D = 6.6 × 10 -11 m β /s α was obtained from fitting the curve of mean square displacement versus t by equation Interested readers can be referred to [29] for additional comparison between the theoretical prediction and the CTRW simulation.Third, the PFG signal attenuation expressions (22 a-c) are also consistent with those obtained by various approximation methods such as the non-Gaussian phase distribution (nGPD) method [26], and instantaneous signal attenuation method [25].The same SGP approximation results, Equation ( 23) can be obtained by all these five different methods: the EPSDE method [17], the instantaneous signal attenuation method [25], the nGPD method [26], the method to observe the signal intensity at the origin, and the modified-Bloch equation method [29].Moreover, Equation ( 27), can be approximately obtained in this paper, which agrees with that obtained from the instantaneous signal attenuation method [25].The numerical evaluation of PFG signal attenuation by Equation ( 20) or ( 30) can be conveniently performed by the direct integration method.Figure 2b shows that the results of the direct integration method agree with that obtained from the Adomian decomposition method.Compared to the Adomian decomposition method, the computing speed of the direct integration method improves by orders of magnitudes.The speed improvement is because in the direct integration method, each in Equation ( 32) is only a single time integration, while in the Adomian decomposition method, is a superposition of many terms of time integration.Additionally, from Equations ( 33) and (34), the direct integration method can be used to calculate Mittag-Leffler type functions and their derivative.Figure 3 shows that the Mittag-Leffler function (MLF) calculated from the direct integration method agree with that calculated by other methods in [47,48].The Pade approximation method in [48] is only for subdiffusion.Because the direct integration method does not cause overflow, which can be a useful method to calculate Mittag-Leffler type functions.The Fortran code for MLF calculation can be obtained from the link: https://github.com/GLin2017/Mittag-Lefflerfunction-calculated-by-Direct-Integration.This direct integration method will be a great help to the application of current theoretical results to PFG anomalous diffusion in NMR and MRI.integration method, the Pade approximation method [48], and the method used in [47].
In the current results, only PFG anomalous diffusion of single quantum coherence is considered.
The results here can be easily extended to handle PFG anomalous diffusion of intramolecular multiple quantum coherence (MQC) [26,38].As the gradient field induced phase evolution of an n order intramolecular MQC is n times faster than the corresponding single quantum coherence, the PFG signal attenuation of intramolecular MQC is the same as that of single quantum coherence with an effective gradient intensity ng [40].Therefore, for intramolecular MQC in PGSE or PGSTE experiments, the PFG signal attenuation equation is ( From Equation (36), it is obvious that the This paper only considers the spin self-diffusion that can be described by the time-space fractional diffusion equation based on the fractional derivative.In general, the two methods used in this paper can be applied to other types of anomalous diffusions such as that described by the timespace diffusion equation based on the fractal derivative [49], etc.Additionally, only the symmetric anomalous diffusion in homogeneous spin systems is studied here.In real applications, the anomalous diffusion can take place in complicated systems such as inhomogeneous system, restricted geometries [27], anisotropic system [42], non-symmetric system, etc. Additonally, the gradient field may be nonlinear [50].The current methods may face challenges in these complicated situations, which reminds us that much effort is required to the studying of PFG anomalous diffusion.
diffusion coefficients for the real space diffusion and the effective phase diffusion respectively, and their units are effective phase diffusion and the real spin particle diffusion belong to the same type diffusion, they should obey the same type diffusion equation.The one-dimensional real space diffusion equation based on the fractional derivative is[5,[30][31][32]

Figure 1 .
Figure 1.(a) PGSTE pulse sequences; (b) PGSE pulse sequence.The gradient pulse width is δ , and the diffusion delay is Δ.The time-dependent behavior of wavenumber

.
Figure 2b shows the PFG signal attenuation of anomalous diffusion obtained by Equations (22a-c) based on the Adomian decomposition method, Equation (32) based on the direct integration method, and the CTRW simulation.Other parameters used in Figure 2 are 2

Figure 2 . 6 .
Figure 2. Comparison of the PFG signal attenuation from theoretical predictions with CTRW simulation: (a) mean square displacement from CTRW simulation, (b) PFG signal attenuation from different methods: the Adomian decomposition method with Equations (22a-c), the direct integration method with Equation (32), the and the CTRW simulation.The parameters used are

Figure 3 .
Figure 3. Comparing the MLF calculated by different methods: Equation (33) from the direct t