Index of /rozprawy2/11451

Pełen tekst

(1)Faculty of Physics and Applied Computer Science. Doctoral thesis Karol Borkowski. Analysis and correction of errors in diusion tensor imaging due to gradient inhomogeneity. Supervisor: dr hab. Artur T. Krzy»ak Faculty of Geology, Geophysics and Environmental Protection, AGH, Cracow. Cracow, November 2018.

(2) Declaration of the author of this dissertation:. Aware of legal responsibility for making untrue statements I hereby declare that I have written this dissertation myself and all the contents of the dissertation have been obtained by legal means. data, podpis autora. Declaration of the thesis Supervisor:. This dissertation is ready to be reviewed. data, podpis promotora rozprawy. 2.

(3) F irst of all, I would like to acknowledge Alexandra. Elbakyan for her indisputable contribution to remove all bariers in the way of science. I would like to sincerely thank the following persons who supported my researches: dr hab. Artur Krzy»ak for giving me the oportunity to realize my researches under his supervision, dr Krzysztof Kªodowski for his invaluable help, especially at the beginning of my PhD stidies, prof. Henryk Figiel for the successfull cooperation during the rst two years of my PhD studies, dr Bogdan Figura for his valuable suggestions to clarify the mathematical formulation of the problem. I would also like to thank for the nancial support from: Marian Smoluchowski Cracow Scientic Consortium, The National Centre of Research and Development within the following projects: PBS2/A2/16/2013, STRATEGMED2/265761/10/NCBR/2015.. 3.

(4) Acronyms BSD-DTI DTI. b-matrix spatial distribution in diffusion tensor imaging. diffusion tensor imaging. DWI diffusion weighted imaging EPI. echo planar imaging. FA. fractional anisotropy. FID. free induction decay. FOV field of view GLR gradient-linearity ratio HARDI. high angular resolution diffusion imaging. MRI magnetic resonance imaging NMR nuclear magnetic resonance PGSE pulsed gradient spin echo RF. radio frequency. ROI. region of interest. sBSD-DTI SD. simplified b-matrix spatial distribution in diffusion tensor imaging. standard deviation. S-DTI standard diffusion tensor imaging SE. spin echo. SNR signal-to-noise ratio uBSD-DTI. uniform b-matrix spatial distribution in diffusion tensor imaging. 4.

(5) Contents Acronyms ................................................................................................................................... 4 Contents ...................................................................................................................................... 5 Streszczenie ................................................................................................................................ 7 Abstract ...................................................................................................................................... 8 I.. Introduction ........................................................................................................................ 9 I.1. Aim and motivation ......................................................................................................... 9 I.2. Organization of the thesis ................................................................................................ 9. II.. Background ...................................................................................................................... 10 II.1. The physics of nuclear magnetic resonance ................................................................. 10 Nuclear Magnetic Moment............................................................................................... 10 Magnetization Vector ....................................................................................................... 11 Time Evolution ................................................................................................................. 11 Sample Excitation ............................................................................................................ 13 Relaxation......................................................................................................................... 14 The link between the quantum spin and classical magnetization vector .......................... 15 II.2. The principles of magnetic resonance imaging ............................................................ 17 Echo Formation ................................................................................................................ 18 Spatial Encoding .............................................................................................................. 18 II.3. Diffusion magnetic resonance imaging ........................................................................ 21 The b-matrix ..................................................................................................................... 23 Applications of DTI ......................................................................................................... 23 II.4. High angular resolution diffusion-weighted imaging ................................................... 24 II.5. Diffusion gradient inhomogeneity ................................................................................ 27 II.5. B-matrix spatial distribution in DTI ............................................................................. 27 Determination of the diffusion properties of the phantom ............................................... 28. III.. Contribution .................................................................................................................. 30. III.1. The generalized Stejskal-Tanner equation for non-uniform diffusion gradients ........ 30 Theory .............................................................................................................................. 30 Discussion ........................................................................................................................ 33 Conclusion ........................................................................................................................ 34 III.2. Theoretical analysis of the phantom rotation in BSD-DTI ......................................... 35 III.3. Computer simulations of BSD-DTI versus standard DTI – overview ........................ 37 The phantoms ................................................................................................................... 37 5.

(6) The spatial pattern functions ............................................................................................ 38 III.4. Theoretical validation of the BSD-DTI and uBSD-DTI approaches .......................... 39 Materials and methods ..................................................................................................... 39 Results .............................................................................................................................. 41 Discussion ........................................................................................................................ 44 Conclusion ........................................................................................................................ 45 III.5. The simplified BSD-DTI calibration using three positions of the phantom ............... 46 Theory .............................................................................................................................. 46 Materials and methods ..................................................................................................... 46 Results .............................................................................................................................. 49 Discussion ........................................................................................................................ 50 Conclusion ........................................................................................................................ 51 III.6. Analysis and correction of errors in DTI-based tractography due to diffusion gradient inhomogeneity ...................................................................................................................... 52 Materials and methods ..................................................................................................... 52 Results .............................................................................................................................. 54 Discussion ........................................................................................................................ 59 Conclusion ........................................................................................................................ 60 III.7. Analysis and correction of the systematic the eigenvalues shift due to diffusion gradient inhomogeneity ...................................................................................................................... 61 Materials and methods ..................................................................................................... 61 Results .............................................................................................................................. 62 Discussion ........................................................................................................................ 74 Conclusion ........................................................................................................................ 76 References ................................................................................................................................ 77. 6.

(7) Streszczenie Obrazowanie tensora dyfuzji (ang. diffusion tensor imaging lub DTI) jest unikalną techniką obrazowanie magnetyczno-rezonansowego (ang. magnetic resonance imaging lub MRI), która pozwala na nieinwazyjne badanie mikrostruktury żywych tkanek. Dostarcza ona nowe rodzaje kontrastu, jak np. anizotropia frakcyjna, które nie są dostępne przy użyciu żadnej innej metody obrazowania. Ponadto, technika DTI pozwala na rekonstrukcję włókien w procesie traktografii, który uważany jest za jedno z najważniejszych osiągnięć ostatnich dekad w dziedzinie neuroobrazowania. Dokładność pomiaru DTI jest ograniczona przez jednorodność zastosowanych gradientów pola magnetycznego. System cewek gradientowych zaprojektowany jest tak, aby wytwarzać pole magnetyczne o wartości zmieniającej się liniowo wzdłuż danego kierunku. Powyższe założenie nie jest jednak spełnione w przypadku praktycznie dowolnego skanera MRI, co prowadzi to błędnego oszacowania tensora dyfuzji. Celem niniejszej rozprawy było zbadanie wpływu niejednorodności gradientów na obrazowanie tensora dyfuzji oraz rozwój niedawno zaproponowanej metody korekcji polegającej na wyznaczeniu rozkładu przestrzennego macierzy b (ang. b-matrix spatial distribution in DTI lub BSD-DTI). Na początku wyprowadzona została matematyczna zależność pomiędzy tensorem dyfuzji a intensywnością sygnału ważonego dyfuzyjnie dla ogólnego przypadku niejednorodnych gradientów. Następnie przeprowadzone zostały symulacje komputerowe eksperymentów DTI zakładając różne rozkłady przestrzenne pola gradientowego. Tensor dyfuzji został obliczony przy założeniu jednorodnych gradientów (metoda standardowa) oraz stosując technikę BSD-DTI w trzech wariantach. Ważniejsze wyniki symulacji zostały potwierdzone eksperymentalnie przy użyciu skanera MRI o indukcji pola magnetycznego wynoszącej 3 tesle. W części teoretycznej udowodniono, że w przypadku echa spinowego, równanie Stejskala-Tannera jest prawidłowe niezależnie od założenia o jednorodności gradientów. W pozostałych przypadkach posiada ono również część urojoną, która informuje o przesunięcia fazy sygnału ważonego dyfuzyjnie. Otrzymane wyniki stanowią teoretyczną podstawę dla metody BSD-DTI oraz metody korekcji opartej na zastosowaniu tzw. tensora cewki (ang. coil tensor). Na podstawie wyników eksperymentów oraz symulacji ustalono, że niejednorodność gradientów dyfuzyjnych w obrazowaniu DTI prowadzi do trzech rodzajów błędów: obniżenia dokładności i precyzji pomiaru, rotacji zmierzonej elipsoidy dyfuzji względem jej rzeczywistej orientacji oraz przesunięcia wartości własnych uśrednionych w obrębie jednorodnego obszaru. Wielkość niektórych z powyższych błędów jest zależna od orientacji włókien wewnątrz skanera. Potencjał metody BSD-DTI w korekcji błędów wynikających z niejednorodności gradientów dyfuzyjnych został teoretycznie zweryfikowany w pełnej oraz dwóch uproszczonych wersjach. Pierwsza z nich (uniform BSD-DTI) polega na założeniu, że własności dyfuzyjne zastosowanego fantomu są jednorodne w obrębie jego całej objętości. Druga uproszczona wersja pozwala uprościć proces kalibracji kosztem precyzji poprzez pomiar fantomu w trzech pozycjach zamiast sześciu. Pełna kalibracja BSD-DTI pozwala osiągnąć najwyższą jakość pomiaru bez względu na stopień niejednorodności gradientów. Z drugiej strony, jest ona najbardziej czasochłonna i najtrudniejsza do przeprowadzenia, więc zastosowanie jednej z uproszczonych wersji może stanowić korzystny kompromis.. 7.

(8) Abstract Diffusion tensor imaging (DTI) is a unique magnetic resonance imaging (MRI) modality that allows to examine the microstructure of tissues in-vivo and noninvasively. It delivers new kinds of contrast, like fractional anisotropy, which are unattainable with any other imaging technique. Moreover, DTI allows to infer the course of the fiber tracts in the tractography process, which is believed to be one of the most notable advances in the field of neuroimaging in recent decades. The reliability of DTI is limited by the homogeneity of the applied magnetic field gradient. The gradient system is meant to produce a magnetic field that varies along a given direction linearly with the distance from magnet isocenter. Nevertheless, this assumption is not fulfilled in the case of virtually any MRI scanner, what leads to the erroneous estimation of the diffusion tensor and potentially, to a misdiagnosis. The aim of this work was to identify the implications of the gradient inhomogeneity in DTI and to develop a recently proposed method of calibration, called b-matrix spatial distribution in DTI (BSD-DTI). To address this problem, first the mathematical relation between the diffusion weighted MRI signal and the diffusion tensor was derived for the general case of the inhomogeneous gradient. Then, on that basis, the computer simulations of DTI experiments were performed assuming various spatial distributions of the gradient field and various levels of the distortion. Thereafter, the diffusion tensors were calculated and compared with the assumed properties of the imaged phantoms. To rely not only on simulated data, the major findings were confirmed by means of an MRI experiment using a 3T scanner. In the theoretical part, it has been proven that in the case of spin echo-based sequences, the Stejskal-Tanner equation is still valid, even if the gradients are non-uniform. In other cases, the inhomogeneity introduces a phase shift of the diffusion weighted signal. The obtained results are the theoretical basis for the BSD-DTI method as well as for the correction method based on so-called coil tensor. On the basis of simulations and experiments, it has been found that the diffusion gradient inhomogeneity causes three types of errors in DTI. Namely, the decrease of the precision and accuracy, the deviation between the real and measured orientation of diffusion ellipsoid and the systematic shift of the diffusion tensor eigenvalues averaged over a homogeneous region and overestimation of anisotropy. Interestingly, some effects are dependent on the orientation of fibers inside the scanner. The potential of the BSD-DTI technique in the correction of the above errors was discussed and theoretically validated for its complete and two simplified versions. The former simplification, called uniform BSD-DTI, consists in assuming a high quality of the applied phantom. The latter one simplifies the calibration procedure at the expense of accuracy. The complete approach ensures the highest fidelity of DTI measurement regardless of the severity of the gradient distortion; however, it is also the most complex and time-consuming. Therefore, in some cases, the simplified variants can be a reasonable trade-off.. 8.

(9) I.. Introduction. I.1. Aim and motivation The spatial encoding as well as the diffusion weighting in magnetic resonance imaging is achieved by means of the magnetic field gradients. The resultant spatially variable magnetic field is assumed to be linear, however, it is not true because of the concomitant field, imperfect construction of the gradient coils, magnetic susceptibility interfaces inside the sample and other sequence-specific factors, like induced eddy currents. The gradient inhomogeneity causes two types of artifacts in DTI: geometrical distortions and erroneous estimation of the diffusion properties of the sample. Although, the former kind of artifacts is well described and routinely corrected, not as much has been reported about the latter one. Moreover, the diffusion tensor is commonly estimated using a relation, that has been derived under the assumption that the gradients are uniform throughout the sample, so technically it is not valid. The aim of this work is to establish a general mathematical relation between the MRI signal and the diffusion tensor, analyze the effect of the gradient nonlinearity in DTI and to theoretically validate and further develop the b-matrix spatial distribution technique – a method of great potential for correction of the errors arisen due to the gradient inhomogeneity.. I.2. Organization of the thesis Just like Gaul, the thesis is divided into three parts. The first one is the introduction. In the second part, the NMR phenomenon and magnetic resonance imaging is described, including diffusion tensor imaging. Thereafter, the causes and methods of correction of the diffusion gradient inhomogeneity are discussed with particular emphasis on the b-matrix distribution technique. The third part contains the author’s contribution. Its first two chapters are theoretical. In the former one, the mathematical relation between the diffusion weighted MRI signal and the diffusion tensor is derived for the general case of the inhomogeneous gradient. In the latter one, the rotation of the phantom inside the scanner during the BSD-DTI calibration is analyzed. In the subsequent chapters, the influence of the diffusion gradient inhomogeneity on the DTI measurements is analyzed by means of computer simulations. In each case, the tensor is calculated assuming a spatially constant b-matrix (standard DTI) and using the BSD-DTI technique. The complete and two simplified versions of BSD-DTI are considered. The results of the preliminary MRI experiments have been also presented to support the most important findings of the simulations.. 9.

(10) II.. Background. II.1. The physics of nuclear magnetic resonance Nuclear magnetic resonance (NMR) can be described either classically or quantum-mechanically. Both approaches are equivalent, what was emphatically indicated by the fact, that it was explained independently and almost simultaneously from the classical point of view by Bloch et. al. [1] and from the quantum point of view by Purcell et. al. [2]. Despite the fact, that this phenomenon is based purely on the quantum-mechanical physics of nuclear spins, the measured property is the macroscopic net magnetization yielded by entangled nuclear magnetic moments, and can be treated as a classical vector. In a similar way, many phenomena involving atoms are not regarded as quantum, although the atom formation cannot be explained classically. Nuclear Magnetic Moment The most fundamental fact underlying the nuclear magnetic resonance phenomenon is that nucleons possess magnetic dipole moment 𝜇⃗. It is directly proportional to their intrinsic ⃗⃗ (spin) and pointed in the same direction. Therefore, the relation between angular momentum 𝐾 them can be written down as: ⃗⃗ , 𝜇⃗ = 𝛾𝐾. (II.1.1). where 𝛾 = 267.513 ∙ 106. 𝑟𝑎𝑑 𝑠∙𝑇. is a coefficient of proportionality known as the gyromagnetic. ratio. Due to the wave nature of matter, the values of energy, momentum and angular momentum are quantized. Considering this fact, the magnitude of the nuclear magnetic moment can be also express as: 𝜇 = 𝛾ℏ𝐼. (II.1.2) ℎ. where ℏ = 2𝜋 = 1.054 ∙ 10−34 𝐽 ∙ 𝑠 is the angular momentum quantum and I is the nuclear angular momentum quantum number, which takes integer or half-integer values. 1. The value of I for nucleons is 2. Spins of even numbers of nucleons of the same kind (protons or neutrons) in nuclei cancel out, therefore only the elements with odd numbers of protons or neutrons can be sensed in the NMR experiment. Since the only nucleus used in diffusion measurements is proton – the nucleus of hydrogen, further considerations will be limited only to this case.. 10.

(11) Magnetization Vector ⃗⃗⃗ is the resultant of the moments of The nuclear magnetic moment per unit volume 𝑀 individual spins. ⃗⃗⃗ = 1 ∑𝑖 𝜇⃗𝑖 . 𝑀 𝑉. (II.1.3). In this work it will be referred to as the magnetization. ⃗⃗ = (0 0 𝐵0 )𝑇 in The magnetization of the sample put into the external magnetic field 𝐵 thermal equilibrium is governed by Curie’s law: ⃗⃗⃗ = 𝑀. 𝜒 𝜇0. ⃗⃗, 𝐵. (II.1.4) 𝐻. where 𝜇0 = 12.57 ∙ 10−7 𝑚 is the magnetic permeability and 𝜒 is the nuclear magnetic susceptibility described by the formula: 𝜒=. 𝛾 2 ℏ 2 𝜇0 𝑁 3𝑘𝐵 𝑇. 𝐼(𝐼 + 1) =. 𝜇 2 𝜇0 𝑁 (𝐼+1) 3𝑘𝐵 𝑇. 𝐼. ,. (II.1.5) 𝐽. where N is the number of spins per unit volume and 𝑘𝐵 = 1.38 ∙ 10−23 𝐾 is the Boltzmann constant. In the considered case of protons, formula II.1.5 simplifies to: 𝜒=. 𝜇 2 𝜇0 𝑁. (II.1.6). 𝑘𝐵 𝑇. For a numerical example, let us take the protons in water at the human body temperature 𝐽. 1. T=310 K and B=3 T. In this case we have 𝜇 = 1.4 ∙ 10−26 𝑇 and 𝑁 = 6,65 ∙ 1027 𝑚3 . The ⃗⃗⃗ , obtained from eq. II.1.4 for the nuclear susceptibility 𝜒 = 3.047 ∙ 10−4 , is magnetization 𝑀 𝐴 ⃗⃗. equal to 1.149 ∙ 10−9 and is pointed in the direction of 𝐵 𝑀. Time Evolution In order to derive the formula describing the time evolution of the magnetization, let us begin with the equation: ⃗⃗⃗ = 𝛾𝐾 ⃗⃗ 𝑀. (II.1.7). ⃗⃗ refers to the resultant angular momentum of nuclei per unit where, in this context, vector 𝐾 ⃗⃗ is: volume. The second Newton’s law for rotation applied to 𝐾 𝑑 𝑑𝑡. ⃗⃗ = 𝑇 ⃗⃗, 𝐾. (II.1.8). ⃗⃗ is torque. On the other hand, the definition of the magnetic dipole moment is: where 𝑇. 11.

(12) ⃗⃗ = 𝑀 ⃗⃗⃗ × 𝐵 ⃗⃗ 𝑇. (II.1.9). Combining equations II.1.8, II.1.9 and the time derivative of eq. II.1.7 we obtain: 𝑑 𝑑𝑡. ⃗⃗⃗ = 𝛾𝑀 ⃗⃗⃗ × 𝐵 ⃗⃗. 𝑀. (II.1.10). To solve this equation, let us calculate the cross product: 𝑖 ⃗⃗⃗ × 𝐵 ⃗⃗ = |𝑀𝑥 𝑀 𝐵𝑥. 𝑗 𝑀𝑦 𝐵𝑦. 𝑘 𝑀𝑧 | = 𝐵𝑧. = 𝛾(𝑀𝑦 𝐵𝑧 − 𝑀𝑧 𝐵𝑦 )𝑖 + 𝛾(𝑀𝑧 𝐵𝑥 − 𝑀𝑥 𝐵𝑧 )𝑗 + 𝛾(𝑀𝑥 𝐵𝑦 − 𝑀𝑦 𝐵𝑥 )𝑘,. (II.1.11). where i, j and k are versors indicating the directions of the x, y and z axis, respectively. Considering the fact that the static magnetic field is applied only along the z-direction, the equation of motion becomes: 𝑑 ⃗⃗⃗ = 𝛾𝐵0 (𝑀𝑦 𝑖 − 𝑀𝑥 𝑗), 𝑀 (II.1.12) 𝑑𝑡. or in a more legible form: 𝑑 𝑑𝑡 𝑑 𝑑𝑡 𝑑 𝑑𝑡. 𝑀𝑥 = 𝛾𝐵0 𝑀𝑦. (II.1.13a). 𝑀𝑦 = 𝛾𝐵0 𝑀𝑥. (II.1.13b). 𝑀𝑧 = 0.. (II.1.13c). Eq. II.1.13a can be solved by substituting 𝑀𝑦 = 𝑀𝑥 ∙ 𝑡𝑎𝑛(𝛾𝐵𝑡), and separating the variables, what yields: 1 𝑀𝑥. 𝑑𝑀𝑥 = 𝛾𝐵0 ∙ 𝑡𝑎𝑛(𝛾𝐵𝑡)𝑑𝑡.. (II.1.14). Integration of eq. 1.14 brings: 𝑀𝑥 (𝑡) = −𝑀𝑥 (0) ∙ cos(𝛾𝐵0 𝑡),. (II.1.15). where 𝑀0 is the amplitude. The initial phase has been chosen to be equal to 0. Analogously, the solution of eq. 1.13b is: 𝑀𝑦 (𝑡) = −𝑀𝑦 (0) ∙ sin(𝛾𝐵0 𝑡). (II.1.16). Eq. II.1.16 could be conveniently written down using a complex number: 𝑀⊥ (𝑡) = 𝑀𝑥 (𝑡) + 𝑖𝑀𝑦 (𝑡) = 𝑀⊥ (0) ∙ 𝑒 −𝑖𝛾𝐵0 𝑡 = 𝑀⊥ (0) ∙ 𝑒 𝑖𝜑 ,. 12. (II.1.17).

(13) where 𝑀⊥ is the projection of the magnetization on the xy-plane (transverse plane). Its evolution induces a voltage in the receiver coil and is observed as the NMR signal. The solution of eq. II.1.13c is straightforward: 𝑀𝑧 (𝑡) = 𝑐𝑜𝑛𝑠𝑡.. (II.1.18). It can be noticed, that the magnetization is stationary only if it lays in parallel to the ⃗⃗. Otherwise, it rotates around 𝐵 ⃗⃗ with the angular velocity given static external magnetic field 𝐵 by: ⃗⃗, 𝜔 ⃗⃗ = −𝛾𝐵. (II.1.19). that is known as the Larmor precession. The word “precession” is used as an analogy to the change in the orientation of the rotational axis of a rotating body possessing the moment of inertia. The minus sign in eq. 1.18 determines the clockwise direction of the precession. Sample Excitation The magnetization can be manipulated using a radio frequencies (RF) electromagnetic ⃗⃗1 rotating in the transverse plane with angular velocity 𝜔 pulse 𝐵 ⃗⃗′. When the RF pulse is applied, the total magnetic field becomes: ⃗⃗𝑡𝑜𝑡 = 𝐵 ⃗⃗ + 𝐵 ⃗⃗1 𝐵. (II.1.20). or according to eq. II.1.19: ⃗⃗⃗⃗ ⃗⃗𝑡𝑜𝑡 = − 𝜔 ⃗⃗1 𝐵 +𝐵 𝛾. (II.1.21). It is convenient to describe the effect of the RF pulse in the frame of reference rotating with ⃗⃗1 is stationary. Then, the velocity of a point determined by the angular velocity 𝜔 ⃗⃗′, in which 𝐵 ⃗⃗ relative to the stationary frame is given by the relation: position vector 𝑉 ⃗⃗ 𝑑𝑉 𝑑𝑡. ⃗⃗ ′ 𝑑𝑉. ⃗⃗ , = ( 𝑑𝑡 ) + 𝜔 ⃗⃗′ × 𝑉. (II.1.22). ⃗⃗⃗ as the where the prime symbol ( ′ ) denotes the movement in the rotating frame. Considering 𝑀 ⃗⃗ and comparing eq. II.1.21 with eq. II.1.10, the equation of motion in the vector quantity 𝑉 rotating frame becomes: ⃗⃗⃗ ′ 𝑑𝑀. ⃗⃗⃗ × (𝐵 ⃗⃗𝑡𝑜𝑡 + ( 𝑑𝑡 ) = 𝛾𝑀. ⃗⃗⃗⃗′ 𝜔 𝛾. ⃗⃗⃗⃗ 𝜔. ⃗⃗⃗ × (− + 𝐵 ⃗⃗1 + ) = 𝛾𝑀 𝛾. ⃗⃗⃗⃗′ 𝜔 𝛾. ),. (II.1.23). or briefly: ⃗⃗⃗ ′. 𝑑𝑀 ⃗⃗⃗ × 𝐵 ⃗⃗𝑒𝑓𝑓 , ( 𝑑𝑡 ) = 𝛾𝑀. (II.1.24). 13.

(14) ⃗⃗𝑒𝑓𝑓 is the effective magnetic field in the rotating frame. where 𝐵 ⃗⃗𝑒𝑓𝑓 becomes According to eq. II.1.23, if 𝜔 ⃗⃗ is equal to 𝜔 ⃗⃗′, the effective magnetic field 𝐵 ⃗⃗1 and the Larmor precession occurs around its direction. This phenomenon is known equal to 𝐵 as the nuclear magnetic resonance. The application of the RF pulse of appropriate strength and duration allows to rotate the magnetization about any angle. In practice, it is usually rotated about 90° and 180° (𝜋/2 and 𝜋 rad). Relaxation Equations II.1.13 were derived the under assumption that the spins do not interact either with the environment or with each other. In reality, the resultant magnetic field “felt” by every single nuclear magnetic moment is affected by the surrounding spins, which slightly accelerate or decelerate their rate of precession. The environment of every spin changes dynamically due to the thermal motions of molecules and, as a result, the spins lose their coherence and the signal declines. This process is called the spin-spin or transverse relaxation. It is stochastic and irreversible. On the other hand, the spins interact also with their surroundings (the “lattice”), which is a thermal reservoir. The state of the lowest energy is achieved if the magnetization is given by eq. II.1.4. It is called the equilibrium magnetization and its magnitude is denoted as 𝑀0 . The energy of the spins can be increased by the RF pulse, as described in the previous paragraph. After the excitation, the magnetization comes back to the equilibrium state transferring its energy into the lattice in the form of heat. This process is called the spin-lattice or longitudinal relaxation. Taking into account the above processes, the time evolution of the magnetization can be described by phenomenological Bloch equations [1]: 𝑑 𝑑𝑡. 𝑀𝑥 = 𝛾𝐵0 𝑀𝑦 −. 𝑑. 𝑀 = 𝛾𝐵0 𝑀𝑥 − 𝑑𝑡 𝑦 𝑑. 𝑀𝑥. (II.1.25a). 𝑇2 𝑀𝑦. (II.1.25b). 𝑇2. 1. 𝑀 = 𝑇 (𝑀0 − 𝑀𝑧 ). 𝑑𝑡 𝑧. (II.1.25c). 1. T1 and T2 are time constants characterizing the rate of the spin-lattice and spin-spin relaxations, respectively. Equation II.1.25a and II.1.25b can be solved in the same way as equations II.1.13a and II.1.13b. Their complex solution in the static frame is: −𝑡. 𝑀⊥ (𝑡) = 𝑀⊥ (0) ∙ 𝑒 𝑖𝛾𝐵0 𝑡 ∙ 𝑒 𝑇2 .. (II.1.26). The solution in the rotating frame becomes: −𝑡. 𝑀⊥ (𝑡) = 𝑀⊥ (0) ∙ 𝑒 𝑖𝜑 ∙ 𝑒 𝑇2 .. (II.1.27). 14.

(15) Equation II.1.25c can be solved using an integrating factor or simply by putting the magnetization and time variables on the opposite sides of the equation and integrate. The result is: −𝑡. −𝑡. 𝑀𝑧 (𝑡) = 𝑀𝑧 (0) ∙ 𝑒 𝑇1 + 𝑀𝑜 ∙ (1 − 𝑒 𝑇1 ).. (II.1.28). The rotating magnetization vector surrounded by the receiving coil induces current according to Faraday’s law. In consequence of the relaxation, the induced signal is an exponentially decaying sinusoid known as the free induction decay (FID). The link between the quantum spin and classical magnetization vector A common misconception pertaining NMR is the image of nuclear magnetic moments occupying two energetic levels associated with two possible spins. This phenomenon is known as the Zeeman effect. According to the Boltzmann distribution, the occupation of the state of lover energy is slightly more probable, what yields the magnetization vector. In this model, the nuclear magnetic resonance is perceived as the absorption of quanta of the electromagnetic field ⃗⃗1, which energies are: 𝐵 𝐸 = 𝜔ℏ,. (II.1.29). This energy corresponds to the difference between the two energy states (eigenstates). As it was shown in the Stern-Gerlach experiment [3], the measured energy of a particular spin is indeed quantized; however, the measured value is the macroscopic magnetization. In such a case, the spins inside the sample are entangled and the magnetic moment of a particular nuclei can be any superposition of the eigenstates. Then, the state vector |𝜓⟩ is given by: |𝜓⟩ = 𝑎| ↑⟩ +b| ↓⟩,. (II.1.30). where 𝑎 and b can be any complex numbers which satisfy the equation |𝑎|2 + |𝑏|2 = 1. The arrows represent the eigenstates associated with two possible values of spin. Therefore, a more appropriate image would be a “hedgehog” of magnetic moments (Fig. 1), which are more densely distributed on the upper side of the reference frame, yielding the net magnetization. The spins are rotated during the RF pulse by purely magnetic interactions without affecting their relative orientations.. 15.

(16) Figure 1. The depiction of the magnetization vector as the superposition of entangled nuclear magnetic moment in the absence a) and presence b) of the external static magnetic field. The illustrations were taken from [4] by permission of Csaba Szántay.. 16.

(17) II.2. The principles of magnetic resonance imaging “Twinkle, twinkle, T2*; How I wonder what you are! XY signal soon decays; Why do the spins go out of phase? Twinkle, twinkle, T2*; Something pulls those spins apart. Spin-spin crosstalk sets T2, But by then T2* is through. A brief duration here is sealed By an inhomogeneous field. Twinkle, twinkle, T2*; Now I know just what you are!” Greg Crowther To obtain an image rather than a collective signal from the whole sample, the magnetic field gradient 𝐺⃗ is applied, such that the z-component of the magnetic field is spatially variable. The gradient is defined by the relation: 𝐵𝑧 (𝑟⃗, 𝑡) = (𝐵0 + 𝑟⃗ ∙ 𝐺⃗ (𝑡)). (II.2.1). where 𝑟⃗ = [𝑥 𝑦 𝑧] and symbol (∙) denotes the inner dot product. The gradient makes the Larmor frequency conditional on the spatial position r: 𝜔(𝑟⃗, 𝑡) = 𝛾(𝐵 + 𝑟⃗ ∙ 𝐺⃗ (𝑡)). (II.2.2). The very first MRI experiment, that was performed by Paul Lauterbur, consisted in registering the FID signal in the presence of a static magnetic field gradient [5]. After each acquisition, the object was rotated about an axis perpendicular to the gradient direction and finally, the image was reconstructed by means of a back-projection algorithm. Nowadays, the FID signal is rarely used in MRI, because it decays very rapidly with a time constant T2*. It happens because the magnetic field is inhomogeneous due to the applied gradient, imperfections in the magnet engineering and distortions caused by the magnetic susceptibility interfaces within the sample. The above factors accelerate the transverse relaxation, so T2* is usually much shorter than T2. Their relationship can be expressed by the following equation: 1 𝑇2∗. 1. = 𝑇 + 𝛾Δ𝐵𝑖𝑛ℎ ,. (II.2.3). 2. where Δ𝐵𝑖𝑛ℎ is the magnetic field inhomogeneity in a given voxel.. 17.

(18) Echo Formation Unlike in the case of spin-spin relaxation, the signal decay caused by the magnetic field inhomogeneities can be reverser by applying two consecutive RF pulses [6]. The first pulse rotates the magnetization by a certain angle (usually 90°). It is called the excitation pulse or 𝜋/2-pulse. The second one is applied after a time period 𝜏 and rotates the magnetization by 180°, what yields a so-called spin echo (SE). It is called a refocusing pulse or 𝜋-pulse. The signal intensity in the SE-based sequences is dependent on the true T2 relaxation time. Alternatively, the echo could be obtained without applying the refocusing pulse if the polarity of the local gradient was reversed. A similar effect can be achieved by applying two consecutive magnetic field gradients G>> Δ𝐵𝑖𝑛ℎ of the opposite directions, which yield cocalled gradient recalled echo (GE or GRE). In the GE-based sequences the relaxation processes due to the B0 inhomogeneities are not reversed, so the signal intensity is dependent on T2* relaxation time, not T2. According to equations II.1.27 and II.1.28, the signal intensity in the spin echo experiment is given by the relation: 𝑆(𝑇𝐸, 𝑇𝑅) = 𝜌 ∙ (1 − 𝑒. −𝑇𝑅 𝑇1. −𝑇𝐸. ) ∙ 𝑒 𝑇2 ,. (II.2.4). where 𝜌 is a constant proportional to the proton (spin) density, TE=2𝜏 (echo time) is the interval between the excitation RF pulse and the maximum intensity of the echo. TR (repetition time) is the interval between subsequent excitation pulses. In the case of the GE experiment, T2 in eq. 2.4. is substituted with T2*. Spatial Encoding Most of today's MRI techniques are based on selective excitation and the Fourier transform rather than the back-projection. A set of three gradient coils generates the spatially variable magnetic field in three orthogonal directions, that define the x, y, and z axis of the cartesian laboratory coordinate frame. The first step to obtain a three-dimensional image is to selectively excite the spins in a chosen slice by applying the RF pulse simultaneously with the so-called slice-select gradient Gz. The gradient introduces a position-dependent spread in the Larmor frequency: Δ𝜔 = 𝛾 ∙ Δz∙Gz. (II.2.5). By applying an RF pulse of the same bandwidth, the spins within a slice of thickness Δz are excited. The slice has a rectangular excitation profile, where the all spins inside are flipped by sin(𝑥). the same angle and all spins outside are unexcited, if the RF pulse is a sinc function (. 𝑥. ) in. the time domain. The bandwidth of the RF pulse can be manipulated by adjusting its duration, according to the uncertainty principle for Fourier transforms. Another gradient Gx is applied continuously during the signal acquisition to spatially encode the spins by changing their frequencies along its direction. Therefore, it is called the 18.

(19) frequency-encode gradient. According to eq. II.2.2, it causes a spatially dependent phase accumulation given by: 𝑡. 𝑡. 𝜙(𝑥, 𝑡)=∫0 𝜔( 𝑥, 𝑡′)𝑑𝑡′ = −𝛾 ∫0 (𝐵 + 𝐺𝑥 (𝑡 ′ ) ∙ 𝑥) 𝑑𝑡.. (II.2.6). The 𝛾𝐵0 term can be omitted, what corresponds to transforming the description of this phenomenon into the rotating frame. Then, the accumulated phase becomes: 𝑡. 𝜙(𝑥, 𝑡)=−𝛾𝑥 ∫0 𝐺𝑥 (𝑡 ′ ) 𝑑𝑡′.. (II.2.7). Assuming that the repetition time is long enough for the magnetization to completely recover (e.g. TR>5T1), a discreate signal element 𝜕𝑆 during the application of the frequency-encode gradient becomes: 𝜕𝑆(𝑡) = 𝜌(𝑥) ∙ 𝑒. −𝑡 𝑇∗2. 𝑡. ∙ 𝑒 −𝑖𝛾𝑥 ∫0 𝐺𝑥 (𝑡. ′ )𝑑𝑡 ′. 𝑑𝑥.. (II.2.8). The last gradient Gy is applied prior to the signal acquisition and speeds up or slows down the spins along its direction. After it is switched off, the spins return to their initial frequency but keep the different phases. Following the application of both gradients, the signal from the element 𝜕𝑆 becomes: 𝜕𝑆(𝑡) = 𝜌(𝑥, 𝑦) ∙ 𝑒. −𝑡 𝑇∗2. 𝑡. ∙ 𝑒 −𝑖𝛾(𝑥 ∫0 𝐺𝑥. 𝑡. (𝑡 ′ )𝑑𝑡 ′ +𝑦 ∫0 𝐺𝑦 (𝑡 ′ )𝑑𝑡 ′ ). 𝑑𝑥𝑑𝑦.. (II.2.9). Defining quantities kx and ky such that: 1. 𝑡. 1. 𝑡. 𝑘𝑥 (𝑡) = 2𝜋 𝛾 ∫0 𝐺𝑥 (𝑡 ′ )𝑑𝑡 ′. (II.2.8a). 𝑘𝑦 (𝑡) = 2𝜋 𝛾 ∫0 𝐺𝑦 (𝑡 ′ )𝑑𝑡 ′ ,. (II.2.8b). equation II.2.9 can be simplified to: −𝑡 ∗. 𝜕𝑆(𝑘𝑥 , 𝑘𝑦 ) = 𝜌(𝑥) ∙ 𝑒 𝑇2 ∙ 𝑒 −𝑖2𝜋(𝑘𝑥 𝑥+𝑘𝑦 𝑦) 𝑑𝑥𝑑𝑦. (II.2.9). The resulting raw data matrix is often called the k-space. Since there is no restriction on the direction along which the gradient can be applied, it is reasonable to generalize eq. II.2.9 to three dimensions. It can be rewritten in the vector notation: −𝑡. ⃗⃗ ) = 𝜌(𝑟⃗) ∙ 𝑒 𝑇∗2 ∙ 𝑒 −𝑖𝑘⃗⃗∙𝑟⃗ 𝑑𝑟⃗, 𝜕𝑆(𝑘. (II.2.10). ⃗⃗ (𝑡) is the k-space vector defined by: where 𝑘 ⃗⃗ (𝑡) = 1 𝛾 ∫𝑡 𝐺⃗ (𝑡 ′ )𝑑𝑡 ′ . 𝑘 0 2𝜋. (II.2.11) 19.

(20) Omitting the relaxation term, the signal integrated over the entire volume is: ⃗⃗ ) = ∭ 𝜌(𝑟⃗) ∙ 𝑒 −𝑖𝑘⃗⃗∙𝑟⃗ 𝑑𝑟⃗ 𝑆(𝑘. (II.2.12). It can be noticed that the measured signal is the Fourier transform of the spin density. Therefore, if the whole k-space is collected, the spatial distribution of the spin density (image) can be obtained by the inverse Fourier transform: ⃗⃗ ) ∙ 𝑒 −𝑖𝑘⃗⃗∙𝑟⃗ 𝑑𝑘 ⃗⃗ 𝜌(𝑟⃗) = ∭ 𝑆(𝑘. (II.2.13). There are plenty of imaging sequences for sampling the k-space. The most routinely used is echo-planar imaging (EPI) [7].. 20.

(21) II.3. Diffusion magnetic resonance imaging To get the NMR experiment sensitive to diffusion, two additional gradient pulses are applied on both sides of the refocusing pulse. Such a modified spin-echo sequence is known as pulsed gradient spin echo (PGSE) [8]. If the spins have drifted after the first diffusion gradient into a region where the effective magnetic field, and thus the Larmor frequency, is different, the phase compensation due to the second gradient is not exact. The resultant signal attenuation can be evaluated by modifying the Bloch’s equations, as it was done by H. Torrey [9]: 𝑑 𝑑𝑡. 𝑀𝑥 = 𝛾𝐵𝑀𝑦 −. 𝑑. 𝑀 = 𝛾𝐵𝑀𝑥 − 𝑑𝑡 𝑦 𝑑 𝑑𝑡. 𝑀𝑥 𝑇2 𝑀𝑦 𝑇2. + ∇ ∙ 𝐷∇(𝑀𝑥 − 𝑀𝑥0 ). (II.3.1a). + ∇ ∙ 𝐷∇(𝑀𝑦 − 𝑀𝑦0 ). (II.3.1b). 1. 𝑀𝑧 = 𝑇 (𝑀0 − 𝑀𝑧 ) + ∇ ∙ 𝐷∇(𝑀𝑧 − 𝑀𝑧0 ),. (II.3.1c). 1. where D is the diffusion coefficient. The diffusion terms in eqs. II.3.1 stem from Fick’s law of diffusion [10]. To describe the transverse magnetization, eq. II.3.1b is multiplied by 𝑖 = √−1 and added to eq. 3.1a. Then, using eq. II.1.17 and omitting the relaxation terms, we get: 𝜕𝑀⊥. = −𝑖𝛾𝐵𝑀⊥ + ∇ ∙ 𝐷∇𝑀⊥. 𝜕𝑡. (II.3.2). Equation 3.2 was the starting point for deriving the relation between the signal attenuation due to the diffusion and the diffusion coefficient, which is known as the Stejskal-Tanner equation [11]: 𝐴(𝑡=2𝜏). 𝑙𝑛 (. 𝐴(0). )=-bD,. (II.3.3). where A(t) is the signal amplitude proportional to 𝑀⊥ and the b-value is given by: 1. 𝑏 = 𝛾 2 𝛿 2 (Δ − 3 𝛿)𝑔2 ,. (II.3.4). where 𝛿 is the duration of the diffusion gradient pulse, 𝑔 is its maximal strength and Δ is the interval between two consecutive pulses. To quantify the diffusion in the anisotropic media, a tensor is needed rather than a scalar ̂ is a 3 x 3 second order matrix, which diagonal elements value [12]. The diffusion tensor 𝐷 correspond to the diffusion coefficients along the axes of the laboratory coordinate frame, while the off-diagonal elements determine the rotation of the eigensystem with respect to the laboratory frame. Its explicit form is: 𝐷𝑥𝑥 ̂ = [𝐷𝑦𝑥 𝐷 𝐷𝑧𝑥. 𝐷𝑥𝑦 𝐷𝑦𝑦 𝐷𝑧𝑦. 𝐷𝑥𝑧 𝐷𝑦𝑧 ]. 𝐷𝑧𝑧. (II.3.5). 21.

(22) In the case of uncharged particles, the diffusion tensor is symmetrical. The Stejskal-Tanner equation for anisotropic diffusion is a quadratic form [8]: 𝐴(𝑡=2𝜏). 𝑙𝑛 (. 𝐴(0). ̂𝐷 ̂, ) = −𝑏:. (II.3.6). where the colon denotes the Frobenius product and 𝑏̂ is a 3 x 3 symmetrical matrix, which in the most general way can be expressed as: 2𝜏 𝑏̂ = 𝛾 2 ∫0 𝐹(𝑡)𝑇 𝐹(𝑡) 𝑑𝑡,. (II.3.7). 𝑡 where 𝐹(𝑡) = ∫0 𝐺𝑥 (𝑡 ′ )𝑑𝑡 ′ is the total field offset caused by the gradient. In this work, 𝑏̂ will. be referred to as the b-matrix. Equation II.3.6 applies also to the magnetic resonance imaging. In that case, it is assumed that the diffusion encoding is the same for the whole k-space and has a value as though, it was derived for its center. Such an approximation is justified because most of the diffusion contrast information is encoded by low frequencies. To estimate the diffusion tensor, the N diffusion weighted signals 𝑆𝑖 (at least six) with N different diffusion gradient directions have to be measured in addition to one reference signal 𝑆0 , to collect vector 𝑆⃗: 𝐴1 𝑆⃗ = [𝑙𝑛 (𝐴 ) 0. 𝑇. 𝐴. … 𝑙𝑛 ( 𝐴𝑁 )] .. (II.3.8). 0. Then, the set of Stejskal-Tanner equations for the acquired signals can be written as: 𝑆⃗ = −𝐵̂ 𝑑⃗ ,. (II.3.9). where 𝑑̂ is a column vector containing six independent elements of the diffusion tensor: 𝑑⃗ = [𝐷𝑥𝑥. 𝐷𝑦𝑦. 𝐷𝑧𝑧. 𝐷𝑥𝑦. 𝐷𝑥𝑧. 𝐷𝑦𝑧 ]𝑇. (II.3.10). and 𝐵̂ is a 6 x N matrix containing independents of the b-matrices for all gradient directions: 𝑏1𝑥𝑥 𝐵̂ = [ ⋮ 𝑏𝑁𝑥𝑥. 𝑏1𝑦𝑦 ⋮ 𝑏𝑁𝑦𝑦. 𝑏1𝑧𝑧 ⋮ 𝑏𝑁𝑧𝑧. 2𝑏1𝑥𝑦 ⋮ 2𝑏𝑁𝑥𝑦. 2𝑏1𝑥𝑧 ⋮ 2𝑏𝑁𝑥𝑧. 𝑇. 𝑏21𝑦𝑧 ⋮ ] . 2𝑏𝑁𝑦𝑧. (II.3.11). The formula for the elements of the diffusion tensor is obtained by rearranging eq. II.3.9: 𝑑⃗ = −𝐵̂ −1 𝑆⃗.. (II.3.12). If more than six directions of the gradient are applied, the tensor can be estimated by means of linear regression.. 22.

(23) The b-matrix Apart from experimental techniques like BSD-DTI, that are discussed in further chapters, there are two common ways of determining the b-matrix. The former one consists in integrating all the gradients present during the imaging sequence (eq. II.3.7), as it was described by Mattiello et al. [13]. This approach will be referred to as standard DTIA or S-DTIA. More popular is the latter one, which uses the b-value (eq. 3.4) together with a dyadic product of the gradient unit vector g: 𝑔𝑥 ̂ 𝑔 𝑏 = 𝑏 [ 𝑦 ] [𝑔𝑥 𝑔𝑧. 𝑔𝑥. 𝑔𝑥 ] = 𝑏𝑔𝑔𝑇 .. (II.3.13). This approach is based on the assumptions that the b-value is the same for each direction and that the off-diagonal elements of the b-matrix can be derived from the diagonal ones. It will be referred to as the standard DTIB or S-DTIB. If the method of delivering the b-matrix is not given explicitly, in this work, it should by assumed that the standard DTIB approach was applied. Applications of DTI Diffusion tensor imaging derives new kinds of contrast, which allows to distinguish tissues on the basis of their anisotropy [14]. It is expressed by a scalar factor called Fractional Anisotropy (FA) [15]: (𝜆1 −𝜆2 )2 +(𝜆2 −𝜆3 )2 +(𝜆3 −𝜆1 )2. 𝐹𝐴 = √. 2(𝜆1 2 +𝜆2 2 +𝜆3 2 ). .. (II.3.14). It allows the characterization of subtle changes associated with aging and various neurodegenerative diseases [16]. In tissues with a well-organized microstructure, such as nerves, muscles, ligaments or tendons, water diffusivity is greatest along the fibers direction [17]. Therefore, knowing the direction of the eigenvector associated with the largest eigenvalue it is possible to infer the course of the fiber tracts in a process called tractography [17]–[19]. The fiber tracking has multiple clinical applications, including the assessment of important white matter tracts during surgery or pathological changes due to ischemic stroke or degenerative diseases [21]. Tractography is believed to be one of the most notable advances in the field of neuroimaging in recent decades. A promising field of application of DTI is also cardiac imaging [22]. It delivers unique information about the microstructure of myocardium in-vivo and noninvasively, what is beyond the reach of any other imaging technique. Nevertheless, it is challenging because of the cardiac movement and its short T2 relaxation time.. 23.

(24) II.4. High angular resolution diffusion-weighted imaging Diffusion tensor imaging suffers from a number of limitations, like deviation from the Gaussian diffusion model or so-called crossing fibers. The second issue is a case in about 90% of white matter voxels [23]. The above limitations can be addressed in two ways: focusing on the diffusion process or using some model of the imaged microstructure. Q-space The former approach demands less assumptions, so it is often called model-free. In such a case, the estimation of the preferred directions of diffusion requires the detailed characterization of so-called q-space, which corresponds to the diffusion weighted signals acquired for various gradient directions and b-values. This approach differs from the PGSE method by assuming very strong and short diffusion gradients of duration 𝛿, so that the molecules displacement due to diffusion occurs only during the time interval ∆ between two consecutive diffusion gradient pulses. The q-space imaging can be described by considering the ensemble magnetization, which is the sum of the spins in a given voxel of volume V: 𝜓(𝐺, ∆) = ∫ 𝑒 𝑖2𝜋𝑞𝑅 𝑑𝑉,. (II.4.1). where G is the gradient amplitude and R is the displacement over time ∆ and q is defined as: 1. 𝑞 = 2𝜋 𝛾𝐺𝛿.. (II.4.2). Defining the probability density function (PDF) P(R, ∆), that informs what fraction of spins within a voxel displaced by R over the time ∆, the q-space signal S is: 𝑆(𝑞, ∆) = ∫ P(R, ∆)𝑒 𝑖2𝜋𝑞𝑅 𝑑𝑉.. (II.4.3). According to the above formula, a 3D Fourier transform of q-space data yields the probability density function. This method provides very rich diffusion information, which allows to resolve a complex structure of fibers. Since it is model-free, it does not violate the Gaussian diffusion model assumption. On the other hand, it is very time consuming. An exemplary imaging protocol with about 500 q-values takes about 30 min. Nevertheless, this time could be reduced using the compressed sensing technique. Since a wide range of the q-space spectrum is covered, this technique is often referred to as the diffusion spectrum imaging (DSI). 24.

(25) The fiber orientation distribution function (ODF) Various information can be extracted from the displacement distribution. The most commonly used is the fiber orientation distribution function (ODF), which allows to perform fiber tracking. It is a metric of how likely the fibers are oriented along a given direction. The probability that molecules move along direction u equals the integrate of PDF values along that direction, weighted by the square distance from the origin: 𝑂𝐷𝐹(𝑢) = ∫ 𝑅 2 P(R, ∆)dV.. (II.4.3). Q-ball A less time-consuming model-free approach is the Q-ball imaging. In this technique, the q-space is traversed only on a sphere. This means that multiple directions of the diffusion gradient are used, but each of the same b-value. The fiber ODF function is obtained using the Funk-Radon transformation [24]. Parametric models Although DSI provides the most complete description of diffusion, this information is usually reduced. For example, only 9 parameters are needed to describe three fiber orientations. Therefore, a most efficient would be model-based approaches, for which the minimal number of measurements equals the number of parameters. Multitensor model The most intuitive diffusion model, that allows to distinguish crossing fibers, is the multitensor model. It consists in fitting multiple ellipsoids (tensors) to the measured data, so only 12 DW measurements are required to resolve two crossing fibers. Nevertheless, this approach becomes unstable for a high number of tensors. Diffusion kurtosis imaging The average value of an arbitrary function F(R) is: < 𝐹(𝑅) >= ∫ 𝑃(𝑅, 𝑡)𝐹(𝑅)𝑑𝑉,. (II.4.4). where P(R,t) is the probability density function. The excess kurtosis K(t) in a given direction defined by vector n is: <(𝑛⋅𝑅)4 >. 𝐾(𝑡) = <(𝑛⋅𝑅)2 >2 − 3.. (II.4.5). 25.

(26) For isotropic Gaussian diffusion, the kurtosis is equal to zero. Therefore, it can be understood as a dimensionless metric of the deviation of a given function from the Gaussian distribution. If K(t) is negative, the given function P is more sharply peaked than a Gaussian. If it is positive, the function is less sharply peaked. The diffusion coefficient is given by: 1. 𝐷(𝑡) = 2𝑡 < (𝑛 ⋅ 𝑅)2 >.. (II.4.6). The relationship between the diffusion weighted signal and the diffusional kurtosis is: 𝑆. 1. 𝑙𝑛 (𝑆 ) = −𝑏 ∙ 𝐴𝐷𝐶 + 6 𝑏 2 𝐴𝐷𝐶 2 𝐾𝑎𝑝𝑝 ,. (II.4.7). 0. where 𝐾𝑎𝑝𝑝 is the apparent kurtosis. The diffusional kurtosis can be estimated by calculating a limit: 𝐾(∆) = lim 𝐾𝑎𝑝𝑝 (∆, 𝛿).. (II.4.8). 𝛿→∞. The diffusion kurtosis tensor has a rank of 4 and 81 components. Nevertheless, it is fully symmetric, so only 15 components are independent and therefore, at least 15 diffusion weighted measurements are required to determine it.. 26.

(27) II.5. Diffusion gradient inhomogeneity Although the magnetic field gradient is assumed to be uniform, is has been shown multiple times that this assumption is not fulfilled in the case of virtually any MRI scanner [25]– [29]. The reason is that, according to Maxwell’s equations [30], a linear variation of a single component of the magnetic field cannot exist in free space. As a consequence, the gradient field is distorted by a so-called concomitant field [31]. Apart from that, the gradients can be distorted by eddy currents induced in the RF coil and other conductive surfaces due to rapidly switched gradient pulses, B0 field inhomogeneities and imperfect construction of the gradient coils. The gradient inhomogeneity causes two types of artifacts in DTI: spatial image distortions and inaccurate estimation of the diffusion tensor. While the former kind of artifacts is well described and routinely corrected [24], [29], [30], there is no gold standard in the correction of the latter one. The diffusion gradient field can be mapped using an isotropic phantom [34], but it allows to correct only the gradient’s magnitude, not its direction. A more robust correction, that was proposed by R. Bammer et al., can be achieved by multiplying the gradient by a so-called gradient coil tensor 𝐿̂(𝑟⃗) [35]: 𝐿𝑥𝑥 (𝑟⃗) 𝐿𝑥𝑦 (𝑟⃗) 𝐿𝑥𝑧 (𝑟⃗) 𝐺⃗𝑎𝑐𝑡 (𝑟⃗) = [𝐿𝑦𝑥 (𝑟⃗) 𝐿𝑦𝑦 (𝑟⃗) 𝐿𝑦𝑧 (𝑟⃗)] 𝐺⃗𝑑𝑒𝑠 =𝐿̂(𝑟⃗)𝐺⃗𝑖𝑛𝑡 , 𝐿𝑧𝑥 (𝑟⃗) 𝐿𝑧𝑦 (𝑟⃗) 𝐿𝑧𝑧 (𝑟⃗). (II.5.1). where 𝐺⃗𝑎𝑐𝑡 (𝑟⃗) is the actual distribution of the gradient and 𝐺⃗𝑖𝑛𝑡 is the intended gradient. The components 𝐿̂𝑖𝑗 (𝑟⃗) are expressed by the formula: 𝐿̂(𝑟⃗) =. ⃗⃗))𝑖 𝜕(𝐵𝑧 (𝑟 𝜕𝑟𝑗. 𝐺𝑖. .. (II.5.2). They represent the real magnetic field gradient generated by the ith coil for a unit gradient intended in the jth direction. Since the magnetic field in eq. 5.2 has only one non-zero component, it can be treated as a scalar function. Bz can be expressed by means of spherical harmonics [26] and indicated experimentally by the magnetic field mapping [36] or on the basis of specification provided by the manufacturer. The downside of the experimental approach is that a different sequence is used for calibration and different for the imaging, so the sequencerelated factors, like eddy currents, are not considered. On the other hand, the specification of the gradient system is a commercially sensitive information. Moreover, the exact functioning of the gradient coils can slightly differ even between scanners of the same model.. II.5. B-matrix spatial distribution in DTI An alternative approach to the correction of the diffusion gradient inhomogeneity is the b-matrix spatial distribution in the DTI technique (BSD-DTI). It allows to indicate the exact form of the b-matrix for each voxel separately using an anisotropic phantom [29]. The real form. 27.

(28) of the b-matrix is mapped regardless of the sources of the distortion and without the need of any prior knowledge about the gradient parameters. This method was introduced in 2008 in the Polish patent application [37], where the full concept of using the b-matrix spatial distribution in DWI and DTI measurements was presented. In the following years, American [38], European [39] and Japanese [40] patents were also granted. The patents were granted in all the markets, where the PCT procedures were initiated [41]. In 2013 a matrix of isotropic and anisotropic phantoms was proposed to indicate the bmatrix spatial distribution. This approach, called the tensor spatial distribution in DTI (TSDDTI), was described in the Polish patent application [42]. To carry out the BSD-DTI calibration an anisotropic phantom is imaged using the same DTI protocol that is to be used in the proper measurement. Then, the phantom is rotated M-1 (M⩾6) times by given sets of Euler angles and imaged after each rotation to collect matrix 𝑆̂(𝑟⃗): 𝐴(𝑟⃗). 𝑆̂ =. …. 𝑙𝑛 ( 𝐴(𝑟⃗)1,𝑁 ). ⋮. ⋱. ⋮. 0,1. 𝐴(𝑟⃗)𝑀,1. 𝑙𝑛 (𝐴(𝑟⃗). [. 0,𝑀. ) ⋯. 𝑇. 𝐴(𝑟⃗). 𝑙𝑛 (𝐴(𝑟⃗)1,1 ). 0,1. 𝐴(𝑟⃗)𝑀,𝑁. 𝑙𝑛 ( 𝐴(𝑟⃗). 0,𝑀. ,. (II.5.3). ) ]. where N is the number of diffusion gradient directions and M is the number of the phantom positions. The rotation angles can be chosen freely, subject that the diffusion tensor of the phantom in the laboratory coordinate frame changes after rotation. In particular, it should not be rotated around its symmetry axis. Since the diffusion properties of the phantom are known in each of the M orientations, they can be collected in a matrix: 𝐷𝑥𝑥,1 ̂ 𝑑=[ ⋮ 𝐷𝑥𝑧,𝑀. 𝐷𝑦𝑦,1 ⋮ 𝐷𝑦𝑦,𝑀. 𝐷𝑧𝑧,1 ⋮ 𝐷𝑧𝑧,𝑀. 𝐷𝑥𝑦,1 ⋮ 𝐷𝑥𝑧,𝑀. 𝐷𝑥𝑧,1 ⋮ 𝐷𝑥𝑧,𝑀. 𝑇. 𝐷𝑦𝑧,1 ⋮ ] . 𝐷𝑦𝑧,𝑀. (II.5.4). Then, the set of the Stejskal-Tanner equations can be written down for the whole experiment in the form of one matrix equation: 𝑆̂(𝑟⃗) = −𝑑̂ 𝐵̂ (𝑟⃗). (II.5.5). Finally, the independent elements of the b-matrix for each gradient direction can be obtained using the formula: 𝐵̂ (𝑟⃗) = −𝑑̂ −1 𝑆̂(𝑟⃗),. (II.5.6). which is the essence of the BSD-DTI technique. Determination of the diffusion properties of the phantom To perform the BSD-DTI calibration the diffusion properties of the phantom must be already well-known. Nevertheless, it may be challenging to ascertain them, since the MRI28.

(29) based techniques are imprecise. A possible solution could be to perform the DTI measurement applying as many gradient directions as reasonably possible (e.g. 512) and appropriately long repetition time, to prevent the phantom from heating up. Moreover, the diffusion tensor of the phantom can be estimated applying the BSD-DTI calibration using an isotropic medium of known diffusion coefficient, such as water. For the anisotropic medium, equation II.3.6 simplifies to: D 0 𝐴 ̂ = −𝑏̂(𝑟⃗): [ 0 𝐷 𝑙𝑛 (𝐴 ) = −𝑏̂(𝑟⃗): 𝐷 0 0 0. 0 0 ] = −(𝑏𝑥𝑥 + 𝑏𝑦𝑦 + 𝑏𝑧𝑧 )D. D. (II.5.7). Hence, the b-matrix trace (the b-value in eq. II.3.4) for a given gradient direction can be calculated by: 𝑏′ =. −𝑙𝑛(. 𝐴 ) 𝐴0. 𝐷. ,. (II.5.8). where the prime symbol denotes the corrected value. Finally, knowing the elements bij of the intended b-matrix, the elements after correction b’ij can be obtained using the formula: 𝑏𝑖𝑗. 𝑏′ 𝑏. ,𝑖 = 𝑗. 𝑏′𝑖𝑗 = { . 𝑏′ ∙𝑏′ 𝑏𝑖𝑗 √ 𝑏𝑖𝑖 ∙𝑏 𝑗𝑗 , 𝑖 ≠ 𝑗 𝑖𝑖. (II.5.9). 𝑗𝑗. The above procedure should be applied to every diffusion weighted measurement separately.. 29.

(30) III. Contribution III.1. The generalized Stejskal-Tanner equation for non-uniform diffusion gradients Since the Stejskal-Tanner equation (eq. II.3.6) has been derived under the assumption that the gradient is homogeneous within the sample, formally it is not valid. Moreover, in the strategies to correct the diffusion gradient inhomogeneity, that were described in the previous chapter, the diffusion tensor is ultimately estimated using this equation. It implies an underlying assumption that the gradients are uniform, what is a direct contradiction. For that reason, in this chapter, the general relation between the NMR signal and the diffusion tensor is derived for the general case of the non-uniform diffusion gradient. The presented proof has been published in the Journal of Magnetic Resonance [43]. Theory The diffusion weighted signal can be modeled either as the diffusion of spin-bearing molecules, or as the diffusion of magnetization. The former approach is based on the cumulant expansion, while the latter one consists in solving the Bloch-Torrey equation (eq. II.3.2). The following discussion represents the latter approach. Replacing the diffusion coefficient in eq. II.3.2 with the diffusion tensor and introducing the notation used by Stejskal, where 𝑀⊥ is replaced by 𝜓(𝑟⃗, 𝑡), the Bloch-Torrey equation becomes: 𝜕𝜓 𝜕𝑡. ̂ ∇𝜓, = −𝑖𝛾𝐵𝜓 + ∇ ∘ 𝐷. (III.1.1). where (∘) is used to indicate the outer dot product and the gradient of a vector, when it is used with nabla. 𝜓 is given by: 𝜓 = 𝑀0 ∙ 𝑒 𝑖𝜑 ,. (III.1.2). but this time, the phase is assumed to be spatially variable, so that 𝜑(𝑟⃗, t): R4⟶R. The spinspin relaxation term has been omitted, as it was done by Stejskal and Tanner. It could be done without the loss of generality, since we consider only the signal attenuation due to diffusion, not the absolute signal intensities. The arguments of ψ and 𝜑 are left off for brevity. In the absence of the diffusion part the decreasing amplitude of the transverse magnetization vector can be counted by replacing the constant M0 with a function of time A(t). Then, eq. III.1.2 becomes: 𝜓 = 𝐴(𝑡) ∙ 𝑒 𝑖𝜑 .. (III.1.3). The time derivative of the above equation is:. 30.

(31) 𝜕𝜓 𝜕𝑡. =. 𝜕(𝐴(𝑡)∙𝑒 −𝑖𝜑 ) 𝜕𝑡. =. 𝜕𝐴(𝑡) 𝜕𝑡. 𝜕𝜑. 𝑒 𝑖𝜑 + 𝑖𝐴(𝑡) 𝜕𝑡 𝑒 𝑖𝜑 .. (III.1.4). Now, let us consider the diffusion part of eq. III.1.1, where the nabla operator is taken with respect to r: ̂ ∇𝜓 = ∇ ∘ 𝐷 ̂ ∇ [𝐴(𝑡) ∙ 𝑒 𝑖𝜑 ] = 𝐴(𝑡)[∇ ∘ 𝐷 ̂ ∇𝑒 𝑖𝜑 ] = 𝐴(𝑡)[∇ ∘ (𝐷 ̂ 𝑒 𝑖𝜑 ∙ 𝑖∇φ)]. ∇∘𝐷. (III.1.5). The chain rule for the gradient of the vector product is: ∇ ∘ (𝑓𝑤) = ∇f ∘ w + f(∇ ∘ w). (III.1.6). Using equation III.1.6, the diffusion part of eq. III.1.1 can be further transformed: ̂ 𝑒 𝑖𝜑 ∙ 𝑖∇φ)] = −𝑖𝐴(𝑡)[∇𝑒 𝑖𝜑 ∘ 𝐷 ̂ ∇𝜑 + 𝑒 𝑖𝜑 ∇ ∘ 𝐷 ̂ ∇𝜑], 𝐴(𝑡)[∇ ∘ (𝐷. (III.1.7). and finally: ̂ ∇𝜓 = −𝐴(𝑡) ∙ 𝑒 −𝑖𝜑 [∇𝜑 ∘ 𝐷 ̂ ∇𝜑 − 𝑖(∇ ∘ 𝐷 ̂ ∇𝜑)]. ∇∘𝐷. (III.1.8). Substituting equations III.1.4 and III.1.8 into eq. III.1.1 yields: 𝜕𝐴(𝑡) 𝜕𝑡. 𝜕𝜑. 𝜕𝜑. ̂ ∇𝜑 − 𝑖(∇ ∘ 𝐷 ̂ ∇𝜑)] 𝑒 𝑖𝜑 + 𝑖𝐴(𝑡) 𝜕𝑡 𝑒 𝑖𝜑 = 𝑖𝐴(𝑡) 𝜕𝑡 𝑒 𝑖𝜑 − 𝐴(𝑡) ∙ 𝑒 −𝑖𝜑 [∇𝜑 ∘ 𝐷. (III.1.9). Clearing terms algebraically and dividing both sides by ei , gives: 𝜕𝐴(𝑡) 𝜕𝑡. ̂ ∇𝜑 − 𝑖(∇ ∘ 𝐷 ̂ ∇𝜑)]. = −𝐴(𝑡)[∇𝜑 ∘ 𝐷. (III.1.10). In the general case of non-uniform gradients, the effective magnetic field in the rotating frame can be defined using a pattern function 𝑝⃗(𝑟⃗) instead of 𝑟⃗: 𝐵𝑧 (𝑟⃗, 𝑡) = (𝐵𝑧 (𝑟⃗))𝑥 + (𝐵𝑧 (𝑟⃗))𝑦 + (𝐵𝑧 (𝑟⃗))𝑧 = 𝑝⃗(𝑟⃗) ∙ 𝐺⃗ (𝑡) = 𝑝𝑥 𝐺𝑥 + 𝑝𝑦 𝐺𝑦 + 𝑝𝑧 𝐺𝑧 . (III.1.11) Such an approach consists in regarding the gradient as linear in a curvilinear coordinate system given by 𝑝⃗(𝑟⃗). A similar expression for Beff was employed e.g. by Zhu in his US Patent Application Publication [44]. (𝐵𝑧 (𝑟⃗))𝑖 indicates the magnetic field generated by the ith gradient coil. Assuming the spatial distribution of the gradient field determined by 𝑝(𝑟⃗), the phase of 𝜓 in the PGSE sequence is given by: 𝜑(𝑟⃗, t) = −𝑝⃗(𝑟⃗) ∙ 𝛾[𝐹(𝑡) − (𝜁 − 1)𝑓] = −𝑝⃗(𝑟⃗) ∙ ⃗⃗⃗⃗ 𝑘′(𝑡), 31. (III.1.12).

(32) where 𝑓 ≡ 𝐹(2𝜏) and 𝜁 = {. 1, 𝑡 < 𝜏 . −1, 𝑦 > 𝜏. Equation III.1.12 considers the fact that the phase is reversed by the refocusing pulse. In such an approach, ⃗⃗⃗⃗ 𝑘′(𝑡) can be regarded as the k-space vector multiplied by 2𝜋 at any time of the experiment. Employing eq. III.1.12, the gradient of the phase in eq. III.1.10 is: ⃗⃗⃗⃗(𝑡). ⃗⃗ ′(𝑡)) = −(∇ ⋄ 𝑝⃗(𝑟⃗))𝑘′ ∇𝜑(𝑟⃗, t) = −∇ (𝑝⃗(𝑟⃗) ∙ 𝑘. (III.1.13). The gradient of the pattern function 𝑝⃗(𝑟⃗) in the above equation has the form:. ∇ ⋄ 𝑝⃗(𝑟⃗) =. 𝜕𝑝𝑥. 𝜕𝑝𝑥. 𝜕𝑝𝑥. 𝜕𝑥 𝜕𝑝𝑦. 𝜕𝑦 𝜕𝑝𝑦. 𝜕𝑧 𝜕𝑝𝑦. 𝜕𝑥 𝜕𝑝𝑧. 𝜕𝑦 𝜕𝑝𝑧. 𝜕𝑧 𝜕𝑝𝑧. 𝜕𝑦. 𝜕𝑧. [ 𝜕𝑥. ≝ 𝐿̂(𝑟⃗).. (III.1.14). ]. It can be regarded as the Jacobian matrix for the coordinates change from Cartesian to those given by 𝑝⃗(𝑟⃗). Substitution of equations III.1.13 and III.1.14 into eq. III.1.10 yields: 𝜕𝐴(𝑡) 𝜕𝑡. ⃗⃗⃗⃗(𝑡) ∘ 𝐷 ⃗⃗⃗⃗(𝑡) − 𝑖 (∇ ∘ 𝐷 ⃗⃗⃗⃗(𝑡))]. ̂ 𝐿(𝑟⃗)𝑘′ ̂ 𝐿(𝑟⃗)𝑘′ = −𝐴(𝑡) [𝐿(𝑟⃗)𝑘′. (III.1.15). Equation 1.15 can be solved by separating variables and integration: 𝐴(𝑡) 𝜕𝐴(𝑡′). ∫𝐴(0). 𝐴(𝑡′). 𝑡. ⃗⃗⃗⃗′ (𝑡′) ∘ 𝐷 ⃗⃗⃗⃗′ (𝑡′) − 𝑖 (∇ ∘ 𝐷 ⃗⃗⃗⃗′ (𝑡′))] 𝜕𝑡′. ̂ 𝐿̂(𝑟⃗)𝑘 ̂ 𝐿̂(𝑟⃗)𝑘 = − ∫0 [𝐿̂(𝑟⃗)𝑘. (III.1.16). Substituting t=2τ the solution is: 2τ 𝐴(2τ) ⃗⃗⃗⃗′ (𝑡) ∘ 𝐷 ⃗⃗⃗⃗′ (𝑡) + 𝑖∇ ∘ 𝐷 ̂ 𝐿̂(𝑟⃗)𝑘 ̂ 𝐿̂(𝑟⃗) ∫2τ ⃗⃗⃗⃗ 𝑙𝑛 ( 𝐴(0) ) = − ∫0 𝐿̂(𝑟⃗)𝑘 𝑘′(𝑡) 𝜕𝑡. 0. (III.1.17). Equation III.1.17 is the most general form of the Stejskal-Tanner equation. In the case of the spin echo-based sequences, it can be further simplified by noticing that, according to eq. 1.12, 2τ the integral ∫0 ⃗⃗⃗⃗ 𝑘′(𝑡) 𝜕𝑡 is equal to zero. Therefore, the equation becomes: 2τ 𝐴(2τ) ⃗⃗⃗⃗′ (𝑡) ∘ 𝐷 ⃗⃗⃗⃗′ (𝑡) 𝜕𝑡 ̂ 𝐿̂(𝑟⃗)𝑘 𝑙𝑛 ( 𝐴(0) ) = − ∫0 𝐿̂(𝑟⃗)𝑘. or in matrix notation:. 32. (III.1.18).

(33) 2τ ′ 𝐴(2τ) ⃗⃗⃗⃗′ (𝑡) 𝜕𝑡. ̂ 𝐿̂(𝑟⃗)𝑘 𝑙𝑛 ( 𝐴(0) ) = − ∫0 ⃗⃗⃗⃗ 𝑘 (𝑡)𝑇 𝐿̂(𝑟⃗)𝑇 𝐷. Using eq. II.3.7 the final form of the generalized Stejskal-Tanner equation is: 𝐴(2τ). ̂. 𝑙𝑛 ( 𝐴(0) ) = −𝑏̂(𝑟⃗): 𝐷. (III.1.19). Knowing the theoretical b-matrix 𝑏̂𝑠𝑡𝑑 , which was calculated assuming a homogeneous gradient, its real form can be obtained for each voxel separately using the tensor 𝐿̂(𝑟⃗): 2τ ′ 𝑏̂(𝑟) = 𝐿̂(𝑟⃗) ∫0 ⃗⃗⃗⃗ 𝑘 (𝑡)𝑇 ⃗⃗⃗⃗ 𝑘 ′ (𝑡) 𝜕𝑡𝐿̂(𝑟⃗)𝑇 = 𝐿̂(𝑟⃗)𝑏̂𝑠𝑡𝑑 𝐿̂(𝑟⃗)𝑇 .. (III.1.20). The diffusion tensor estimated using the standard DTI approach can be corrected in a similar way: ̂𝑟𝑒𝑎𝑙 (𝑟⃗) = 𝐿̂(𝑟⃗)𝑇 𝐷 ̂ 𝐿̂(𝑟⃗), 𝐷. (III.1.21). ̂𝑟𝑒𝑎𝑙 is the corrected diffusion tensor. where 𝐷 According to equations III.1.11 and III.1.14, the elements of the tensor 𝐿̂(𝑟⃗) can be calculated using the formula: 𝐿̂(𝑟⃗) =. ⃗⃗))𝑖 𝜕(𝐵𝑧 (𝑟 𝜕𝑟𝑗. 𝐺𝑖. ,. (III.1.22). which is identical with equation II.5.2. Therefore, the above considerations also provide the theoretical justification of the gradient coil tensor-based approach (eq. II.5.1). Discussion The Stejskal-Tanner equation is still valid in the case of the PGSE experiment even if the assumption of the gradient uniformity is not fulfilled. In such a case, the b-matrix is spatially variant according to eq. III.1.20. Its real form can be established knowing the matrix L(𝑟⃗) or directly by using the b-matrix spatial distribution in the DTI (BSD-DTI) technique [28]. The former approach has been proposed by Bammer et al. [35]. In that paper, the matrix 𝐿̂(𝑟⃗) was called the “gradient coil tensor” and introduced without derivation. 𝐿̂(𝑟⃗) could be calculated on the basis of the specification provided by the manufacturer or indicated experimentally [36]. Correction using the coil tensor was adapted by D. Malyarenko et al. [45], [46] and E. Tan et al. [34]. The coil tensor described by Bammer et al. is characteristic for a given coil system. It is assumed that the magnetic field in a given voxel is a linear function of the intended gradient, what is not necessarily true in general case. Moreover, it does not take into account such sequence-dependent factors as eddy currents, so a given tensor 𝐿̂(𝑟⃗) is valid only for a given pattern p(𝑟⃗), which in turn, depends on the imaging protocol. Therefore, the robustness of the 33.

(34) field mapping-based methods is limited, since a different sequence is used for calibration and different is for the imaging. A prominent alternative could be the BSD-DTI technique, which is able to map the real form of the b-matrix for a particular imaging protocol. The real magnetic field gradient is not necessarily symmetric. In other words, the absolute values of the effective magnetic field are not necessarily exactly the same, when 2. opposite gradients are intended. In such a case, the integral  k (t )dt for gradient echo-based 0. sequences may not be equal to zero. Therefore, according to equation III.1.17, the signal phase is shifted proportionally to the diffusion coefficient and to the convexity of pattern p(𝑟⃗). According to eq. III.1.13, if the diffusion gradient is assumed to be linear, gradient ∇ϕ is equal to –k’(t) and its divergence is equal to zero. Thus, the classic Stejskal-Tanner equation can be regarded as the special case of eq. III.1.17. Conclusion In the case of spin echo-based sequences, the Stejskal-Tanner equation is valid, even if the assumption of homogeneous gradient field is not fulfilled. In other cases, the equation has an additional imaginary part, with correspond to the signal component shifted in phase by 90°. For the reliable estimation of the diffusion tensor the b-matrix should be derived for each voxel separately, what could be done either by indicating the coil tensor 𝐿̂(𝑟⃗) or by applying the calibration using an anisotropic phantom.. 34.

(35) III.2. Theoretical analysis of the phantom rotation in BSD-DTI Equation II.5.6 is valid under the assumption that the diffusion properties of the phantom are the same across its entire volume. Such high-quality phantoms would be expensive and difficult to manufacture; however, the BSD-DTI method is not limited only to such cases. The diffusion properties of the phantom can be found in the laboratory coordinate frame in three steps: First of all, the phantom is placed inside the scanner in such a way, that the axes of the reference frame associated with the phantom (r, s, t) are aligned with the axes of the laboratory frame (x, y, z). Then, the DTI measurement is performed to establish the spatial distribution of ̂ (𝑞⃗), where 𝑞⃗ = (𝑟, 𝑠, 𝑡) is a vector in the coordinate frame associated with the diffusion tensor 𝐷 the phantom. The application of multiple diffusion gradient directions and additional calibration using an isotropic phantom (e.g. water-ice phantom) of well-known diffusion coefficient, would ensure the highest reliability of this procedure. The second step is to physically rotate the phantom by a given set of Euler angles about its center, and the vector 𝑟⃗ ′ is calculated by multiplying it by an appropriate rotation matrix: 1 0 0 cos (𝛼) 𝑅𝑥 (𝛼) = [ 0 sin(𝛼) cos(𝛽) 0 𝑅𝑦 (𝛽) = [ −𝑠𝑖𝑛(𝛽). 0 −𝑠𝑖𝑛(𝛼)], 𝑐𝑜𝑠(𝛼). (III.2.1a). 0 𝑠𝑖𝑛(𝛽) 1 0 ], 0 𝑐𝑜𝑠(𝛽). (III.2.1b). cos(𝛾) −𝑠𝑖𝑛(𝛾) 0 𝑅𝑧 (𝛾) = [𝑠𝑖𝑛(𝛾) cos(𝛽) 0], 0 0 1. (III.2.1c). for the rotations about the x, y and z axes, by 𝛼, 𝛽 and 𝛾 angles, respectively. It is convenient to describe a rotation about an arbitrary axis as the effect of two or more consecutive rotations about x, y or z axis. In such a case, the order of multiplication corresponds to the order of rotations: 𝑟⃗ 𝑇 = 𝑅𝑥 (𝛼)𝑅𝑦 (𝛽)𝑅𝑧 (𝛾) ∙ 𝑞⃗ 𝑇 = 𝑅(𝛼𝛽𝛾) ∙ 𝑞⃗ 𝑇 .. (III.2.2). ̂ (𝑟⃗) in a voxel of already known The last step is to calculate the diffusion tensor 𝐷 ̂ (𝑞⃗), which was obtained in the first step. It can be done by a coordinates knowing the tensor 𝐷 rotation transformation: ̂ (𝑟⃗) = 𝑅(𝛼𝛽𝛾) ∙ 𝐷 ̂ (𝑞⃗) ∙ 𝑅(𝛼𝛽𝛾)𝑇 . 𝐷 This transformation corresponds to the rotation of diffusion ellipsoid from initial orientation using a set of the Euler angles.. 35.

(36) Considering the phantom’s nonuniformities, eq. II.5.6, that is the basis of BSD-DTI, becomes: 𝐵̂ (𝑟⃗) = −𝑑̂ (𝑟⃗)−1 𝑆̂(𝑟⃗).. (III.2.3). In this work, the BSD calibration involving the above procedure will be referred to as the complete BSD-DTI, while the calibration assuming perfectly homogeneous phantom will be referred to as the uniform BSD-DTI or uBSD-DTI. The analysis described in this section was presented at the 37th Annual International Conference of the IEEE EMBS in Milan and published in the proceedings [47].. 36.

(37) III.3. Computer simulations of BSD-DTI versus standard DTI – overview In the following chapters, the influence of gradient inhomogeneity on the estimated diffusion tensor is analyzed using computer simulations. Although the purpose of each simulated experiment was different, they were all conducted according to a similar scheme. Firstly, the spatial distortion of the diffusion gradient was defined by a pattern function 𝑓⃗(𝑟⃗). Secondly, the BSD-DTI calibration was simulated to indicate the b-matrix distribution 𝑏̂(𝑟⃗). Thereafter, the DTI measurement of an isotropic or anisotropic phantom was simulated and finally, the diffusion tensor was calculated twofold: assuming a spatially constant b-matrix 𝑏̂𝑠𝑡𝑑 (standard DTI or S-DTI) and using one of the BSD-DTI modalities. The general workflow of the simulations is presented in Figure 2.. Figure 2. The general workflow of the simulations.. The idea of analyzing the influence of the diffusion gradient inhomogeneity by means of computer simulations of the standard DTI vs. BSD-DTI was presented at the 34th Annual Scientific Meeting of the ESMRMB in Barcelona [48]. The phantoms The first virtual phantom (P1) used in the simulations was assumed to be a cube made with thin glass plates separated by thin layers of water. The diffusivity in parallel and in perpendicular to the plates was equal to 0.002 mm2/s and 0.0005 mm2/s, respectively. It was used to perform the BSD-DTI calibration and as the object measured during the proper imaging. 37.