creating and evaluating synthetic MR images

(1)

Implementation and verification of a quantitative MRI method for

creating and evaluating synthetic MR images

Implementering och verifiering av en kvantitativ MR-metod för att skapa och utvärdera syntetiska MR-bilder

Aleksander Blagoiev

Faculty of Science, Health and Technology

Degree Project for Master of Science in Engineering Physics 30 HP (ECTS)

Supervisors: Marcus Berg & Jakob Heydorn Lagerlöf Examiner: Lars Johansson

(2)

Abstract

The purpose of this thesis was to implement and quantitatively test a quantitative MRI (qMRI) method, from which synthetic MR images are created and also evaluated. The parameter maps of T₁, T₂^∗, and effective proton density (PD^∗) were tested with reference tubes containing different relaxation times, and concentrations of water (H₂O) and heavy water (D₂O). Two normal volunteers were also used to test qMRI method, by performing regional analysis on the parameter maps of the volunteers.

The synthetic FLASH MR images were evaluated by: using the relative standard deviation of a region of interest (ROI) as a measure for the signal-to-noise ratio (SNR), implanting artificial multiple sclerosis (MS) lesions in the parameter maps used to create the synthetic images, and an MRI radiologist opinion of the images. All MRI measurements were conducted on a 3.0 Tesla scanner (Siemens MAGNETOM Skyra^fit). The results from reference tube testing, shows that the implementation was reasonably successful, although the T₂^∗ maps can not display values on voxels which have T₂ exceeding 100 ms. In vivo parameter map ROI values were consistent between volunteers. The SNR and contrast-to-noise ratio of synthetic images are comparable to their measured counterparts depending on T_E. The artificial MS lesions were distinguishable from normal appearing tissue in a T₁-weighted synthetic FLASH. The radiologist thought the a synthetic T₂^∗-weighted FLASH was somewhat promising for clinical use after further research and development, however a synthetic T₁-weighted FLASH had clinical value.

(3)

Sammanfattning

Syftet med detta arbete var att implementera och kvantitativt undersöka en kvantitativ MR (qMRI) metod, för att sedan skapa och utvärdera syntetiska MR-bilder. qMRI-metodens parameterkartor (T₁, T₂^∗ och effektiv proton densitets PD^∗) undersöktes med olika typer av referensprover. Dessa prover innehöll skilda relaxationstider, samt olika koncentrationer av vatten (H₂O) och tungt vatten (D₂O). In vivo parameterkartor från frivilliga granskades genom att jämföra T₁, T₂^∗ och PD^∗ värdena på intresseområden (ROIs) mellan frivilliga och publicerade värden. Syntetiska FLASH MR-bilder utvärderades genom att: använda relativa standardavvikelsen av ett intresseområde (ROI) som ett mått på signal-brusförhållande (SNR), implantera artificiell multipel skleros (MS) lesioner i de frivilligas parameterkartor för att se ifall dessa kan identifieras i de syntetiska MR-bilder, och slutligen utvärderade en MR-radiolog bilderna. MR-mätningarna utfördes på 3.0 Tesla MR-kamera (Siemens MAGNETOM Skyra^fit). Resultaten från referensproverna visar att implementeringen var rimligen framgångsrik, även om beräknade T₂^∗ för voxlar som har T₂ över 100 ms inte är pålitliga. Frivilligas parameterkartor visade på bra överensstämmelse, dessvärre inte med publicerade. SNR och kontrast-till-brus-förhållandet (CNR) för syntetiska bilder är jämförbara med deras uppmätta motsvarigheter, beroende på T_E. De artificiella MS-lesionerna kunde tydligt skiljas från normal omgivande vävnad i en T₁-viktad syntetisk FLASH. Radiologen tyckte att en syntetisk T₂^∗-viktad FLASH var något lovande för klinisk användning efter ytterligare förbättringar, medan en syntetisk T₁-viktad FLASH hade kliniskt värde.

(4)

Acknowledgments

First and foremost, I want to thank Jakob Heydorn Lagerlöf for giving me the opportunity to do this thesis work; for the conversations related and non-related to the thesis, and the time you have spent helping me perform measurements.

I like to thank Marcus Berg for the valuable discussions and input to the thesis. I have immensely enjoyed the physics anecdotes you have told during this project.

I want to thank the people at the department of medical physics at Örebro university hospital, for allowing us to borrow the phantom and reference tubes used in this thesis.

Last but not least, I like to thank my family for their support during my studies.

(5)

Chapter 1 Introduction

Magnetic resonance (MR) imaging or simply MRI, has superior soft tissue contrast compared to most medical imaging modalities. It is nonhazardous when used properly, and does not expose the patient to ionizing radiation. These are some of the facts that make MRI highly attractive for imaging pathologies and abnormalities in human beings. However, an MRI examination can take 30+ minutes, the scanner itself is expensive and there is also a shortage of trained MRI personal in general; the bottom line being, there is often a deficiency of MRI time in hospitals. There exist different ways of solving this, one in particular is to reduce the MRI scan time, making it possible to examine more patients per day. Today technology exist which decreases scan time, for example parallel imaging, but there are also new ones on the horizon one of which is synthetic MRI. Synthetic MRI is based on quantitative MRI (qMRI), there one measures the physical parameters of tissue which gives a conventional MR image its contrast. Knowing how these parameters are distributed spatially across the patient is called a parameter map. With the parameter maps and an appropriate equation, the signal from an MRI scan can be simulated. This has various advantages, primarily no extra time has to be spent on measuring images with different contrast, since they can be synthesized with similar quality. Synthetic MRI can at least be dated back to 1985 when Bobman et al. published the article

”Cerebral Magnetic Resonance Image Synthesis” [1]. The subject of synthetic MRI has in recent times received new energy, Elster writes: ”The idea [of synthetic MRI] lay relatively dormant until the mid 2000’s when new rapid parameter estimation sequences were developed and integrated into user-friendly commercial products.” [2]. There exist also an open-source qMRI software for data handling and computation of parameter maps, it is called ”hMRI-toolbox” and it mainly uses the pulse sequence FLASH for parameter map estimation [3].

1.1 Aim

The aim of this project was to create synthetic MR human brain images from scratch which are reliable to a reasonable degree, and have comparable image quality to conventional MR images; given constraints to use standard scanner built-in pulse sequences, and limited access to an MRI scanner.

The purpose of this is not to compete with existing software, but rather to see how good the synthetic images can be, given the constraints mentioned.

1.2 Objective

• Based on research by other people, find existing pulse sequences on the MRI scanner for parameter map estimation.

• Create the necessary algorithm for data handling/processing and calculation of parameter maps.

• Test the parameter maps with reference samples to estimate the uncertainty of the implemented qMRI method.

(8)

• Evaluate the parameter map estimation on in vivo brain tissue, by comparing with published values and between volunteers.

• Create synthetic images of the volunteers and evaluate them by: measuring the signal-to-noise, contrast-to-noise ratio and implanting of artificial pathologies.

• Show synthetic images to MRI radiologist for assessment of clinical value.

(9)

Chapter 2 Theory

The discoveries and technological advances that have provided the technology of magnetic resonance (MR) imaging, or simply MRI, are many and to give a comprehensive overview is beyond the scope of this report. Here I will mainly give the minimum which is necessary to understand and criticize the results in this report. The chapter starts out with the basic parts of an MRI scanner, then the phenomenon of nuclear magnetic resonance (NMR), spatial localization of the NMR signal, and lastly miscellaneous topics specific to the thesis.

2.1 Brief Comparison of X-ray Imaging and MRI

Before going into MRI, it is fruitful to know what distinguishes it from other medical imaging principles. Most of us are familiar with x-ray imaging (formerly called projectional radiography):

there ionizing radiation is directed towards a sample¹ and behind it, is a plate which records the transmitted radiation. This creates a shadow image of the electron density of the sample, called a radiograph. The image is 2D although the object is of course 3D. Using computers and multiple radiographs from different angles of the sample, cross-sectional images and 3D reconstruction of the sample can be created; this is called computed tomography (CT). A sample in projection radiography is ” [...] merely reflecting or attenuating the radiation, which is directed onto it [...]” [4, p. 4]. MRI, instead uses non-ionizing radiation (radio frequencies) in a fundamentally different way, explained by Forshult: ”MRI is a true three-dimensional imaging technique, where the object itself is active, responding in various ways to the radiation sent into it – not only onto it.” [4, p. 4-5].

2.2 Basic MRI Hardware

An MRI scanner consists of many parts, but the central ones are perhaps: superconducting coil, gradient coils and radio frequency coils. The superconducting coil generates a static magnetic field, B0, which usually has a magnitude of a few Tesla (T). There exist |B₀| of 0.2 T to 7 T in clinical use, and in research laboratories even 10.5 T for human imaging [5]. The superconductor material is often niobium-titanium (NbTi) alloy which is wound into a coil [6]. For the superconducting coil to be superconducting, it has to be subjected to cryogenic temperatures. Often liquid helium is used as coolant, which boils at roughly 4 Kelvin (K).

B₀ is highly uniform in the region were the patient is imaged, called bore (see Figure 2.1). The reasons for this are given in the report. The location within the bore, which has the most uniform magnetic field is termed isocenter, and it is the preferred position for placing patients when imaging [7, p. 106]. B₀ is not turned off after the MRI scanner is installed [8], if it is desired to turn off B₀ quickly in an emergency, quenching has to be initiated. Then superconductivity is lost and the energy stored in B₀ field is dissipated by boiling off liquid helium.

1Sample, object and patient is used interchangeably in this chapter.

(10)

To increase the homogeneity of B₀ at operation site special hardware and software can be used.

The installation of such hardware and software is termed shimming, and can be both active and passive. The MRI scanner is placed in a room with special lining inside the walls, for instance copper is used as a lining. This means that the MRI scanner is inside a Faraday cage. This is necessary to prevent extraneous radio frequency sources to introduce artifacts on the MR images.

There are three sets of gradient coils, one for each spatial direction (x, y, z).² These coils generate a magnetic field strength in the range mT per meter (mT/m) [10, p. 144], which superposes on B₀. The gradient field is denoted G, and it purposefully results in a nonuniform magnetic field strength inside the bore. With this gradient, one can determine from which spatial part of the patient the NMR signal is departing.

Radio frequency (RF) coils are what emit the RF magnetic field pulses and receive the NMR signal from the sample. There exist coils that perform both functions and coils which only do one of each [9, ch. 9]. The transmitted RF magnetic field is denoted B₁. The primary transmitting coil is the body coil, it encircles the entire bore. There exist plenty of different coils which are used for specific body parts on patient, for instance: head, spine and extremities coils [9, ch. 2].

(a)

Superconducting coil

Gradient coils

RF coil Scanner bore Service connections and

piping for quenching Transverse cross section

x y

(b)

Figure 2.1: (a) Image of a 1.5 T scanner. (b) Typical regions in the scanner which contain the main coils, note that these regions might not be true for the scanner in image (a).

2.3 Nuclear Magnetic Resonance

Often the term magnetic resonance is heard in the context of quantum mechanics, but it is also a phenomenon that can be seen/experienced in the macroscopic world. For instance, take a compass needle that is aligned with Earth’s magnetic field and place an electromagnet perpendicular to it.

Turning the electromagnet on and off at the right time will induce oscillation of the compass needle, that increases with time (assuming no friction). If the compass needle had angular momentum about its axis, then it would be magnetic resonance. A distinct definition is given by Slichter: ”Magnetic resonance is a phenomenon found in magnetic systems that possess both magnetic moments and angular momentum. [...] the term resonance implies that we are in tune with a natural frequency of the magnetic system” [11, p. 1].

The underpinning of MRI is the manifestation of magnetic resonance in the quantum mechanical

2For the z-axis gradient a Maxwell coil is used and for the x and y-axis Golay coils [9, ch. 9].

(11)

regime. One of the most important equations used in NMR and MRI is the Bloch equation d

dtM = M × γB −M_xe_x+ M_ye_y

T₂ −(M_z− M₀)e_z

T₁ .

In Subsection 2.3.7 I explain what the terms mean, but until that, I will derive one part of the Bloch equation (which I call ”simplified Bloch equation”), namely

d

dtM = M × γB,

with quantum mechanics (assuming B = B_ze_z). I also introduce/explain important aspects related to MRI. The simplified Bloch equation can as well be derived with classical physics, but since quantum theory governs atoms and molecules I thought it be most appropriate.

2.3.1 Nuclear Magnetic Dipole Moment and Zeeman Splitting

Atomic nuclei with nonzero spin can experience the phenomenon of nuclear magnetic resonance (NMR). The net spin results in a magnetic field around the nucleus. The first term in the monopole expansion of the magnetic field, magnetic dipole, is concluded from experiment to be [10, p. 26]:

µ = γJ , (2.1)

where J is the net angular momentum of the nucleus, and γ is the gyromagnetic ratio or magnetogyric ratio³. Assuming the nucleus is at rest, then the only contribution to J is the spin (which is also called intrinsic angular momentum) S of the nucleus, meaning J = S. The operator representing S is denoted

S = ˆˆ S_xe_x+ ˆS_ye_y + ˆS_ze_z, (2.2) where e_x, e_y and e_z are Cartesian unit vectors.

Imagine we have a stationary magnetic dipole moment in a homogeneous magnetic field B, the energy of the system is

E = −µ · B, (2.3)

which is a standard expression from electromagnetism. We choose the z-direction along the B:

E = −µ_zB_z. (2.4)

The operator representing the total energy in quantum mechanics is the Hamiltonian operator ˆH.

Because the total energy of the system given by equation (2.4) the Hamiltonian is

H = −ˆˆ µzBz = −γ ˆJzBz= −γ ˆSzBz. (2.5) To keep going, I need to use the eigenfunctions and eigenvalues of the spin operator. For this reason I briefly review some spin physics now: the commutations relations for the spin operators ˆS are [12, ch. 6]

[ ˆSx, ˆSy] = i~ ˆSz, [ ˆSz, ˆSx] = i~ ˆSy, [ ˆSy, ˆSz] = i~ ˆSx. (2.6) From these three relations the eigenfunctions of ˆSz(|s, m_si in bra–ket notation) and eigenvalues (~ms) are

Sˆ_z|s, m_si = ~ms|s, m_si , (2.7) s is the spin of a particle, ~ is Planck’s constant divided by 2π and ms has the allowed values

−s, −s + 1, −s + 2, . . . , +s. If a particle has s = a then the following is true:

|S|²= a(a + 1)~² and Sˆ_i has (2a + 1) eigenstates for i ∈ {x, y, z}.

3See Subsection 2.3.9 for more details on γ.

(12)

Now back to the derivation: solving the time independent Schrodinger equation ( ˆH |E_`i = E_`|E_`i) for the system:

−γB_zSˆ_z|E_`i = E_`|E_`i (2.8)

where γB_z is a constant, we get

|E_`i = |s, m_si and E_`= −γB_z~ms. (2.9) from comparing equations (2.7) and (2.8). Taking an example of s = 1/2 system, the two states m_s = ±1/2 no longer have the same energy – they split. Splitting of energy levels due to a magnetic field is called the Zeeman effect. The state m_s = +1/2 is lower in energy then m_s = −1/2, were the difference is

∆E = −γB_z~

−1

2 − γB_z~ 1

2 = γB_z~. (2.10)

A transition between the two states can happen if the right amount of energy is supplied. When it comes to a photon, its energy has to satisfy ~ω = ∆E, solving for the angular frequency of the photon, ω, gives

ω = γBz (2.11)

This equation is very important in NMR and MRI; the angular frequency is called Larmor frequency after the physicist Joseph Larmor (1857-1942). I will often refer to (2.11) as the Larmor equation.

We will use (2.11) many times, especially when discussing spatial encoding of the NMR signal. The Larmor frequency will also be rederived in a different context. This subsection as been important in many respects, and having established that there is a split in energy of the spin eigenstates, will be important latter when discussing magnetization.

2.3.2 Precession of hSi in a Magnetic Field

Assume a particle with s = 1/2 at rest in a uniform magnetic field B = B_ze_z. Using the standard formula in quantum mechanics, which describes the time rate of the expectation value for an observable [12, ch. 5],

d

dthOi = i

~h[ ˆH, ˆO] i +

*∂ ˆO

∂t +

, (2.12)

where O is an arbitrary observable and ˆH is the Hamiltonian operator of the system. Substituting for O the spin S, gives

d

dthSi = i

~

h[ ˆH, ˆS] i +

*∂ ˆS

∂t +

. (2.13)

The last term in (2.13) is zero since the spin operators have no explicit time dependence. From this we have in component form:

d

dthS_xi = i

~

h[ ˆH, ˆS_x] i , d

dthS_yi = i

~

h[ ˆH, ˆSy] i , d

dthS_zi = i

~

h[ ˆH, ˆSz] i .

(2.14a) (2.14b) (2.14c) From (2.5) we know ˆH = γB_zSˆ_z. Using this with the commutations relations (2.6) results in:

d

dthS_xi = −γB_zi

~

h[ ˆS_z, ˆS_x] i = γB_zhS_yi , d

dthS_yi = −γB_zi

~

h[ ˆS_z, ˆS_y] i = −γB_zhS_xi , d

dthS_zi = −γB_zi

~

h[ ˆSz, ˆSz] i = 0.

(2.15a) (2.15b) (2.15c)

(13)

Writing these set of equations in vector notation d

dthSi = hSi × γB. (2.16)

Equation (2.16) seems to have the same ”structure” as the simplified Bloch equation. In fact, (2.16) is latter used to derive the simplified Bloch equation.

Instructive Solution to Equation _dt^d hSi = hSi × γB

A simple solution to (2.16) is helpful in understanding what the equation describes. We see the expectation value of the z-component is constant, but not x and y, they are coupled. Shuffling the equations (2.15) into the following form:

d²

dt² hS_xi = −γ²B_z²hS_xi , d²

dt²hS_yi = −γ²B_z²hS_yi , d

dthS_zi = 0.

(2.17a) (2.17b) (2.17c) One possible solution for hS_xi, hS_yi and hS_zi is

hS_xi = cos(γB_zt), hS_yi = sin(γB_zt), hS_zi = 1.

(2.18a) (2.18b) (2.18c) This solution is definitely not the most general, however, it still has the main attributes from more general solutions. We see the x and y-component of hSi are oscillating with an angular frequency ω = γB_z. This type of motion, can be described as precession about the magnetic field pointing in the positive z-direction. As I said before, we have derived the Larmor frequency (or Larmor equation) again but in a conceptually different approach. The solution (2.18) does not yields all of the necessary insights needed for NMR. It is possible to retrieve them from equations (2.17), but this is a tedious undertaking. Hence we take a different route which is presented in the next subsection.

2.3.3 Quantum Mechanical Motivation for Flip Angle

The derivation given here is adopted from Brown et al. [10, section 5.3], it differs by using matrix computation. Once again: assume we have a s = 1/2 particle at rest, subjected to a static magnetic field B = B_ze_z. The most general state of this particle is given by a superposition of two possible eigenstates of ˆSz [10, p. 77] (in principle I could have chosen eigenstates of ˆSx or ˆSy):

Ψ(t) = C₊ 1 2,1

2

e⁻~ⁱE+t+ C₋ 1 2,−1

2

e⁻~ⁱE−t. (2.19) C+ and C₋ are complex constants, E_±= ∓γB_z~/2 coming from (2.8). The particle being stationary is the reason for no spatial variables in Ψ(t). Computing the expectation value of S is the current objective. For computational ease matrices are used. The eigenfunctions of ˆS_z are employed as a

”vector” basis, that is to say

|a statei = a 1 2 ,1

2

+ b

1 2 ,−1

2

˙

=

"

a b

#

, (2.20)

where a and b are complex constants. Writing the operators ˆS_x, ˆS_y, ˆS_z in the chosen basis, results in the Pauli matrices:

Sˆ_x= ~ 2

"

0 1 1 0

#

, Sˆ_y = ~ 2

"

0 −i i 0

#

and Sˆ_z = ~ 2

"

1 0

0 −1

#

. (2.21)

(14)

Defining the constants C₊ and C₋ to be

C₊= cos(θ/2)e^iα⁺, C−= sin(θ/2)e^iα⁻ (2.22) where θ, α₊ and α₋ are real. Writing Ψ(t) with the new constants in vector form gives

Ψ(t) =





cos(θ/2)e^−i(α⁺⁻^E+~ t)

sin(θ/2)e^−i(α⁻⁻^E−~ t)



≡ Ψ. (2.23)

Finally, we can compute the expectation value of the spin:

hS_xi = hΨ(t)| ˆS_x|Ψ(t)i = Ψ^†~ 2

"

0 1 1 0

# Ψ = ~

2sin(θ) cos(α₋− α₊− γB_zt), hS_yi = hΨ(t)| ˆS_y|Ψ(t)i = Ψ^†~

2

"

0 −i i 0

# Ψ = ~

2sin(θ) sin(α₋− α₊− γB_zt), hS_zi = hΨ(t)| ˆS_z|Ψ(t)i = Ψ^†~

2

"

1 0

0 −1

# Ψ = ~

2cos(θ).

(2.24a)

(2.24b)

(2.24c) Where † is the conjugate transpose (adjoint). Defining φ = α− − α₊ and ω = γB_z gives a more pleasant expression, namely:

hS_xi = ~

2sin(θ) cos(φ − ωt), hS_yi = ~

2sin(θ) sin(φ − ωt), hS_zi = ~

2cos(θ).

(2.25a) (2.25b) (2.25c) The new things here, compared to equations (2.18), are the phase factor φ and time independent sinusoidal prefactor with θ as argument. We clearly see θ to be the angle that hSi makes with the z-axis, and φ the initial azimuthal angle at t = 0. A schematic of this is shown in Figure 2.2. θ is a helpful aid when discussing an important quantity in MRI, namely flip angle (also called tip angle).

For this reason, have I chosen to denote flip angle with θ in this report.

x

y z

θ

φ hS(t)i

hS(t = 0)i B

Figure 2.2: Precession of the spin expectation value.

(15)

2.3.4 Oscillating Magnetic Field Effect on Spin

Now we will see how the expectation value of the spin, hSi, behaves when an oscillating magnetic field is applied to it. The full details of the derivation is given by Brown et al. [10, section 5.4]. I only set up the derivation and show the main result from [10, section 5.4].

If you have forgotten: assume we have a s = 1/2 particle at rest, subjected to the time dependent magnetic field

B(t) = B₁(t) + B₀ = B₁(cos(ωt)e_x− sin(ωt)e_y) + B₀e_z. (2.26) From equation (2.5) we observe the Hamiltonian is time dependent:

H(t) = − ˆˆ µ · B(t) = −γ ˆS · B(t) = γSˆxB1cos(ωt) − ˆSyB1sin(ωt) + ˆSzB0

. (2.27) Using matrices and introducing

Sˆx⁰ = ˆSxcos(ωt) − ˆSysin(ωt) = ~ 2

"

0 1 1 0

#

cos(ωt) −~ 2

"

0 −i i 0

#

sin(ωt) = ~ 2

"

0 e^iωt e^iωt 0

#

, (2.28) gives

H(t) = γBˆ ₁~ 2

"

0 e^iωt e^iωt 0

#

+ γB₀~ 2

"

1 0

0 −1

#

= ~ 2

"

ω0 ω1e^iωt ω₁e^iωt ω₀

#

, (2.29)

where ω₀= γB₀ and ω₁= γB₁. We have to solve the Schrodinger equation for the system, i~d

dtΨ(t) = ~ 2

"

ω₀ ω₁e^iωt ω1e^iωt ω0

#

Ψ(t), (2.30)

to get the wave function, in order to determine how hSi behaves. By making the ansatz

Ψ(t) = ψ₁⁰(t)e^iω⁰^t/2 1 2,1

2

+ ψ₂⁰(t)e^−iω⁰^t/2 1 2,−1

2

=

"

ψ₁⁰(t)e^iω⁰^t/2 ψ₂⁰(t)e^−iω⁰^t/2

#

(2.31) and assuming B₁ is in resonance with the system, meaning ω = ω₀, the wave function Ψ(t) can be determined from the Schrodinger equation. When we have the solution and define

Sˆ_y⁰ = ˆS_xsin(ωt) + ˆS_ycos(ωt), (2.32) we can compute hS_x⁰i,S_y⁰and hS_zi. The results of this is [10, p. 82]

hS_x⁰i = hΨ(t)| ˆS_x⁰|Ψ(t)i = hS_x⁰(t = 0)i ,

hS_y⁰i = hΨ(t)| ˆS_y⁰|Ψ(t)i = hS_y⁰(t = 0)i cos(ω₁t) + hSz(t = 0)i sin(ω₁t), hS_zi = hΨ(t)| ˆS_z|Ψ(t)i = − hS_y⁰(t = 0)i sin(ω₁t) + hS_z(t = 0)i cos(ω₁t).

(2.33a) (2.33b) (2.33c) Where the initial conditions are

hS_x⁰(t = 0)i = ~

2sin(θ) cos(φ), hS_y⁰(t = 0)i = ~

2sin(θ) sin(φ), hS_z(t = 0)i = ~

2cos(θ).

(2.34a) (2.34b) (2.34c) The initial conditions describes a point on a sphere with radius ~/2, polar angle θ and azimuthal angle φ. From the definition of ˆSx⁰ and ˆSy⁰, we see that they are representing S components in a coordinate system that is rotating about the z-axis at angular frequency ω. This means S_x⁰, S_y⁰ and S_z are the components of S in the coordinate system (x⁰, y⁰, z) which rotates about the z-axis of the stationary

(16)

coordinate system (x, y, z). To drive the point home: the difference between S_x⁰, S_y⁰ and S_x, S_y is that the former is the transverse component of S in a rotating coordinate system, while the latter is the transverse component of S in a static system. It is common to denote (x, y, z) as the laboratory frame and (x⁰, y⁰, z) as rotating frame.

Equations (2.33) express that hSi has a constant magnitude, and (2.33a) means the x⁰ component of the spin is unchanging, while the y⁰ and z are sinusoidally changing with frequency ω₁. This motion of hSi is called precession about the x⁰-axis. The point of this subsection was to give a quantum mechanical motivation for ”spin tipping” (which is important in both MRI and NMR), where the expectation value of the spin is tipped from equilibrium with an RF pulse consisting of the resonance frequency of the spin system.

2.3.5 Magnetic Dipole Moment and Magnetization

I hope that in the two previous sections, I have convinced the reader that the expectation value of the spin, hSi, behaves rather classically. All of this could in principle have been omitted by stating the Ehrenfest theorem: expectation values of observables follow classical laws. I believe spin to be such a strange quantity in physics that the analysis can be useful for intuition.

From equation (2.1) we see µ behaves similar as S in a static and oscillating magnetic field (even if γ is negative, the characteristic behavior shown above still applies). We can substitute µ/γ for S in equation (2.16) giving

d

dthµi = hµi × γB. (2.35)

We have almost derived the simple Bloch equation, to get there we have to define a macroscopic quantity, namely the magnetization of a volume element⁴ containing many non-interacting dipole moments:

M = 1

∆V

*_∆V X

i

µ_i +

= 1

∆V

X

i

hµ_ii , (2.36)

where ∆V is the volume of the volume element and the sum includes all µ_i in ∆V . Because each term in M obeys (2.35) then M follows it:

d

dtM = M × γB. (2.37)

The magnitude of M depends on how many µ are pointing in the same direction. This can be analyzed with studying how the two states m_s= ±1/2 are populated by protons in a volume element.

2.3.6 How the Spin 1/2 Eigenstates are Populated

Suppose we have many s = 1/2 non interacting protons at rest which are in a uniform magnetic field B₀ = B₀e_z. There are two distinct energy states for S_z, namely m_s = ±1/2 with energy E± = ∓γ~B0 respectively. Now the question is: how are these two states populated in equilibrium?

In this case it would be most logical to use the Fermi-Dirac distribution, since it describes the statistics of these systems. However, in the classical limit the Fermi-Dirac function can be approximated by the Maxwell-Boltzmann distribution. Using some classical statistical mechanics arguments [13] the equation

N_−1/2

N_+1/2 = e^−γ~B⁰^/k^B^T, (2.38)

can be derived [10, section 6.2]. Where N_±1/2 is the amount of protons in the states m_s = ±1/2, respectively. The equation assumes equilibrium has been reached. Substituting numbers for the constants:

N_−1/2

N_+1/2 = exp −267 · 10⁶ Hz/T · 1.0 · 10⁻³⁴Js · 3 T 1.38 · 10⁻²³J/K · 310 K

!

≈ 0.99998. (2.39)

4Volume element and voxel are used interchangeably in this report. The ”vo” in voxel stands for ”volume”, and ”el”

for ”element”.

(17)

Meaning

N_+1/2≈ 1

0.99998N−1/2≈ 1.000019N_−1/2. (2.40) There are roughly 19 parts per million (ppm) more protons in state m_s = +1/2 than m_s = −1/2.

In other words, for every million of hydrogen nuclei there are 19 more in state m_s = +1/2 than ms = −1/2. There would be no magnetization if there was an equal number of protons in the states m_s = ±1/2. But from above we have see there is a net amount of protons in state m_s = 1/2, resulting in a magnetization at nonzero temperatures.

2.3.7 Bloch Equation

MRI is based on the macroscopic magnetization M = M_xex + M_yey + M_zez of a volume element of an object, resulting from B₀. The phenomenological equation⁵ governing M is named the Bloch equation:⁶

d

dtM = M × γB −Mxex+ M_yey

T₂ −(M_z− M₀)e_z

T₁ , (2.41)

where M₀ is the equilibrium magnetization due to B₀. B is the external magnetic field environment of the volume element. The last two terms in the Bloch equation are relaxation terms, they describe the return of M to equilibrium (assuming B only has a z-component) after M has been tipped away from it [10, p. 60]. Figure 2.3 illustrates this. The tipping is done with an RF pulse that contains the resonance frequency of the volume element. The M polar angle after the RF pulse is called flip angle, and I denote it with θ. Usually the intended flip angle is given as a prefix to the RF pulses, for example: 30^◦ RF pulse. T₁ is called longitudinal (or spin-lattice) relaxation time, T₂ transverse (spin-spin) relaxation time. If both T₁ and T₂ → ∞, then the Bloch equation reduces to (2.37) which we derived with quantum mechanics.

The Bloch equation was first proposed by Felix Bloch (1905-1983) in the article ”Nuclear Induction”

[14].⁷ Slichter gives the following comment about the Bloch equations: ”Although they have some limitations, they have nevertheless played a most important role in understanding resonance phenomena, since they provided a very simple way of introducing relaxation effects.” [11, p. 35].

x

y z

M

B

θ

Figure 2.3: Recovery of the magnetization to equilibrium where the blue curve illustrates the trajectory;

assuming B = Be_z.

5This means the equation is empirically motivated.

6Also known as Bloch equations, referring to each vector component of the equation.

7In the article Bloch uses the H-field instead of B-field.

(18)

T₁ (R₁)

From the Bloch equation it follows that the T₁ relaxation time represents the M_z desire to return to equilibrium M₀ (given B = Be_z). A quantum mechanical motivation of T₁ relaxation is given by Brown et al [10, section 6.2] and Slichter [11, ch. 1]. The main idea of the motivation goes along these lines: for a transition between the states m_s= ±1/2 to take place, there has to be some sort of energy source or sink. It provides/takes energy, making the transition energetically possible. In solids the energy source/sink is the lattice, and for liquids it is molecules. This means T₁ relaxation, physically is an energy transfer process with the surrounding chemical environment. The inverse of T₁ is denoted R₁.

T₂ (R₂) and T₂^∗ (R^∗₂)

From the Bloch equation it follows that the T₂relaxation time represents the decay of M xy-component, given B = Be_z. T₂ relaxation stems mainly from fluctuations in the magnetic field within a volume element, stemming from Brownian motion of the constituent atoms/molecules and from nuclei who are interacting via their magnetic fields (given they have nonzero spin) [10, p. 620]. These local field fluctuations change the protons Larmor frequency, and if the individual hµi where in phase at some time, then they would dephase because of this.

Inhomogeneity of B₀ and other external fields (such as magnetic fields from surrounding tissues of patients) will also contribute to changes in Larmor frequency. The parameter T₂^∗ (reduced transverse relaxation time) combines the last mentioned complication and T₂relaxation [10, p. 9]. T₂^∗is dependent on experimental conditions, this makes it less useful than T₁ and T₂; nevertheless T₂^∗ is relevant in MRI. T₂^∗ is less then or equal to T₂, and the inverse of T₂^∗ is R^∗₂; likewise R₂ ≡ 1/T₂.

Table 2.1: Approximate T₁ and T₂ values of human tissue at 1.5 T and 37 ^◦C [10, p. 56].

Human Body Tissue T₁ (ms) T₂ (ms) Brain gray matter (GM) 950 100 Brain white matter (WM) 600 80 Cerebrospinal fluid (CSF) 4500 2200

Muscle 900 50

Fat 250 60

Oxygenated blood 1200 200

Deoxygenated blood 1200 100

2.3.8 Measuring the NMR Signal

How do we measure the NMR signal from a volume element of an object? It is based on Faraday’s law

∇ × E = −∂

∂tB, (2.42)

where E is the electric field. From Faraday’s law we can derive V = −d

dtΦ, (2.43)

assuming the area related to Φ is time independent. V is the electromotive force (emf) and Φ is the magnetic flux through a surface, defined by

Φ = Z Z

B · ds. (2.44)

The magnetizations M trajectory after it has experienced an RF pulse (which contains the resonance frequency of the volume element) is a spiraling motion towards equilibrium. The motion implies a

(19)

changing magnetic field, and having a coil in the vicinity (perpendicular to the equilibrium direction, ideally) means we have a change in magnetic flux through the coil; resulting in an emf that can be measured, Figure 2.4 illustrates this. From Figure 2.4 it is clear that only the transverse (xy) component of M produces a changing flux. This means only the xy-component results in signal. It is important to have enough signal to overcome the noise, and from subsection 2.3.6 we know that only few ppm more occupy the eigenstate m_s = +1/2 than m_s= −1/2. To get a sufficient signal then, we need many nuclei and the high abundance of hydrogen in biological tissue makes it the first choice for imaging humans. As a short side note: Bloch’s article [14] received its name ”Nuclear Induction”, due to the induced emf one can measure [15].

x

y z

M

M B

time Voltage

Figure 2.4: Illustration of how the NMR signal is generated.

2.3.9 Gyromagnetic Ratio

Here I present a rough formula for the gyromagnetic ratio γ, it will give a qualitative understanding of how mass affects γ. We will arrive at the formula by using hand-waving arguments.

Imagine an electron with charge −e, mass m and zero spin orbiting a stationary nucleus at radius r. The time it takes for the electron to complete one revolution is T . The resulting angular momentum about the nucleus is

L = r × mv = re_r× m2πr

T e_φ= 2πmr²

T e_z. (2.45)

The magnetic dipole moment for the system is [16, p. 243]

µ = IAe_r× e_φ= −e

T πr²e_z. (2.46)

Equation (2.1) for the system we have currently is

µ = γL. (2.47)

Substituting the expression for µ and L in the equation, and solving for γ yields:

γ = −e

2m. (2.48)

Expression (2.48) should only be taken as a qualitative formula for the gyromagnetic ratio of stationary charged particles, since the underlying assumptions are not applicable. Although for a stationary electron with spin 1/2 the equation (2.48) is off by a factor of 2 roughly (called electron g-factor).

For a stationary proton it is even more off compared to experiment. Still the tendency of lower γ for

(20)

larger mass is correctly given by (2.48). This tendency can be seen in Table 2.2 where experimentally determined gyromagnetic ratios are shown. The fact that different nuclei possess different γ makes it possible to distinguish distinct chemical elements in NMR spectroscopy. It also makes it in principle possible to mix heavy water (D₂O) with water (H₂O) and only get signal from one species. This is the idea behind one type of reference tube used in this thesis.

A short side note: relativistic quantum theory predicts the electron g-factor to 2 of a free electron [17]. This is consistent with experiment up to the one thousandth decimal place. However, using quantum electrodynamics (QED) many of the decimals have been correctly determined. In fact to such an extent, that it is classified as the most accurate predication of any physics theory ever made (at the time of writing) [17]. The exact physical cause for the dipole moment of electrons and protons seems to be unknown today. There are some theoretical models that try to predict the protons dipole moment by studying currents of quarks within the proton [15].

With what we have uncovered here and the Larmor equation, we see a way to identify atomic nuclei with nonzero spin: put them in a uniform magnetic field, excite the nuclei and measure the frequencies of light they emit. Sometimes the chemical environment (electrons) of the nuclei can reduce B₀ in the region of the nucleus, which shifts the Larmor frequency. This shift can be measured and used to determine the chemical environment of the nucleus. Both elemental and structural identification of a sample is performed in NMR spectroscopy by these principles.

Table 2.2: The gyromagnetic ratio for some nuclei. Negative sign for γ/2π means that the magnetic moment is anti-parallel to the spin vector. Data from source 13 within [4].

Nucleus Spin s γ/2π (MHz/T)

Hydrogen ¹H 1/2 42.58

Deuterium D (²H) 1 6.54

Carbon ¹³C 1/2 10.71

Nitrogen¹⁴N 1 3.08

Oxygen¹⁷O 5/2 −5.77

2.4 Tissue Parameters

The tissue parameters⁸ which have been used in this thesis are: T₁, T₂^∗ and PD^∗ (also denoted ρ^∗).

The two former have been introduced in Subsection 2.3.7. Here I explain what PD and PD^∗ are, also a short description of hydrogen in biological tissue relevant for MRI.

Proton density (PD) or spin density is a parameter that measures the NMR active hydrogen nuclei in tissue. It is a dimensionless quantity given in percent units (p.u. or pu), where 100 pu is the PD value for H₂O or cerebrospinal fluid (CSF). PD^∗ denotes that data which is used to computed PD, has been affected by T₂^∗ effects. This means that PD and PD^∗ are not identical, but related.

2.4.1 Hydrogen

The hydrogen nucleus consists only of a single proton⁹, and it is the nucleus of choice for clinical MRI due to: most biological tissues have hydrogen as their main constituent [4, p. 3];in other words, there are plenty of hydrogen nuclei (which is needed for a sufficiently strong NMR signal because only few ppm of protons contribute).

Hydrogen atoms in biological tissues exist in different molecules, for instance H₂O, proteins and fat. A hydrogen atom can perhaps rotate or move more in molecule X then in molecule Y, this has some effect on the relaxation times. For this reason, not all hydrogen nuclei contribute equally to the

8A note to prevent confusion: tissue T1, T2 and PD values are often referred to as tissue parameters. In principle we can change them, but most of the time they are thought as being intrinsic (not changing).

9I use the words hydrogen nucleus and proton interchangeably throughout this report.

(21)

NMR signal. Intermolecular forces further complicate the situation, for instance: an H₂O molecule can be tightly bound to a macromolecule, which affects the H₂O molecules ability to move and rotate [18]. For this reason, the notion of free and bound H₂O exist [18]. It turns out that the nuclei in free and loosely bound H₂O are the main contributors to most types of NMR signal from biological tissues [19].

2.5 k-space in MRI

Fourier space or k-space is a concept which can be found in many areas of study: from digital image processing to solid state physics, and MRI is not an exception. k-space is the mathematical space where a Fourier transformed function ”lives”. The Fourier transform is perhaps mostly thought of as a mathematical transform which decomposes frequencies (wavenumbers) of a signal (image) in time (real space) domain to Fourier space. A useful way of viewing a MRI scanner, is by thinking of it as measuring an objects Fourier transform. The justification for this point of view is presented below.

Some important properties of k-space are mentioned and concepts relating k-space with real space are also explained in this section.

2.5.1 Mathematical Motivation for Using k-space in MRI

The NMR active proton density (PD or ρ) is generally a function of r = xe_x+ye_y+ze_zin 3D, denoted ρ(r). Now, think of the magnetization vector being tipped into the transverse plane. There it is precessing at the Larmor frequency, assuming infinite relaxation times. We can approximate the signal received in the receiver coil from a infinitesimal volume element d³r at position r₀ = x₀e_x+y₀e_y+z₀e_z as

dS(r₀, t) = ρ(r₀) cos(φ(r₀, t)) = ρ(r₀) cos(−φ(r₀, t)) = ρ(r₀) Reⁿe^−iφ(r⁰^,t)^o. (2.49) Here φ(r₀, t) is the phase of the precessing magnetization. I drop the Re{. . . } notation, but it is implicitly understood that only the real part of the phasor, ρ(r)e^iφ(r,t), is meant.

By definition

dφ(r, t)

dt = ω(r, t) ⇐⇒ φ(r, t) = Z t

t0

ω(r, t⁰)dt⁰. (2.50) Where t₀ in the integral is the time from which you start measuring the phase. I set t₀ = 0 since I only want to measure the phase after the RF excitation, which I define to be t = 0. From the Larmor equation we know

ω(r, t) = γ|B(r, t)| = γB(r, t), (2.51)

where B(r, t) is the external magnetic field of the MRI scanner, stemming from the gradient and superconducting coils. Assuming the gradient coils can create a linearly varying B(r, t), we get

B(r, t) = B0+ G_x(t)x + G_y(t)y + G_z(t)z = B₀+ G(t) · r. (2.52) Writing the phase in terms of magnetic fields:

φ(r, t) = γ Z t

0

B0+ G(t⁰) · rdt⁰ = γB₀t + γ Z t

0

G(t⁰)dt⁰· r = ω₀t + 2πk · r, (2.53) where

ω0= γB₀ and k = γ 2π

Z t 0

G(t⁰)dt⁰. (2.54)

Using the expression for φ(r, t) in (2.53), and write the signal from an arbitrary volume element:

dS(r, t) = ρ(r)e^−i(ω⁰^t+2πk·r). (2.55)

creating and evaluating synthetic MR images

Implementation and verification of a quantitative MRI method for