Practical computation of the diffusion MRI signal of realistic neurons based on Laplace eigenfunctions

Abstract

The complex transverse water proton magnetization subject to diffusion‐encoding magnetic field gradient pulses in a heterogeneous medium such as brain tissue can be modeled by the Bloch‐Torrey partial differential equation. The spatial integral of the solution of this equation in realistic geometry provides a gold‐standard reference model for the diffusion MRI signal arising from different tissue micro‐structures of interest. A closed form representation of this reference diffusion MRI signal called matrix formalism, which makes explicit the link between the Laplace eigenvalues and eigenfunctions of the biological cell and its diffusion MRI signal, was derived 20 years ago. In addition, once the Laplace eigendecomposition has been computed and saved, the diffusion MRI signal can be calculated for arbitrary diffusion‐encoding sequences and b‐values at negligible additional cost. Up to now, this representation, though mathematically elegant, has not been often used as a practical model of the diffusion MRI signal, due to the difficulties of calculating the Laplace eigendecomposition in complicated geometries. In this paper, we present a simulation framework that we have implemented inside the MATLAB‐based diffusion MRI simulator SpinDoctor that efficiently computes the matrix formalism representation for realistic neurons using the finite element method. We show that the matrix formalism representation requires a few hundred eigenmodes to match the reference signal computed by solving the Bloch‐Torrey equation when the cell geometry originates from realistic neurons. As expected, the number of eigenmodes required to match the reference signal increases with smaller diffusion time and higher b‐values. We also convert the eigenvalues to a length scale and illustrate the link between the length scale and the oscillation frequency of the eigenmode in the cell geometry. We give the transformation that links the Laplace eigenfunctions to the eigenfunctions of the Bloch‐Torrey operator and compute the Bloch‐Torrey eigenfunctions and eigenvalues. This work is another step in bringing advanced mathematical tools and numerical method development to the simulation and modeling of diffusion MRI.

Practical computation of the diffusion MRI signal based on Laplace eigenfunctions: permeable interfaces

“…In our previous work, 37 we depended on the Partial Differential Equation Toolbox of MATLAB to solve the generalized eigenvalue problem in Equation (29), but this is no longer needed in our current implementation. The MATLAB syntax [P, lambda] = eigs(S + Q, M, Neig, ‘‘smallestreal’', ‘‘IsSymmetricDefinite’', true) directly computes the smallest eigenvalues and eigenfunction nodal coordinates using an iterative algorithm, given the mass, stiffness, and flux matrices M , S , and Q .…”

“…We note that in our previous work 37 our code contained a numerical implementation of J x , J y , and J z that had a slight error. We have since corrected this error and we now describe the correct implementation.…”

“…The calculation of the Laplace eigendecomposition in complicated geometries using Monte Carlo based simulations would be essentially impossible due to computational time and memory limitations. In a recent work 37 using the FEM, we implemented inside SpinDoctor the computation of the matrix formalism representation for biological cells subject to impermeable boundary conditions and showed its usefulness in understanding the diffusion MRI signal from realistic neuron geometries. In this work, we extend our simulation framework to include geometries that contain permeable cell membranes.…”

Fourier Representation of the Diffusion MRI Signal Using Layer Potentials

“…The Bloch Torrey PDEbased methods solve the Bloch-Torrey partial differential equation, which describes the evolution of the complex transverse water proton magnetization under the influence of diffusion-encoding magnetic field gradient pulses. The predominant numerical methods to solve this PDE include the finite difference method [14], the finite element method [17,21,2], and the Matrix Formalism method [6,7,18]. Although some recent works in the diffusion MRI community [28,34] still utilize Monte-Carlo simulations, the Bloch-Torrey PDE-based methods have recently demonstrated their potential, including in high-performance computing settings [23,24,26,17,18] and in manifolds settings for thin-layer and thin-tube geometries [25].…”

“…The predominant numerical methods to solve this PDE include the finite difference method [14], the finite element method [17,21,2], and the Matrix Formalism method [6,7,18]. Although some recent works in the diffusion MRI community [28,34] still utilize Monte-Carlo simulations, the Bloch-Torrey PDE-based methods have recently demonstrated their potential, including in high-performance computing settings [23,24,26,17,18] and in manifolds settings for thin-layer and thin-tube geometries [25]. In addition to numerical efficiency, some Bloch-Torrey PDE-based methods allow for a better understanding of the diffusion mechanism.…”