A numerical method for the fractional Fitzhugh–Nagumo monodomain model

A fractional FitzHugh–Nagumo monodomain model with zero Dirichlet boundary conditions is presented, generalising the standard monodomain model that describes the propagation of the electrical potential in heterogeneous cardiac tissue. The model consists of a coupled fractional Riesz space nonlinear reaction-diffusion model and a system of ordinary differential equations, describing the ionic fluxes as a function of the membrane potential. We solve this model by decoupling the space-fractional partial differential equation and the system of ordinary differential equations at each time step. Thus, this means treating the fractional Riesz space nonlinear reaction-diffusion model as if the nonlinear source term is only locally Lipschitz. The fractional Riesz space nonlinear reaction-diffusion model is solved using an implicit numerical method with the shifted Grunwald–Letnikov approximation, and the stability and convergence are discussed in detail in the context of the local Lipschitz property. Some numerical examples are given to show the consistency of our computational approach. References B. Baeumer, M. Kovaly, and M. M. Meerschaert, Fractional reproduction-dispersal equations and heavy tail dispersal kernels, Bulletin of Mathematical Biology 69:2281–2297, 2007. doi:10.1007/s11538-007-9220-2 B. Baeumer, M. Kovaly, and M. M. Meerschaert, Numerical solutions for fractional reaction-diffusion equations, Computers and Mathematics with Applications 55:2212–2226, 2008. doi:10.1016/j.camwa.2007.11.012 N. Badie and N. Bursac, Novel micropatterned cardiac cell cultures with realistic ventricular microstructure, Biophys J 96:3873–3885, 2009. doi:10.1016/j.bpj.2009.02.019 A. Bueno-Orovio, D. Kay, K. Burrage, Fourier spectral methods for fractional-in-space reaction-diffusion equations, Technical report, University of Oxford, 2013. A. Bueno-Orovioy, D. Kay, V. Grau, B. Rodriguez and K. Burrage, Fractional dffusion models of electrical propagation in cardiac tissue: electrotonic effects and the modulated dispersion of repolarization, Technical report, University of Oxford, 2013. K. F. Decker, J. Heijman, J. R. Silva, T. J. Hund and Y. Rudy, Properties and ionic mechanisms of action potential adaptation, restitution, and accommodation in canine epicardium, Am. J. Physiol Heart Circ. Physiol. 296:H1017–H1026, 2009. doi:10.1152/ajpheart.01216.2008 J. S. Frank and G. A. Langer, The myocardial interstitium: its structure and its role in ionic exchange, J Cell Biol 60:586–601, 1974. doi:10.1083/jcb.60.3.586 A. L. Hodgkin and A. F. Huxley, A quantitative description of membrane current and its application to conduction and excitation in nerve, J. Physiol. (Lond), 117:500–544, 1952. http://jp.physoc.org/content/117/4/500.html R. FitzHugh, Impulses and Physiological States in Theoretical Models of Nerve Membrane, Biophys. J., 1:445–466, 1961. doi:10.1016/S0006-3495(61)86902-6 D. Kay, I. W. Turner, N. Cusimano and K. Burrage, Reflections from a boundary: reflecting boundary conditions for space-fractional partial differential equations on bounded domains, Technical report, University of Oxford, 2013. . F. Liu, V. Anh and I. Turner, Numerical solution of space fractional Fokker-Planck equation. J. Comp. and Appl. Math., 166:209–219, 2004. doi:10.1016/j.cam.2003.09.028 F. Liu, P. Zhuang, V. Anh and I. Turner and K. Burrage, Stability and convergence of the difference methods for the space-time fractional advection-diffusion equation. Appl. Math. Comp., 191:12–20, 2007. doi:10.1016/j.amc.2006.08.162 R. Magin, O. Abdullah, D. Baleanu and X. Zhou, Anomalous diffusion expressed through fractional order differential operators in the Bloch–Torrey equation, Journal of Magnetic Resonance 190:255–270, 2008. doi:10.1016/j.jmr.2007.11.007 M. M. Meerschaert, J. Mortensenb and S. W. Wheatcraft, Fractional vector calculus for fractional advection-dispersion, Physica A, 367:181–190, 2006. doi:10.1016/j.physa.2005.11.015 L. C. McSpadden, R. D. Kirkton and N. Bursac, Electrotonic loading of anisotropic cardiac monolayers by unexcitable cells depends on connexin type and expression level, Am. J. Physiol. Cell Physiol. 297:C339–C351, 2009. doi:10.1152/ajpcell.00024.2009 J. Nagumo, S. Arimoto, and S. Yoshizawa, An active pulse transmission line simulating nerve axon, Proceedings of the IRE, 50:2061–2070, 1962. doi:10.1109/JRPROC.1962.288235 S. F. Roberts, J. G. Stinstra and C. S. Henriquez, Effect of nonuniform interstitial space properties on impulse propagation: a discrete multidomain model, Biophys J 95:3724–3737, 2008. doi:10.1529/biophysj.108.137349 J. Sundnes, G. T. Lines, X. Cai, B. F. Nielsen, K. A. Mardal and A. Tveitio, Computing the electrical activity in the heart, Springer-Verlag, 2006. G. D. Smith, Numerical Solution of Partial Differential Equations: Finite Difference Methods, Clarendon Press, Oxford, 1985. F. J. Valdes-Parada, J. A. Ochoa-Tapia and J. Alvarez-Ramirez, Effective medium equations for fractional Fick law in porous media, Physica A, 373:339–353, 2007. doi:10.1016/j.physa.2006.06.007 Q. Yang, F. Liu and I. Turner, Stability and convergence of an effective numerical method for the time-space fractional Fokker-Planck equation with a nonlinear source term, International Journal of Differential Equations, 2010:464321, 2010, doi:10.1155/2010/464321 W. Ying, A multilevel adaptive approach for computational cardiology, PhD thesis, Duke University, 2005. Q. Yu, F. Liu, I. Turner and K. Burrage, A computationally effective alternating direction method for the space and time fractional Bloch-Torrey equation in 3-D, Appl. Math. Comp., 219:4082–4095, 2012. doi:10.1016/j.amc.2012.10.056 Q. Yu, F. Liu, I. Turner and K. Burrage, Stability and convergence of an implicit numerical method for the space and time fractional Bloch-Torrey equation, the special issue of Fractional Calculus and Its Applications in-Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 371:20120150, 2013. doi:10.1098/rsta.2012.0150 Q. Yu, F. Liu, I. Turner and K. Burrage, Numerical simulation of the fractional Bloch equations, J. Comp. Appl. Math., 255:635–651, 2014. doi:10.1016/j.cam.2013.06.027 P. Zhuang, F. Liu, V. Anh and I. Turner, Numerical methods for the variable-order fractional advection-diffusion equation with a nonlinear source term, SIAM J. Num. Anal., 47:1760–1781, 2009. doi:10.1137/08073059

High‐order ADI orthogonal spline collocation method for a new 2D fractional integro‐differential problem

This is the peer reviewed version of the following article: Qiao L, Xu D, Yan Y. (2020). High-order ADI orthogonal spline collocation method for a new 2D fractional integro-differential problem. Mathematical Methods in the Applied Sciences, 1-17., which has been published in final form at https://doi.org/10.1002/mma.6258. This article may be used for non-commercial purposes in accordance with Wiley Terms and Conditions for Use of Self-Archived Versions.We use the generalized L1 approximation for the Caputo fractional deriva-tive, the second-order fractional quadrature rule approximation for the inte-gral term, and a classical Crank-Nicolson alternating direction implicit (ADI)scheme for the time discretization of a new two-dimensional (2D) fractionalintegro-differential equation, in combination with a space discretization by anarbitrary-order orthogonal spline collocation (OSC) method. The stability of aCrank-Nicolson ADI OSC scheme is rigourously established, and error estimateis also derived. Finally, some numerical tests are give

Can anomalous diffusion models in magnetic resonance imaging be used to characterise white matter tissue microstructure?

During the time window of diffusion weighted magnetic resonance imaging experiments (DW-MRI), water diffusion in tissue appears to be anomalous as a transient effect, with a mean squared displacement that is not a linear function of time. A number of statistical models have been proposed to describe water diffusion in tissue, and parameters describing anomalous as well as Gaussian diffusion have previously been related to measures of tissue microstructure such as mean axon radius. We analysed the relationship between white matter tissue characteristics and parameters of existing statistical diffusion models.A white matter tissue model (ActiveAx) was used to generate multiple b-value diffusion-weighted magnetic resonance imaging signals. The following models were evaluated to fit the diffusion signal: 1) Gaussian models - 1a) mono-exponential decay and 1b) bi-exponential decay; 2) Anomalous diffusion models - 2a) stretched exponential, 2b) continuous time random walk and 2c) space fractional Bloch-Torrey equation. We identified the best candidate model based on the relationship between the diffusion-derived parameters and mean axon radius and axial diffusivity, and applied it to the in vivo DW-MRI data acquired at 7.0 T from five healthy participants to estimate the same selected tissue characteristics. Differences between simulation parameters and fitted parameters were used to assess accuracy and in vivo findings were compared to previously reported observations.The space fractional Bloch-Torrey model was found to be the best candidate in characterising white matter on the base of the ActiveAx simulated DW-MRI data. Moreover, parameters of the space fractional Bloch-Torrey model were sensitive to mean axon radius and axial diffusivity and exhibited low noise sensitivity based on simulations. We also found spatial variations in the model parameter β to reflect changes in mean axon radius across the mid-sagittal plane of the corpus callosum.Simulations have been used to define how the parameters of the most common statistical magnetic resonance imaging diffusion models relate to axon radius and diffusivity. The space fractional Bloch-Torrey equation was identified as the best model for the characterisation of axon radius and diffusivity. This model allows changes in mean axon radius and diffusivity to be inferred from spatially resolved maps of model parameters

Computational modelling of diffusion magnetic resonance imaging based on cardiac histology

The exact relationship between changes in myocardial microstructure as a result of heart disease and the signal measured using diffusion tensor cardiovascular magnetic resonance (DT-CMR) is currently not well understood. Computational modelling of diffusion in combination with realistic numerical phantoms offers the unique opportunity to study effects of pathologies or the efficacy of improvements to acquisition protocols in a controlled in-silico environment. In this work, Monte Carlo random walk (MCRW) methods are used to simulate diffusion in a histology-based 3D model of the myocardium. Sensitivity of typical DT-CMR sequences to changes in tissue properties is assessed. First, myocardial tissue is analysed to identify important geometric features and diffusion parameters. A two-compartment model is considered where intra-cellular compartments with a reduced bulk diffusion coefficient are separated from extra-cellular space by permeable membranes. Secondary structures like groups of cardiomyocyte (sheetlets) must also be included, and different methods are developed to automatically generate realistic histology-based substrates. Next, in-silico simulation of DT-CMR is reviewed and a tool to generate idealised versions of common pulse sequences is discussed. An efficient GPU-based numerical scheme for obtaining a continuum solution to the Bloch--Torrey equations is presented and applied to domains directly extracted from histology images. In order to verify the numerical methods used throughout this work, an analytical solution to the diffusion equation in 1D is described. It relies on spectral analysis of the diffusion operator and requires that all roots of a complex transcendental equation are found. To facilitate a fast and reliable solution, a novel root finding algorithm based on Chebyshev polynomial interpolation is proposed. To simulate realistic 3D geometries MCRW methods are employed. A parallel simulator for both grid-based and surface mesh--based geometries is presented. The presence of permeable membranes requires special treatment. For this, a commonly used transit model is analysed. Finally, the methods above are applied to study the effect of various model and sequence parameters on DT-CMR results. Simulations with impermeable membranes reveal sequence-specific sensitivity to extra-cellular volume fraction and diffusion coefficients. By including membrane permeability, DT-CMR results further approach values expected in vivo.Open Acces

Development of Methodologies for Diffusion-weighted Magnetic Resonance Imaging at High Field Strength

Diffusion-weighted imaging of small animals at high field strengths is a challenging prospect due to its extreme sensitivity to motion. Periodically rotated overlapping parallel lines with enhanced reconstruction (PROPELLER) was introduced at 9.4T as an imaging method that is robust to motion and distortion. Proton density (PD)-weighted and T2-weighted PROPELLER data were generally superior to that acquired with single-shot, Cartesian and echo planar imaging-based methods in terms of signal-to-noise ratio (SNR), contrast-to-noise ratio and resistance to artifacts. Simulations and experiments revealed that PROPELLER image quality was dependent on the field strength and echo times specified. In particular, PD-weighted imaging at high field led to artifacts that reduced image contrast. In PROPELLER, data are acquired in progressively rotated blades in k-space and combined on a Cartesian grid. PROPELLER with echo truncation at low spatial frequencies (PETALS) was conceived as a postprocessing method that improved contrast by reducing the overlap of k-space data from different blades with different echo times. Where the addition of diffusion weighting gradients typically leads to catastrophic motion artifacts in multi-shot sequences, diffusion-weighted PROPELLER enabled the acquisition of high quality, motion-robust data. Applications in the healthy mouse brain and abdomen at 9.4T and in stroke patients at 3T are presented. PROPELLER increases the minimum scan time by approximately 50%. Consequently, methods were explored to reduce the acquisition time. Two k-space undersampling regimes were investigated by examining image fidelity as a function of degree of undersampling. Undersampling by acquiring fewer k-space blades was shown to be more robust to motion and artifacts than undersampling by expanding the distance between successive phase encoding steps. To improve the consistency of undersampled data, the non-uniform fast Fourier transform was employed. It was found that acceleration factors of up to two could be used with minimal visual impact on image fidelity. To reduce the number of scans required for isotropic diffusion weighting, the use of rotating diffusion gradients was investigated, exploiting the rotational symmetry of the PROPELLER acquisition. Fixing the diffusion weighting direction to the individual rotating blades yielded geometry and anisotropy-dependent diffusion measurements. However, alternating the orientations of diffusion weighting with successive blades led to more accurate measurements of the apparent diffusion coefficient while halving the overall acquisition time. Optimized strategies are proposed for the use of PROPELLER in rapid high resolution imaging at high field strength

Monte Carlo Simulation of Diffusion Magnetic Resonance Imaging

The goal of this thesis is to describe, implement and analyse Monte Carlo (MC) algorithms for simulating the mechanism of diffusion magnetic resonance imaging (dMRI). As the inverse problem of mapping the sub-voxel micro-structure remains challenging, MC methods provide an important numerical approach for creating ground-truth data. The main idea of such simulations is first generating a large sample of independent random trajectories in a prescribed geometry and then synthesizing the imaging signals according to given imaging sequences. The thesis starts by providing a concise introduction of the mathematical background for understanding dMRI. It then proceeds to describe the workflow and implementation of the most basic Monte Carlo method with experiments performed on simple geometries. A theoretical framework for error analysis is introduced, which to the best of the author's knowledge, has been absent in the literature. In an effort to mitigate the costly nature of MC algorithms, the geometrically adaptive fast random walk algorithm (GAFRW) is implemented, first invented by D.Grebenkov. Additional mathematical justification is provided in the appendix should the reader find details in the original paper by Grebenkov lacking. The result suggests that the GAFRW algorithm only provides moderate accuracy improvement over the crude MC method in the geometry modeled after white matter fibers. Overall, both approaches are shown to be flexible for a variety of geometries and pulse sequences

Spatially Regularized Spherical Reconstruction: A Cross-Domain Filtering Approach for HARDI Signals

Despite the immense advances of science and medicine in recent years, several aspects regarding the physiology and the anatomy of the human brain are yet to be discovered and understood. A particularly challenging area in the study of human brain anatomy is that of brain connectivity, which describes the intricate means by which different regions of the brain interact with each other. The study of brain connectivity is deeply dependent on understanding the organization of white matter. The latter is predominantly comprised of bundles of myelinated axons, which serve as connecting pathways between approximately 10¹¹ neurons in the brain. Consequently, the delineation of fine anatomical details of white matter represents a highly challenging objective, and it is still an active area of research in the fields of neuroimaging and neuroscience, in general. Recent advances in medical imaging have resulted in a quantum leap in our understanding of brain anatomy and functionality. In particular, the advent of diffusion magnetic resonance imaging (dMRI) has provided researchers with a non-invasive means to infer information about the connectivity of the human brain. In a nutshell, dMRI is a set of imaging tools which aim at quantifying the process of water diffusion within the human brain to delineate the complex structural configurations of the white matter. Among the existing tools of dMRI high angular resolution diffusion imaging (HARDI) offers a desirable trade-off between its reconstruction accuracy and practical feasibility. In particular, HARDI excels in its ability to delineate complex directional patterns of the neural pathways throughout the brain, while remaining feasible for many clinical applications. Unfortunately, HARDI presents a fundamental trade-off between its ability to discriminate crossings of neural fiber tracts (i.e., its angular resolution) and the signal-to-noise ratio (SNR) of its associated images. Consequently, given that the angular resolution is of fundamental importance in the context of dMRI reconstruction, there is a need for effective algorithms for de-noising HARDI data. In this regard, the most effective de-noising approaches have been observed to be those which exploit both the angular and the spatial-domain regularity of HARDI signals. Accordingly, in this thesis, we propose a formulation of the problem of reconstruction of HARDI signals which incorporates regularization assumptions on both their angular and their spatial domains, while leading to a particularly simple numerical implementation. Experimental evidence suggests that the resulting cross-domain regularization procedure outperforms many other state of the art HARDI de-noising methods. Moreover, the proposed implementation of the algorithm supersedes the original reconstruction problem by a sequence of efficient filters which can be executed in parallel, suggesting its computational advantages over alternative implementations

Doctor of Philosophy

dissertationDiffusion tensor MRI (DT-MRI or DTI) has been proven useful for characterizing biological tissue microstructure, with the majority of DTI studies having been performed previously in the brain. Other studies have shown that changes in DTI parameters are detectable in the presence of cardiac pathology, recovery, and development, and provide insight into the microstructural mechanisms of these processes. However, the technical challenges of implementing cardiac DTI in vivo, including prohibitive scan times inherent to DTI and measuring small-scale diffusion in the beating heart, have limited its widespread usage. This research aims to address these technical challenges by: (1) formulating a model-based reconstruction algorithm to accurately estimate DTI parameters directly from fewer MRI measurements and (2) designing novel diffusion encoding MRI pulse sequences that compensate for the higher-order motion of the beating heart. The model-based reconstruction method was tested on undersampled DTI data and its performance was compared against other state-of-the-art reconstruction algorithms. Model-based reconstruction was shown to produce DTI parameter maps with less blurring and noise and to estimate global DTI parameters more accurately than alternative methods. Through numerical simulations and experimental demonstrations in live rats, higher-order motion compensated diffusion-encoding was shown to successfully eliminate signal loss due to motion, which in turn produced data of sufficient quality to accurately estimate DTI parameters, such as fiber helix angle. Ultimately, the model-based reconstruction and higher-order motion compensation methods were combined to characterize changes in the cardiac microstructure in a rat model with inducible arterial hypertension in order to demonstrate the ability of cardiac DTI to detect pathological changes in living myocardium