Close
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected
Invert selection
Deselect all
Deselect all
Click here to refresh results
Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Improved brain dynamic contrast enhanced MRI using modelbased reconstruction
(USC Thesis Other)
Improved brain dynamic contrast enhanced MRI using modelbased reconstruction
PDF
Download
Share
Open document
Flip pages
Contact Us
Contact Us
Copy asset link
Request this asset
Transcript (if available)
Content
Improved Brain Dynamic Contrast Enhanced MRI
using ModelBased Reconstruction
by
Yi Guo
A Dissertation Presented to the
FACULTY OF THE USC GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFORNIA
In Partial Fulfillment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
(ELECTRICAL ENGINEERING)
August 2017
Dedicated to my parents and my wife
ii
Acknowledgments
A great many people have helped me and guided me in completing this thesis, without
them it is never possible for me to achieve what I have done today, and I would hope to
express my deepest gratitude to these people with this opportunity.
I would like to offer my deepest and most sincere gratitude to my advisor and mentor,
Professor Krishna S. Nayak. He is truly an exceptional researcher and an excellent teacher,
who taught me many things in research and in personal development. He is always
optimistic, funny, and helpful, and his guidance and encouragement help me overcome
many difficulties and achieve multiple milestones throughout my PhD. I am forever
indebted to him for his help and support in completing my graduate studies.
I would also like to thank Professor Justin P. Haldar for serving on my doctoral
committee. His rigorous mathematical knowledge and insight has helped me greatly in
times I am stuck with algorithm design and performance. I also learned a great many MRI
imaging and reconstruction knowledge from his MRI class and vector space class.
I am also thankful to Professor Hao Li for serving on my doctoral committee. His
valuable comments and feedback offered me insightful thoughts on my project, and helped
me think the fundamentals of my research topic and the potential impact.
I am thankful to Dr Meng Law, who helped me with scheduling and scanning patients
in USC hospital. I am also grateful to learn many clinical knowledge and radiology practice
from the helpful discussion with him.
I am also thankful to Dr Mark S. Shiroishi. He is kind and patient in helping me rate
image quality of our proposed MRI technique. I also learned a lot about the clinical side of
MRI imaging through the collaboration with him in scanning patients.
I am also very grateful to Dr R. Marc Lebel, who helped and collaborated with me
throughout my entire PhD study. I have inherited his postdoc project at USC for a smooth
start of my PhD journey, and he offered tremendous help in many aspects of my research
with his solid knowledge in both clinical and mathematical part of MRI.
I would also like to give my special thanks to Dr Sajan Goud Lingala, who collaborated
closely with me in the DCEMRI project. His solid knowledge in MRI image
reconstruction helped me immensely in developing and debugging my algorithms. He is
also very kind, patient and helpful, making every discussion easy and fruitful.
I am also privileged to be part of the lively and friendly MREL family. I am grateful to
every past and present member. I would like to thank Weiyi Chen, Ahsan Javed, Vanessa
Landes, Yannick Bliesener, Yongwan Lim, Ziyue (Brian) Wu, Hung Phi Do and Eamon
Doyle for their continuous help and support. I also hope to give my special thanks to
Terrence Jao, Yinghua Zhu, and Xin Miao, who helped me greatly in both my research and
personal life. I would also like to thank my friends outside the lab for helping me overcome
various difficulties during my stay here at USC.
iii
Most importantly, I am deeply grateful to my parents and my wife for their unfailing
love, support and the sacrifice they made to help me fulfill my PhD study and realize my
dream.
iv
Table of Contents
Acknowledgments .......................................................................................................... ii
List of Publications ....................................................................................................... vi
Abbreviation .................................................................................................................. ix
Abstract .......................................................................................................................... x
Chapter 1. Introduction ................................................................................................ 1
Motivation .......................................................................................................... 1
Dynamic Contrast Enhanced Magnetic Resonance Imaging (DCEMRI) ........ 3
Organization of the Dissertation ........................................................................ 5
Chapter 2. Background ................................................................................................ 6
MRI Physics ....................................................................................................... 6
2.1.1 MR Signal Generation ................................................................................. 6
2.1.2 T1 and T2 Relaxation ................................................................................... 7
2.1.3 MRI Signal Space........................................................................................ 9
2.1.4 Parallel Imaging ........................................................................................ 10
DCEMRI Basics ............................................................................................. 11
2.2.1 Calibration Scans....................................................................................... 12
2.2.2 DCEMRI Scan ......................................................................................... 15
2.2.3 Preprocessing for TK modelling .............................................................. 16
2.2.4 Trackerkinetic (TK) Modeling ................................................................. 18
2.2.4.1 Common TK Models for Brain Tumor Evaluation ............................. 19
2.2.4.2 Model Selection ................................................................................... 20
Sparse Sampling and Constrained Reconstruction .......................................... 20
2.3.1 Sparse Sampling ........................................................................................ 21
2.3.2 Constrained reconstruction ........................................................................ 22
2.3.2.1 l2 Norm Constraint ............................................................................... 23
2.3.2.2 Compressed Sensing (l0 and l1 Constraint) .......................................... 24
Chapter 3. Clinical Evaluation of Constrained Reconstruction in Brain DCEMRI . 27
Introduction ...................................................................................................... 27
Materials and Methods ..................................................................................... 29
3.2.1 HighResolution WholeBrain DCEMRI ................................................ 29
3.2.2 Experimental Methods .............................................................................. 31
3.2.3 Comparisons / Evaluation ......................................................................... 35
Results .............................................................................................................. 36
3.3.1 Experimental Results................................................................................. 36
3.3.2 Quantitative Assessment ........................................................................... 40
3.3.3 Qualitative Assessment ............................................................................. 41
v
Discussion ........................................................................................................ 43
Conclusion ....................................................................................................... 45
Chapter 4. Modelbased Direct Reconstruction for DCEMRI ................................. 47
Introduction ...................................................................................................... 47
Theory .............................................................................................................. 48
4.2.1 Direct TK Mapping ................................................................................... 48
4.2.2 Indirect TK Mapping ................................................................................. 51
Methods ............................................................................................................ 51
4.3.1 Digital Phantom......................................................................................... 51
4.3.2 InVivo Retrospective Evaluation ............................................................. 52
4.3.3 InVivo Prospective Evaluation ................................................................ 54
Results .............................................................................................................. 54
Discussion ........................................................................................................ 62
Conclusion ....................................................................................................... 65
Appendix 4A ................................................................................................... 65
Chapter 5. Joint Arterial Input Function and Tracker Kinetic Parameter Estimation
from UnderSampled DCEMRI using a Model Consistency Constraint ........................ 67
Introduction ...................................................................................................... 67
Theory .............................................................................................................. 68
5.2.1 Model Consistency Constraint .................................................................. 68
5.2.2 Joint AIF and TK Parameter Estimation ................................................... 69
5.2.3 Theoretical Benefits .................................................................................. 69
Methods ............................................................................................................ 70
5.3.1 Data Sources .............................................................................................. 70
5.3.2 Demonstration of TK Solver Flexibility ................................................... 71
5.3.3 Demonstration of TK Model Flexibility ................................................... 71
5.3.4 Demonstration of Joint AIF and TK estimation ........................................ 71
5.3.5 Demonstration with Prospectively Undersampled Data ........................... 72
Results .............................................................................................................. 72
Discussion ........................................................................................................ 80
Conclusion ....................................................................................................... 80
Appendix 5A ................................................................................................... 81
Chapter 6. Conclusion ............................................................................................... 83
References .................................................................................................................... 85
vi
List of Publications
Journal Articles
1. Yi Guo, Sajan Goud Lingala, Yannick Bliesener, R. Marc Lebel, Yinghua Zhu,
Krishna S. Nayak. Joint arterial input function and tracker kinetic parameter
estimation from undersampled DCEMRI using a model consistency constraint.
Magnetic Resonance in Medicine (in revision)
2. Sajan Goud Lingala, Yi Guo, Yinghua Zhu, Robert Marc Lebel, Meng Law, Krishna
S. Nayak. Tracer Kinetic Models as Temporal Constraints during DCEMRI
reconstruction, Medical Physics (in revision)
3. Yi Guo, Sajan Goud Lingala, Yinghua Zhu, R Marc Lebel, Krishna S Nayak. Direct
Estimation of TracerKinetic Parameter Maps from Highly Undersampled Brain
DCEMRI, Magnetic Resonance in Medicine DOI: 10.1002/mrm.26540
4. Yi Guo, R Marc Lebel, Yinghua Zhu, Sajan Goud Lingala, Mark S. Shiroishi, Meng
Law, Krishna S. Nayak. Highresolution Wholebrain DCEMRI Using Constrained
Reconstruction: Prospective Clinical Evaluation in Brain Tumor Patients, Medical
Physics, 43 (2016): 20132023
5. Y Zhu, Y Guo, SG Lingala, RM Lebel, M Law, KS Nayak. GOCART: GOlden
angle CArtesian Randomized Timeresolved 3D MRI. Magnetic Resonance
Imaging 34.7 (2016): 940950
6. Xin Miao, Yi Guo, Terrence Jao, Muhammad Usman, Claudia Prieto, Sajan Goud
Lingala, Krishna S Nayak. Accelerated Cardiac CINE MRI Using Locally Low
Rand and Finite Difference Constraints, Magnetic Resonance Imaging 34.6 (2016):
707714.
Invited papers
1. Yi Guo, Yinghua Zhu, Sajan Goud Lingala, R. Marc Lebel, Mark Shiroishi, Meng
Law and Krishna Nayak. Highresolution wholebrain dynamic contrastenhanced
MRI using compressed sensing, 3 August 2015, SPIE Newsroom. DOI:
10.1117/2.1201507.006016
Conference Papers
1. Yi Guo, Sajan Goud Lingala, R. Marc Lebel, Krishna S. Nayak. Joint estimation of
arterial input function and tracer kinetic parameters from undersampled DCEMRI,
ISMRM 2017, p3751 (Oral Power Pitch Presentation)
vii
2. Yi Guo, Sajan Goud Lingala, Krishna S. Nayak. Reconstruction of DCE tracer
kinetic parameters from undersampled data with a flexible model consistency
constraint, ISMRM 2017, p2301 (Oral Presentation)
3. Yi Guo, Yinghua Zhu, Sajan Goud Lingala, R. Marc Lebel, and Krishna S. Nayak,
Direct Reconstruction of Kinetic Parameter Maps in Accelerated Brain DCEMRI
using the ExtendedTofts Model, ISMRM 2016, p0868 (Oral presentation)
4. Sajan Goud Lingala, Yi Guo, Yinghua Zhu, Naren Nallapareddy, R. Marc Lebel,
Meng Law, and Krishna Nayak. Accelerated brain DCEMRI using Contrast Agent
Kinetic Models as Temporal Constraints, ISMRM 2016, p0651 (Oral presentation)
5. Yi Guo, Yinghua Zhu, Sajan Goud Lingala, R. Marc Lebel, and Krishna S. Nayak,
Direct Reconstruction of TracerKinetic Parameter Maps from Prospective Highly
Undersampled DCEMRI, ISMRM data sampling & image reconstruction
workshop, 2016 (Oral presentation)
6. Sajan Goud Lingala, Yi Guo, Yinghua Zhu, Naren Nallapareddy, R. Marc Lebel,
Meng Law, Krishna S. Nayak, Accelerated DCE MRI using contrast agent kinetic
models as temporal constraints, ISMRM data sampling & image reconstruction
workshop, 2016 (Oral presentation)
7. S.G.Lingala, Y.Mohsin, S.Bhave, X.Miao, Y.Guo, K.S.Nayak, E.DiBella, M.Jacob,
Dataadaptive reconstruction algorithms for accelerated dynamic MRI: an open
source MATLAB package, ISMRM reconstruction workshop, 2016.
8. KS Nayak, Y Guo, Y Zhu, SG Lingala, RM Lebel, N Nallapareddy, MS Shiroishi,
M Law. "Improved clinical DCEMRI pipeline for high resolution, whole brain
imaging: application to brain tumor patients." RSNA 2015, Chicago, Oral
Presentation: #SSA1804. (Oral presentation)
9. Yi Guo, Yinghua Zhu, Sajan Goud Lingala, R. Marc Lebel, and Krishna S. Nayak,
Highly Accelerated Brain DCE MRI with Direct Estimation of Pharmacokinetic
Parameter Maps, ISMRM 2015, p0573 (Oral presentation)
10. Yi Guo, R. Marc Lebel, Yinghua Zhu, Mark S. Shiroishi, Meng Law, and Krishna
S. Nayak, Highresolution Wholebrain DCE MRI of Brain Tumor using
Constrained Reconstruction: Prospective Clinical Evaluation, ISMRM 2015, p3050
11. Xin Miao, Sajan Goud Lingala, Yi Guo, Terrence Jao, and Krishna S. Nayak,
Accelerated Cardiac Cine Using Locally Low Rank and Total Variation Constraints,
ISMRM 2015, p0571 (Oral presentation)
12. Yinghua Zhu, Yi Guo, Sajan Goud Lingala, R. Marc Lebel, Meng Law, Krishna
Nayak, Evaluation of GLACIER sampling for 3D DCEMRI, ISMRM 2015, p2535
viii
13. Yinghua Zhu, Yi Guo, Sajan Goud Lingala, Samuel Barnes, R. Marc Lebel, Meng
Law, Krishna Nayak, Evaluation of DCEMRI Data Sampling, Reconstruction and
Model Fitting Using Digital Brain Phantom, ISMRM 2015, p3052
14. Sajan Lingala, Yi Guo, Yinghua Zhu, Samuel Barnes, Robert Marc Lebel, Krishna
Nayak, Accelerated DCE MRI Using Constrained Reconstruction Based on
PharmacoKinetic model dictionaries, ISMRM 2015, p0196
15. Weiyi Chen, Yi Guo, Ziyue Wu, Krishna S. Nayak, MRI Constrained
Reconstruction without Tuning Parameters Using ADMM and Morozov's
Discrepency Principle, ISMRM 2015, p3410
16. R. Marc Lebel, Yi Guo, Yinghua Zhu, Sajan Goud Lingala, Richard Frayne, Linda
B Andersen, Jacob Easaw, Krishna S Nayak, The comprehensive contrastenhanced
neuro exam, ISMRM 2015, p3705
17. Yinghua Zhu, Yi Guo, R. Marc Lebel, Meng Law, Krishna S Nayak, Randomized
Golden Ratio Sampling For Highly Accelerated Dynamic Imaging, ISMRM 2014,
p4365
18. Yi Guo, Xiaoke Wang, Sheng Fang, Kui Ying, Hua Guo, Shi Wang, Enhanced
motion correction combining PROPELLER and Parallel Imaging, ISMRM 2012,
p4862
Patents
1. Krishna S. Nayak, Yi Guo, R. Marc Lebel, Sajan Goud Lingala, Yinghua Zhu,
“Method for Improved Dynamic Contrast Enhanced Imaging Using TracerKinetic
Models as Constraints”, United States Patent Application number: 62/336/033, filed
May 2016.
2. Yinghua Zhu, Yi Guo, Krishna S. Nayak, Robert Marc Lebel, “DYNAMIC 3D MRI
DATA SAMPLING” United States Patent Application 15/075,716, filed March,
2016
ix
Abbreviation
ADMM Alternating Direction Methods of Multipliers
AIF Arterial Input Function
B1 Radiofrequency magnetic field strength
B1
+
Transmit radiofrequency magnetic field
BBB Blood Brain Barrier
CT Computed Tomography
DCE Dynamic ContrastEnhanced
DICOM Digital Imaging and Communications in Medicine
eTofts ExtendedTofts modelling
FSPGR Fast spoiled gradient echo sequence
Kep Rate constant from interstitium to plasma
K
trans
Volume transfer constant from plasma to interstitium
lBFGS Limited memory Broyden–Fletcher–Goldfarb–Shannon
algorithm
MRI Magnetic Resonance Imaging
mSSIM Mean Structural Similarity Index
NLCG Nonlinear conjugate gradient algorithm
PACS A picture archiving and communication system
Patlak Linearized PK modelling
PET Positron Emission Tomography
TK TracerKinetic
R1 1/T1, Longitudinal relaxation rate
r1 Contrast agent relaxivity
rMSE Root Mean Squared Error
ROI Region of Interest
SPECT Singlephoton emission computed tomography
T1 Longitudinal relaxation time constant
ve Volume fraction of extravascular extracellular space
vp Volume fraction of blood plasma
x
Abstract
Dynamic Contrast Enhanced Magnetic Resonance Imaging (DCEMRI) is a non
invasive technique that provides information about the delivery of blood and agents within
the blood to bodily organs. In oncology imaging, DCEMRI is used to generate tracer
kinetic (TK) parameter maps to evaluate tumor severity, differentiate malignant from
benign, and to evaluate response to therapy. Current DCEMRI methods have poor
reproducibility, low spatial resolution, and limited spatial coverage.
The overall goal of my dissertation is to improve brain DCEMRI by using constrained
reconstruction and modelbased reconstruction from undersampled raw data (sampled
below the Nyquist rate). Such techniques can provide substantially higher spatiotemporal
resolution, better coverage, and improved image quality with the same scan time and
contrast agent dose. Using modelbased reconstruction, TK parameter maps, the endpoint
product of DCEMRI, can be directly estimated from undersampled data with improved
quality and reproducibility.
I evaluated a specially tailored constrained reconstruction technique for DCEMRI in
brain tumor patients. The proposed technique is able to provide high spatial resolution and
wholebrain coverage via undersampling and constrained reconstruction with
multiple sparsity constraints. Conventional fullysampled DCEMRI and the proposed
wholebrain DCEMRI were performed on the same brain tumor patients for evaluation.
The proposed method is able to provide wholebrain highresolution DCEMRI with
improved image quality compared to conventional DCEMRI, using only 1/30
th
of the raw
data. These advantages may allow comprehensive permeability mapping in the brain,
which is especially valuable in the setting of large lesions or multiple lesions spread
throughout the brain.
I also developed a modelbased direct reconstruction technique for DCEMRI, where
the TK parameter maps can be directly estimated from rawdata. By comparison with a
stateoftheart indirect constrained reconstruction technique, the proposed direct approach
provides improved TK map fidelity, and enables much higher acceleration up to 100× .
With the prospective study, this method is shown to be clinically feasible and provide high
quality wholebrain TK maps.
To address the limitation of the direct reconstruction, I proposed a model consistency
constrained reconstruction for DCEMRI. The proposed reconstruction poses the tracer
kinetic (TK) model as a model consistency constraint, enabling the inclusion of different
TK solvers and the joint estimation of the arterial input function (AIF) from highly under
xi
sampled data. Good quality TK maps and patient specific AIF can be reconstructed at high
undersampling rates up to 100× .
1
Chapter 1. Introduction
Motivation
Cancer is the second leading cause of death in the United States, accounting for
approximately 600,000 deaths in 2013 (1). Brain cancer refers to cancerous tumors that
form inside human brain tissue. Cancerous tumors can be divided into primary tumors that
start within the brain, and secondary tumors that have spread from somewhere else, also
known as brain metastases. Brain cancer has a high mortality rate, depending on tumor
type, tumor size, and how early the tumor is detected. The mortality rate of brain tumor
can be greatly reduced if small lesions can be detected and treated early.
A variety of examinations are available for brain tumor evaluation and diagnosis, among
which Magnetic Resonance Imaging (MRI) and Computed Tomography (CT) are the most
commonly used (2). In clinical practice, multiple modalities are used jointly to accurately
diagnose the tumor. Reproduced from (3), Figure 1.1 shows representative CT and MRI
images.
2
Figure 1.1 Multimodality medical imaging scan of the brain of a middleaged man with a
primary central neural system lymphoma (A) A CT scan shows a spontaneous mild
hyperdense lesion. (B) A fluid attenuated inversion recovery (FLAIR)weighted MRI of a
peripheral oedema. (C and D) A T1weighted spinecho MRI without (C) and with (D)
contrast with homogeneous enhancing of two lesions in contact with the ventricles.
Contrastenhanced T1w MRI provide high detectability of the tumor (see arrows) (E) A
diffusionweighted MRI with a high signal of the lesion. (F) An apparent diffusion
coefficient map of a diffusionweighted MRI that shows a restriction of intratumoral
diffusion (see arrow), which is suggestive of a high cellular tumor. Reproduced from figure
7 in Ricard, D., Idbaih, A., Ducray, F., Lahutte, M., HoangXuan, K., & Delattre, J. Y.
(2012). Primary brain tumours in adults. The Lancet, 379(9830), 1984–1996.
MRI is a noninvasive imaging tool that provides excellent softtissue contrast, arbitrary
reformatting, and 3D imaging. Comparing to CT or PET that expose patients to radiation,
MRI is safe and noninvasive, making it ideal to longitudinally monitor tumor progression
and treatment effectiveness.
As shown in Figure 1.1, MRI provides multiple different contrasts. For brain tumor
evaluation, efforts have been made to standardize the imaging protocol in clinical trials,
and in (4), the consensus recommends the following key elements for a brain tumor
3
imaging protocol: (i) a precontrast, 3dimensional, isotropic, inversionrecovery T1
weighted gradient echo (IRGRE) sequence; (ii) an axial, 2dimensional T2weighted fluid
attenuated inversion recovery (FLAIR) sequence obtained using a turbospinecho (TSE)
readout; (iii) an axial, 2dimensional, 3directional (isotropic) diffusionweighted imaging
(DWI) sequence obtained using echoplanar (EPI) or radial acquisition; (iv) an axial, 2
dimensional T2weighted TSE sequence; and (v) a postcontrast, 3dimensional, isotropic,
T1weighted IRGRE sequence with matching acquisition parameters to precontrast T1
weighted images.
The key element for a brain tumor imaging protocol is the pre and postcontrast T1
weighted images obtained from a contrast injection. Apart from the aforementioned
sequences for clinical standardization, emerging MRI techniques like dynamic contrast
enhanced (DCE) or dynamic susceptibility contrast (DSC) perfusion MRI have been
proposed to track the contrast agent dynamics during the contrast administration.
Comparing to static T1 or T2weighted anatomic images, these dynamic imaging
techniques can be used for quantitative assessment of specific pathophysiologic parameters,
more accurate grading of intracranial tumors, and differentiation of tumors from normal
tissue (5).
Dynamic Contrast Enhanced Magnetic Resonance Imaging (DCEMRI)
Current brain tumor evaluation criteria, Response Assessment in NeuroOncology Criteria
(RANO) (6, 7), can only consider lesions that are larger than 1cm in maximal diameter.
This is due to the poor resolution of conventional MRI techniques. DCEMRI is not part
of the RANO criteria yet, but it is increasingly of interest. It can provide important
information about neurovascular system, so it may be a powerful tool for evaluating
effectiveness of therapies.
DCEMRI is a noninvasive imaging tool that tracks the delivery of blood to plasma,
capillary, and bodily organs by injecting an exogenous contrast agent. In T1 weighted DCE
MRI, a Gadoliniumbased contrast agent will shorten the T1 value, a key parameter in MRI
imaging that will result in change of signal intensity. Thus the delivery of blood can be
traced by the brightened signal intensity, and lesions and abnormal regions can be labeled.
As shown in Figure 1.2, the blood–brain barrier (BBB) separates the parenchyma of the
central nervous system from the blood. By tracking the contrast agent dynamics from blood
to abnormal tissues through BBB break down, DCEMRI is able to assess many brain
pathologies that cause an opening of the BBB, such as tumors, multiple sclerosis and acute
ischemic strokes. While these diseases show relatively large abnormalities in BBB
functionality, there is also growing interest in the application of DCEMRI to pathologies
4
associated with more subtle and chronic BBB disruption, such as cerebral small vessel
disease, diabetes and Alzheimer’s disease (8).
Figure 1.2 Schematic drawing of the neurovascular unit. The BBB is formed by endothelial
cells that line brain capillaries and are sealed by tight junctions. Astrocytes, pericytes,
microglial cells and basement membranes interact with the endothelium of the BBB,
providing functional and structural support. Reproduced from A. K. Heye, et al.,
“Assessment of blood–brain barrier disruption using dynamic contrastenhanced MRI. A
systematic review,” NeuroImage: Clinical, vol. 6, pp. 262–274, 2014.
From the dynamic signal changes in DCEMRI, trackerkinetic (TK, or termed
pharmacokinetic in previous literatures) parameter maps (K
trans
, Kep and vp etc.) can be
derived. K
trans
is an important TK parameter that indicates BBB permeability surface and
blood flow. It can be used for direct evaluation of the tumor severity. These enable DCE
MRI to be a powerful tool to characterize tumor biology and treatment response (9–11).
However, DCEMRI currently is not the standard of care in many centers conducting
clinical trials in oncology (12). This is because 1) low resolution and limited coverage of
the dynamic scans cannot provide adequate information for the wholebrain; 2) accuracy
and reproducibility of quantitative DCEMRI kinetic modelling is difficult to trust and
evaluate due to different TK modelling strategies and reconstruction techniques. The
clinical translation of DCEMRI is considered a major challenge. Improvement of
resolution, coverage, and modeling reproducibility can greatly help clinical translation of
DCEMRI, and is the main focus of my research.
5
Organization of the Dissertation
• Chapter 2 provides background and review of related technology that is helpful to
understand later chapters. It gives highlevel description of MRI physics, sampling
schemes and reconstruction techniques that are related to the author’s work. It also
provides detailed description of a DCEMRI pipeline including sequence, sampling,
preprocessing, correction, reconstruction and kinetic modeling.
• Chapter 3 describes a prospective evaluation of brain DCEMRI with constrained
reconstruction. A novel constrained reconstruction with multiple constraints
specifically tailored to DCE was applied to brain tumor patients to achieve a whole
brain highresolution DCEMRI. Double injection on the same patients allowed us
to compare directly to the current best practices. In fifteen brain tumor patients that
were scanned, this new approach has shown the potential to vastly improve
visualization and characterization of brain lesions with DCEMRI.
• Chapter 4 describes a novel and efficient reconstruction scheme to directly estimate
tracerkinetic parameter maps from highly undersampled DCEMRI data. By
comparison with a stateoftheart indirect compressed sensing method, it
demonstrates that the proposed direct approach provides improved TK map fidelity,
and enables much higher acceleration. With the prospective study, this method is
shown to be clinically feasible and provide highquality wholebrain TK maps.
• Chapter 5 describes a joint AIF and TK estimation framework using a model
consistency constraint. By posing the TK model as a model consistency constraint
instead of directly forcing the model. The framework can allow for model deviation
and provide better quality TK maps. It decouples the direct reconstruction problem
into two welldefined subproblems, thus greatly reducing the complexity of the
algorithm. This formulation also has the benefit of easy inclusion of different TK
solver, and joint estimation of AIF from undersampled data.
6
Chapter 2. Background
MRI Physics
This section will introduce basic MRI physics on MR signal generation, T1 and T2
relaxation, and concept of kspace and (k,t)space. It will also introduce Parallel imaging,
one of the acceleration techniques that is routinely used in MRI scanner now.
2.1.1 MR Signal Generation
Because of its presence in water and fat in human body, MRI typically utilizes hydrogen
protons for imaging. The hydrogen proton naturally has a magnetic moment, or spin, with
its positive charge and angular momentum. Naturally the spins have random directions,
and will not have any net magnetization in macroscopic scale as shown in Figure 2.1 (a).
In a MRI system, a strong constant magnetic field (3 Tesla for the commonly used
commercial machines) is present in the bore of the machine. With this external magnetic
force (denoted as B0), the hydrogen spins will accumulate a net magnetization in the
direction (Mz) of the external magnetic force, as shown in Figure 2.1 (b). The angular
frequency of nuclear precession will be:
00
B
(2.1)
This is also known as Larmor frequency. ϒ is a constant that is nucleusdependent
(42.58 MHz/T for Hydrogen).
If a transverse oscillating magnetic field (denoted as B1) is introduced with the same
frequency as Larmor frequency, it will tip down the magnetization by the angle determined
by the strength and duration of B1, as shown in Figure 2.1 (c) and (d). This angle is called
flip angle (FA), this frequency is called resonance frequency, and this oscillating magnetic
field is called RF (radiofrequency) pulse, since B1 oscillates in the radiofrequency range.
After RF excitation, the magnetization rotates in xyplane. The vector component of the
magnetization in this plane is called transverse magnetization (denoted as M xy). This
rotating transverse magnetization will then give rise to a signal in the receiver coil.
7
Figure 2.1. With no external magnetic field present, spins rotate about their axes in random
direction (a). In the presence of a magnetic field, slightly more spins align parallel to the
main magnetic field, B0, and thus produce longitudinal magnetization, Mz (b). An RF pulse
(c) tips the magnetization vector by exactly 90° , causing the entire longitudinal
magnetization to flip over and rotate into transverse magnetization, Mxy (d). Reproduced
from (13).
2.1.2 T1 and T2 Relaxation
The MR signal will rapidly fades and return to stable state due to spinlattice interaction
and spinspin interaction. The spinlattice interaction is called T1 relaxation, and spinspin
interaction is called T2 relaxation.
8
Figure 2.2 T1 and T2 relaxation over time for different tissue type. T1 relaxation will restore
the longitudinal magnetization (a), and T2 relaxation will decrease the transverse
magnetization (b).
T1 relaxation will cause the magnetization to return to stable state, or equilibrium in Mz
direction, as shown in Figure 2.2 (a). T2 relaxation will decrease the transverse
magnetization, and eventually destroy the MR signal, as shown in Figure 2.2 (b). Different
tissue types will have different T1 and T2 relaxation rate, and utilization this property can
create MR images of different contrast between tissues. T1 value of a specific tissue is
defined as the time required to restore 63% of the original longitudinal magnetization (M0),
and T2 value is defined as the time required to lose 63% of the transverse magnetization
(Mxy). Typical T1 and T2 values of different tissues at 3T in literatures (14–16) are listed in
Table 2.1.
Table 2.1 Typical T1 and T2 values of different tissues at 3T.
Tissue T1 value at 3T (ms) T2 value at 3T (ms)
White Matter 1110 ± 45 69 ± 3
Gray Matter 1470 ± 50 99 ± 7
Blood 1550 275 ± 50
Skeletal Muscle 1420 ± 38 50 ± 4
Cerebral Spinal Fluid (CSF) 3120 150
If the contrast of a MRI image is weighted primarily by the signal difference caused by
T1 difference, it is called T1weighted imaging. Similarly, T2weightied imaging is for
9
contrast primarily caused by T2 difference. Typical T1 and T2weighted images of the brain
is shown in Figure 2.3. In T1weighted images, tissues (e.g. CSF) with longer T1 will be
darker, since the signal is not recovered yet. In T2weighted images, tissues with longer T2
will be brighter, since the signal is still preserved.
Figure 2.3 Typical T1 and T2 weighted images of the brain. In T1weighted image, gray
matter and CSF are darker because of long T1 values. In T2 weighted image, gray matter
and CSF are brighter because of long T2 values, in contrast to T1 weighted image.
2.1.3 MRI Signal Space
In additional to the constant magnetic field B0 that is used to create net magnetization,
and oscillating magnetic field B1 that is to excite the signal, another gradient magnetic field
that gradually (usually linearly) changes the magnetic field strength based on spatial
location is applied to the magnetic core. As a result, the MR signal is coded with spatial
information, as shown in the equation below (3D case):
2 ( )
( , , ) ( , , ) d d d
x y z
i k x k y k z
x y z
S k k k I x y z e x y z
(2.2)
, where I is the 3d image being scanned with spatial coordinate x, y and z, kx, ky and kz are
determined by the gradient strength and localization, and S(kx,ky,kz) is the acquired MR
signal. Depending on how the gradient is applied for localization, kx direction is usually
called frequency encoding direction, and ky, kz direction is called phase encoding direction.
In this form, the signal is also in the Fourier transform domain of the image, therefore the
MR signal space is called kspace. Normally we use Nyquist rate to sample the kspace at
10
Cartesian grid, and 3D inverse Fourier transform will reconstruct the images from kspace,
as shown in Figure 2.4 for one slice in 3D. In the context of dynamic imaging, where a
whole (or partial) kspace is acquired at some time interval, the kspace plus time dimension
is also called (k, t)space.
Figure 2.4 Relationship between image and Cartesian Nyquist sampled kspace. In kspace,
Δkx will determine the fieldofview (FOV) of the image, and kx
max
will determine the
resolution of the image. Reproduced using images from Karla Miller, FMRIB, University
of Oxford.
2.1.4 Parallel Imaging
MRI is intrinsically slow using Cartesian Nyquist sampling, where we need to acquire
every point in kspace to avoid aliasing in the images. However, if we can utilize multiple
coils that are receiving the same signal with different coil sensitivity (weighting of the
signal), we are able to restore images by acquiring only part of the kspace. This technique
is called parallel imaging.
As shown in Figure 2.5, if multiple coils are used to acquire only half of the kspace by
sample every other line of kspace, the resulting individual coil images will be wrapped
around. However, by utilizing the different coil weighting (coil sensitivity), the original
image can be reconstructed from these multiple coil images. This image space technique is
called SENSE (Sensitivity encoding) (17). It targets to solve the image x from under
sampled kspace y using the following equation:
1 1 1 1
()
HH
x E E E F y
(2.3)
11
, where E is the coil sensitivity encoding maps, F is Fourier transform, and Ψ is the noise
covariance matrix.
Besides SENSE, other variants of parallel imaging exist. The coil sensitivity
information can be implicitly used in kspace directly to reconstruct fullysampled kspace
from undersampled kspace, and such techniques include GRAPPA (18), SPIRiT (19),
ESPIRiT (20), LORAKS (21) , etc.
Figure 2.5 In parallel imaging, an array of receiver coils simultaneously collects the MR
signals. Scan time is shortened by reducing the number of phaseencoding steps. As a result,
the individual images are obtained with a smaller FOV and show the typical wraparound
artifacts. A complete image without wraparound artifacts is reconstructed by combining
the individual images. Reproduced from Weishaupt D, Kö chli VD, Marincek B: How Does
MRI Work? 2008.
DCEMRI Basics
This thesis focuses on dynamic contrast enhanced MRI (DCEMRI). DCEMRI is one
of the dynamic MRI techniques that employs continuous T1weighted imaging before,
during and after a bolus injection of a gadoliniumbased contrast agent (GBCA). With the
help of additional mapping sequences before DCEMRI scan, changes in GBCA
concentration is derived from changes in signal intensity, then regressed to quantify tracer
kinetic (TK) parameters such as K
trans
(volume transfer constant), vp (fractional plasma
volume), and ve (fractional extravascular extracellular space volume) (9, 22) . Figure 2.6
from (8) illustrated typical T1 weighted DCEMRI images, the signal enhancement of
tumor tissue due to contrast injection, and the derived K
trans
maps.
12
Figure 2.6 Illustration of DCEMRI. The repeated acquisition of T1weighted images after
contrast agent injection allows the calculation of signal enhancement as a function of time
(middle) when compared to the precontrast signal intensity (left). These curves can be
used to calculate maps of quantitative pharmacokinetic parameters (e.g. K
trans
, right).
Reproduced from A. K. Heye, et al., “Assessment of blood–brain barrier disruption using
dynamic contrastenhanced MRI. A systematic review,” NeuroImage: Clinical, vol. 6, pp.
262–274, 2014.
DCEMRI is used for quantitative assessment of brain tumors (23, 24), multiple
sclerosis lesions (25), and Alzheimer’s disease (26) and other neurological disorders that
involve bloodbrain barrier (BBB) disruption. This section will introduce each step of a
complete DCEMRI procedure.
2.2.1 Calibration Scans
B1+ Mapping
B1+ represents the inhomogeneity in the transmitted RF filed, which will affect the
actual flip angle applied to the MRI bore. Flip angle will in turn affect the conversion of
DCEMRI signal intensity to contrast concentration curves, and impact the final TK
parameter value accuracy. Figure 2.7 shows typical B1+ maps in head phantom and human
brain. It shows that the actual flip angle inside the MRI bore is often different from what
we prescribe on the scanner console. Ignoring this variation will lead to inaccurate TK
mapping results. (12)
13
Figure 2.7 Illustration of a typical B1+ maps of a head phantom (left) and human brain
(right) in normalized unit, where 1 means the actual flip angle is equal to the prescribed
flip angle. The inhomogeneity of the B1+ maps can be clearly seen.
The clinically available 3D doubleangle gradientecho method (DAM) (27) is usually
performed to measure B1+ maps. This method acquires two gradientecho images (I1 and
I2) with prescribed flip angles of α and 2α. The actual flip angle (hence the B1+ maps) can
be calculated as:
2
1
arccos
2
I
I
(2.4)
Apart from DAM, a faster and more efficient approach called BlochSiegert shift (28)
was also proposed recently. This is also the standard B1+ mapping protocol on GE MRI
scanners.
T1 Mapping
Precontrast T1 and M0 values are also important parameters that affect the accuracy in
conversion of signal to contrast concentration curves. Figure 2.8 shows a typical M0 and
T1 maps of human brain. Mapping the variation of M0 and T1 maps is important for accurate
TK parameter mapping.
14
Figure 2.8 Typical M0 (left) and T1 (right) maps of human brain. M0 map is in normalized
unit, and T1 map unit is millisecond. T1 values are greatly different for different type of
tissues in the brain.
Various T1 mapping techniques exist in MRI, and in DCEMRI, DESPOT1(29) method
is commonly used for the speed and efficiency. DESPOT1 method acquires multiple
spoiled gradient echo (SPGR) images (S) with different flip angles (α), and fits the T1 and
M0 maps using following signal equation:
1
1
/
0
/
sin (1 )
()
1 cos
TR T
TR T
Me
S
e
(2.5)
, where TR is the repetition time in MRI. The fitting is performed by converting the
problem to a linear matrix inversion problem to solve for e
TR/T1
using the acquired multiple
images S(α), then T1 can be computed easily from e
TR/T1
. Figure 2.9 shows SPGR signal
intensity versus flip angle at different T1 value. Small flip angles of 2⁰, 5⁰, 10⁰ are often
used for T1 mapping utilizing the linear region of the fitting equation.
15
Figure 2.9 SPGR signal intensity versus different flip angle at different T1 value. Based on
these curves, small flip angles of 2⁰, 5⁰, 10⁰ are close to linear region, and are used for the
fitting the T1 values.
2.2.2 DCEMRI Scan
Dynamic contrast enhanced MRI uses a fast spoiled gradient echo sequence (SPGR) to
acquire a series of images during the administration of a Gadoliniumbased contrast agent.
Taking account of the intensity of the signal and the contrast between different tissues, a
flip angle of 10⁰ to 30⁰ is appropriate to get images with good SNR and tissue contrast. In
practice, a flip angle of 15⁰ is usually used for DCEMRI.
Conventionally, data is sampled at Nyquist rate on the Cartesian grid for proper
restoration of the MRI images. SPGR signal intensity follows the same Equation (2.5), and
for dynamic images (S), the signal equation is:
()
0
()
sin (1 )
()
1 cos
TR R t
TR R t
Me
St
e
(2.6) (2.6)
where R(t) is R1 (reciprocal of T1) over time. This equation assumes steasdy state for the
spoiled gradient echo sequence, and not applicable to inflow enhancement. The time
dimension t is neglected in the following text for simplicity.
16
2.2.3 Preprocessing for TK modelling
Tracerkinetic (TK) modelling in DCEMRI mainly involves fitting the TK parameters
using the contrast concentration in tissue (C) and the contrast concentration in artery (also
called artery input function, or AIF). This section describes the method to estimate these
two variables from the acquired data.
Signal to Concentration Conversion
MRI samples in kspace, which is Fourier transform space of the images. Basic MRI
reconstruction involves data correction, inverse Fourier transform, and coil combination
(if multiple coils are used). If multiple coils are used and some undersampling are involved,
we often solve the target images (S) by the following leastsquare problem:
2
2
arg min  ( ) 
S
S UFE S y
(2.7)
where U is the undersampling matrix, F is the multidimensional Fourier transform, E
is the coil sensitivity encoding, and y is the acquired kspace data.
This is one of the methods to solve SENSE (17, 30), a parallel imaging technique
introduced in section 2.1.4. Iterative algorithms like Conjugate Gradient (CG) algorithm
can be used to solve this leastsquare problem.
After we acquire the dynamic images (S), which is in 4d (3d spatial plus time) in our
case, we can get the image differences by subtracting the precontrast first image (S0):
0
S S S (2.8)
Then the contrast concentration (C) is converted from ΔS based on the SPGR signal
equation (2.6) (31):
0
0
0
01
1
1
sin 1 cos
1
ln ,
1
1 cos
sin 1 cos
(  ) /
TR R
Sm
Mm
R m e
TR
Sm
Mm
C R R r
(2.9)
where R is the R1 (reciprocal of T1) values over time, TR is the repetition time, α is the
flip angle, r1 is the contrast agent relaxivity. R0 and M0 are the precontrast R1 and the
equilibrium longitudinal magnetization that are estimated from a T1 mapping sequence as
described in section 2.2.1.
17
AIF Extraction
Before the TK modelling, we need to extract the arterial input function (AIF), which is
an essential parameter in the modelling part. AIF is the concentration curve in an artery,
which can be determined manually by a radiologist. This ROI can also be determined
automatically by clustering or modelbased techniques (32, 33). This is called patient
specific AIF, or patAIF. Temporal resolution is critical in obtaining a highquality AIF,
and a temporal resolution of one second is preferred, which is difficult to achieve using
conventional MRI techniques (34, 35).
In the literature, a fixed populationaveraged AIF (popAIF) was also proposed (36).
PopAIF is used when the papAIF is difficult to obtain from the data. It is robust to noise
and interpatient variability, but may lead to inaccurate TK values since patientspecific
information is not taken account for. PopAIF takes a parameterized function form (11
parameters in Parker’s popAIF model), with a few parameters adjusted by the actual
patient data. Other parameterized AIF was also proposed, with the delay and dispersion
parameter to be adjusted by the actual data (37). Figure 2.10 shows a default popAIF as
proposed in (36), and two deformed shapes by adjusting the parameters.
Figure 2.10 Shapes of populationaveraged AIF using different parameters. “Default”
uses the parameters provided in the literature.
Figure 2.11 illustrates how to extract patAIF from a ROI in patientdata, and the
resulting TK parameter maps. It shows that using popAIF may lead to significantly
different TK maps.
18
Figure 2.11 Illustration of extracting patAIF from a manually selected ROI, and
comparison to popAIF (left). Patlak model parameter (K
trans
and vp) from these two
different AIFs (right).
2.2.4 Trackerkinetic (TK) Modeling
After contrast agent concentration (C) are estimated from the data, TK parameter maps
can be derived from these contrast concentration curves. C(t) and Cp(t) are used to denote
concentration in tissue and artery (i.e. AIF) in this section. TK parameter maps are the end
point product of DCEMRI scan, these maps show important physiological information,
and are of more importance to radiologists for diagnosis and treatment response. For
example, the K
trans
value, which is the forward transfer rate constant from plasma to
intersititium, is correlated to the BBB permeability. This is considered a direct parameter
to show tumor grade and severity (22, 23, 38). The plasma volume, vp, is related to the
blood volume in the plasma, and shows the normal or abnormal blood flow in patients (9,
22). Kep is the interstitiumtoplasma rate constant, and is related to the reabsorption of
contrast agent.
For the TK parameter modelling, different models exist to convert the concentration
curves to TK parameters. Each of these models assume some approximation about the
perfusion and flow condition. Table 2.2 lists 4 representative TK models of different
complexity and their corresponding parameters and assumptions.
19
Table 2.2 List of 4 representative TK models and conditions
Model Parameters Condition
0. Null 0 insufficient vascular filling
1. Steadystate model vp
vascular filling with no
microvascular leakage
2. Patlak model K
trans
, vp
leakage without vascular
reabsorption
3. extended Tofts model K
trans
, K ep, vp leakage with reabsorption
In the brain tumor region, it is shown that three parameters of the extendedmodel are
sufficient to fit dynamic contrast enhanced data for glioblastoma. And in slow leakage
cases, two parameters of Patlak model provide reasonable approximation (39). The section
below will introduce these two widely used TK models.
2.2.4.1 Common TK Models for Brain Tumor Evaluation
ExtendedTofts
A widely used TK model for brain tumor is the extendedTofts model (eTofts) (9). This
model is valid in highly perfused tissues and in weakly vascularized tissues with a well
mixed interstitium (40). The model is defined as the following fitting equation:
()
0
( ) ( ) ( )
ep
t
K t u
trans
p p p
C t K C u e du v C t
(2.10)
The model fits three parameters (K
trans
, Kep, vp). Though accuracy is sacrificed
comparing to more general models, eTofts still provides meaningful TK values with high
precision, and is useful in many clinical applications (40).
Patlak
The Patlak model is a linear simplification of eTofts mode, where backflux from
interstitium is neglected. The model is defined as:
0
( ) ( ) ( )
t
trans
p p p
C t K C u du v C t
(2.11)
Although Patlak assumes many extreme approximations, it is still widely used in many
centers because it is linear and robust.
20
2.2.4.2 Model Selection
For this variety of TK models, no one is considered “correct” and goldstandard for
DCEMRI. Like quote from George E. P. Box, “all models are wrong, but some are useful”.
The models with more parameters to be fit like eTofts will provide better fitting for the
data, but increase the variability of the estimated TK parameters (40–42). There lacks
ground truth for TK parameters in DCEMRI, and the injection of contrast agent to human
body makes it difficult to test the precision and reproducibility of different modeling and
imaging strategies. Many studies used simulated data to test the fitting performance of each
model and reconstruction techniques, but this provides limited information about TK
modelling in real patient data. If the fullysampled data is available, a Ftest based on the
fitting errors of each model can be used to determine which model to use (39). However,
this method is difficult to implement in undersampling reconstruction scenario.
Sparse Sampling and Constrained Reconstruction
MRI imaging needs multiple phase encoding to encode the spatial information of the
MR signal. This phase encoding step takes at least 35 ms each time, making a high
resolution 3D image up to 12 minute to be acquired. For DCEMRI, to capture the
dynamic changes during a contrast injection, tradeoffs have to be made between spatial
resolution, spatial coverage, and temporal resolution. Therefore, many techniques have
been proposed to accelerate DCEMRI scans. Most of these techniques focus on
reconstructing proper images from undersampled data. This section will review some of
these techniques that are relevant to the author’s main contribution.
Parallel imaging is a MRI technique that utilizes the redundancy in coil dimension to
accelerate MRI scans. These techniques include SMASH (43), SENSE (17), GRAPPA (18),
ESPIRIT (20) etc. As introduced in section 2.1.4 and 2.2.3, SENSE is already routinely
used in clinical DCEMRI to achieve 2x acceleration.
Compressed sensing (CS) is another technique that aims to reconstruct images from
sparselyundersampled data (44). CS is particularly well suited for dynamic imaging,
which can exploit the redundancy of information in the temporal dimension, either through
dictionarylearning (45, 46), or highpass filtering (47, 48). A combination of PI and CS
has been shown to greatly accelerate the data acquisition, while achieving significantly
higher spatiotemporal resolution and large spatial coverage with only slight image quality
penalties (19, 20, 49).
Several groups have employed PICS techniques to improve DCE imaging. Zhang et
al. (50) employed a locally lowrank constraint in combination with parallel imaging to
achieve up to 19× acceleration rate in pediatric patients. Wang et al. (51) used a reference
imagebased compressed sensing and achieved acceleration factor of 10× without
21
degrading spatial resolution. Feng et al. (49) used compressed sensing, parallel imaging
and goldenangle radial sampling to achieve fast and flexible DCEMRI. Rosenkrantz et
al. (52) examined 20 prostate cancer patients using a similar scheme to evaluate the results
from constrained reconstruction against conventional DCEMRI.
This section will give a general introduction on sparse sampling and constrained
reconstruction for dynamic MRI imaging.
2.3.1 Sparse Sampling
Compressed sensing theory relies on random sampling of kspace, and a variety of
sampling patterns have been proposed. Various nonCartesian sampling schemes have
been proposed over the years, and the literatures have shown that nonCartesian sampling
have multiple benefits over Cartesian sampling in terms of sampling efficiency, robustness
to motion and suppression of offresonance effects, etc. (30, 53, 54) However, non
Cartesian sampling is also expensive to implement on the MRI machine, and takes
additional computational steps and effort for image reconstruction. In the context of 3D
imaging, which is used in DCEMRI (3D plus time), it is easy to sample the Cartesian k
space grid following a nonCartesian style pattern. This pseudoCartesian sampling can
benefit from both Cartesian and nonCartesian sampling. Here we will focus on these
pseudoCartesian sampling patterns for 3D imaging.
Figure 2.12 illustrates three representative pseudoCartesian sampling patterns:
PoissonDisc (PD), GoldenAngle Radial (GA), and Randomized GoldenAngle Radial
(RGA). The undersampling rate (comparing to Nyquist fully sampling) is 30 folds (R=30x).
PD is a randomized sampling pattern that is initially proposed for compressed sensing,
and has the benefit of creating incoherent artifacts that is easily removable by l1 norm
constraint (55). However, PD is initially proposed for static imaging, and in dynamic 3D
imaging, different patterns have to be generated beforehand for each time frame with a
fixed temporal window, providing less flexibility for dynamic imaging (56, 57).
GA is derived from conventional radial sampling with a goldenangle of 111.246°
between each radial spokes to achieve an efficient coverage of kspace (53). When used in
dynamic 3D MRI, GA can be used to continuously sample the (k,t)space, and a group of
spokes can be binned retrospectively to achieve arbitrary temporal resolution. Because of
the property of goldenangle, each time frame will still have optimal kspace coverage.
RGA is a newly proposed sampling pattern that is tailored for dynamic 3D MRI, where
both the benefits of PD and GA is inherited (57). RGA is based on GA Cartesian sampling,
with random sampling of the kykz phase encode locations along each Cartesian radial spoke.
This sampling pattern will have the benefit of incoherent sampling from PD, and arbitrary
selection of temporal window from GA. Note that in the published literature, RGA is
22
acronymed as GOCART for “GOldenangle CArtesian Randomized Timeresolved 3D
MRI”, here we still use RGA for easy understanding of the connection to original GA
sampling.
Figure 2.12 Illustration of PoissonDisc (PD), GoldenAngle Radial (GA), and
Randomized GoldenAngle Radial (RGA) sampling patterns at R=30x.
In this work, we mainly used traditional (GA) and randomized goldenangle (RGA)
radial sampling patterns.
2.3.2 Constrained reconstruction
After undersampled kspace data is acquired, we will use certain prior knowledge as
constraints to restore the images from these undersampled data. The reconstruction
process mostly involves solving following optimization problem:
2
2
ˆ
arg min   ( )
S
S UFES y R S
(2.12)
Where the first term is called data consistency, and it enforces that the image to be
reconstructed (S) remains consistent with acquired kspace data (y). As introduced in
section 2.3.2.2 and 2.2.3, U is the undersampling matrix, F is the Fourier transform, and
E is the coil sensitivity encoding. The second term is the constraints based on some prior
knowledge of the images.
In the earlier literatures, l2 norm constraints are proposed first to improve SNR and
convergence speed. With the emergence of compressed sensing, l1 or lp norm that
approximate l0 norm condition have been proposed to enforce the sparsity of certain
23
transform domain of the images. The sections below will discuss different variants of
constrained reconstruction, and various algorithms to solve them.
2.3.2.1 l2 Norm Constraint
As introduced in section 2.1.4, in the initial publication on SENSE (17), it proposed to
solve the original problem using matrix inversion (or pseudoinverse). For simple linear
sampling pattern and small data size, it is the most intuitive and straightforward way to
solve the problem. However, if complicated sampling pattern is used, or if the data size is
large, commonly we target to solve a leastsquare problem using iterative algorithm (e.g.
conjugate gradient). And l2 norm is used to provide stable convergence and better
performance for a leastsquare problem. For example, Tikhonov regularization is a
commonly used l2 norm constraint (58):
22
2 0 2
ˆ argmin    ( ) 
x
x UFEx y L x x
(2.13)
, where L is a positive semidefinite linear transformation, x0 denotes the prior information
about the solution x. The regularization parameter λ determines the relative weights with
which these two estimates of error combine to form a cost function. In optimization
problem, an Lcurve, as illustrated in Figure 2.13, is often used to determine an appropriate
value of the regularization parameter that provides a good balance between the data
consistency and the prior knowledge.
24
Figure 2.13 An Lcurve illustrates the two costs during reconstruction of the aliased images
from an array. Using distinct regularization, the reconstruction biases toward minimizing
either the prior error or the model error. A tradeoff between these two error metrics is the
use of regularization at the “corner” of the Lcurve. Reproduced from Lin FH, Wang FN,
Ahlfors SP, Hä mä lä inen MS, Belliveau JW: Parallel MRI reconstruction using variance
partitioning regularization. Magn Reson Med 2007; 58:735–744.
For parallel imaging techniques, l2 norm is mainly helpful for convergence speed and
stability, to achieve higher undersampling rates or better reconstruction image quality,
more constraints that better utilize the prior information of the image have been proposed.
2.3.2.2 Compressed Sensing (l0 and l1 Constraint)
In additional to SENSE, Compressed Sensing (CS) is a newly proposed undersampling
reconstruction technique that allows for higher acceleration in MRI (44, 59).
As shown in Figure 2.14, MRI images will be sparse in certain transform domain, where
most of the signals are zeros. If we just use a portion of the largest values in these transform
domain, we can still restore reasonably good quality images. CS technique utilizes this
property of MRI image, and target to restore the images from vastly undersampled k
space, with a constraint that the image should be sparse in certain transform domain.
25
Figure 2.14 Transform sparsity of MR images. (a) Fully sampled images are mapped by a
sparsifying transform to a (b) transform domain; the several largest coefficients are
preserved while all others are set to zero; the transform is inverted forming a (c)
reconstructed image. Reproduced from (59). 59. Lustig M, Donoho DL, Santos JM, Pauly
JM: Compressed Sensing MRI. (March 2008):72–82.
Specifically, CS targets to reconstruct the image x from undersampled kspace data y
by solving the following constrained optimization problem:
minimize
0
  x
s.t.
2
  UFEx y
(2.14)
, where Ψ is some sparsifying transform, U is the undersampling matrix, F is Fourier
transform and E is the coil sensitivity encoding. The l0 norm is to enforce the sparsity of
the transform domain coefficient. However, directly solving a l0 norm constrained problem
can be challenging, and some approximation approaches have been proposed.
26
In the initial publication on CS (55), l1 norm is proposed to approximate l0 norm, and
the problem can be converted to an unconstrained form:
2
21
ˆ argmin    
u
x
x F Ex y x
(2.15)
Where λ is the regularization parameter, or constraint penalty. Many iterative algorithms
can be used to solve this optimization problem. For example, nonlinear conjugate gradient
algorithm is proposed to solve the problem by calculating the gradient of the cost function.
And for l1 norm, the gradient is calculated by relaxing the l1 norm to an approximating l2
norm.
Another important algorithm for solving constrained reconstruction problem is
augmented Lagrangian methods, where auxiliary variables are introduced to split the
original problem into easy subproblems. For Equation (2.15), auxiliary variables u and v
can be introduced to split the problem into:
2
21
u,
min     . ,
v
y UFu v s t v Tx u Ex
(2.16)
And convert into an unconstrained problem:
12
2 2 2
2 1 1 1 2 2 2 2
u, , , ,
min            
t
t
v S e e
y UFu v TS v e Ex u e
(2.17)
And during iteration, each variable can be solved alternatively while keeping other
variables constant. This method is also called alternating direction method of multipliers
(ADMM). Detailed algorithms how to solve each step can be found in Chapter 3.
27
Chapter 3. Clinical Evaluation of Constrained Reconstruction in Brain
DCEMRI
Introduction
T1weighted DCEMRI is a valuable albeit still evolving technique for mapping the
spatial distribution of brain vascular parameters such as perfusion, permeability, and blood
volume (8, 60). It employs serial T1weighted imaging before, during and after a bolus
injection of a gadoliniumbased contrast agent (GBCA). Changes in GBCA concentration
is derived from changes in signal intensity, then regressed to quantify tracker kinetic (TK)
parameters such as K
trans
(volume transfer constant), vp (fractional plasma volume), and ve
(fractional extravascular extracellular space volume) (9, 22). DCEMRI is used for
quantitative assessment of brain tumors(10, 24, 61), multiple sclerosis lesions (25), and
Alzheimer’s disease (26) and other neurological disorders that involve bloodbrain barrier
(BBB) disruption. DCEMRI is also used in clinical oncologic imaging for assessment of
breast (62) and prostate (63) cancer. In brain tumor evaluation, BBB permeability is
typically characterized by K
trans
(9). While its interpretation may be complex, with some
dependency on blood flow, K
trans
correlates with tumor severity and may be a useful
biomarker for monitoring therapeutic response and outcome (9, 25, 38, 64–66) .
Despite its usefulness, conventional clinical DCEMRI is limited by suboptimal image
acquisition that results in low spatial/temporal resolution and insufficient spatial coverage.
Low temporal resolution has also been linked to poor reproducibility of TK parameters
(23). A typical clinical DCE scan provides 230 s temporal resolution to detect signal
intensity changes resulting from contrast agent perfusion (8, 52). As a result, the inplane
voxel dimensions usually exceed 1 mm
2
and the slices are often greater than 5 mm thick.
The spatial coverage is often inadequate to cover the known pathology, such as in the
setting of multiple metastatic lesions.
Recently, compressed sensing (CS) theory has inspired a wide array of new data
acquisition and constrained reconstruction strategies that aim to reconstruct images from
sparselyundersampled data (44). CS is particularly well suited for dynamic imaging,
which can exploit the redundancy of information in the temporal dimension, either through
dictionarylearning (45, 46), or highpass filtering (47, 48). A combination of parallel
imaging (17) (PI) and CS has been shown to greatly accelerate the data acquisition, while
achieving significantly higher spatiotemporal resolution and large spatial coverage with
only slight image quality penalties (19, 20, 49).
Several groups, including ours, have employed PICS techniques to improve DCE
imaging. Zhang et al. (50) employed a locally lowrank constraint in combination with
28
parallel imaging to achieve up to 19× acceleration rate in pediatric patients. Wang et al.
(51) used a reference imagebased compressed sensing and achieved acceleration factor of
10× without degrading spatial resolution. Feng et al. (49) used compressed sensing, parallel
imaging and goldenangle radial sampling to achieve fast and flexible DCEMRI.
Rosenkrantz et al. (52) examined 20 prostate cancer patients using a similar scheme to
evaluate the results from constrained reconstruction against conventional DCEMRI.
Several groups have also utilized undersampling and constrained reconstruction
techniques to accelerate ContrastEnhanced (CE) Magnetic Resonance Angiography
(MRA). CEMRA is particularly amenable to this approach because subtraction
angiograms are sparse in the image domain. Barger et al. (67) used undersampled 3D
projection reconstruction trajectories and a “tornado” viewsharing scheme to achieve
isotropic resolution, broad coverage and 4s temporal resolution. The typical aliasing when
using undersampling is mitigated by the highcontrast properties in MRA. Haider et al.
(68) used a Cartesian radial technique in combination with 2D SENSE, partial Fourier and
view sharing to achieve 12 mm isotropic resolution and subsecond temporal resolution.
With emerging CS techniques, Trzasko et al. (69) demonstrated reduced noise and artifacts
level in Cartesian radial sampling MRA by utilizing a sparsitydriven nonconvex CS
method, and Lee et al.(70) achieved 1mm isotropic resolution, 1.1s frame rate
(corresponding to an acceleration factor of >100× ) with a CS based GraDeS algorithm. In
CEMRA, high spatiotemporal resolution and broad coverage is achieved by exploiting
high image contrast and a high degree of image domain sparsity.
Accelerating DCEMRI for the purpose of pharmacokinetic modeling is more
challenging than for CEMRA since the reconstruction is not spatially sparse and TK
modeling is performed based on every signal containing voxel. Furthermore,
reconstruction fidelity must be very high for accurate modeling whereas moderate error in
visually assessed angiograms is often tolerable. Finally, the dynamic range of contrast
induced signal change is smaller in tissue than in vessels, enabling subtle compression
artifacts to translate into noticeable errors in tissue parameters. The novelty of our proposed
approach is that multiple sparsity constraints are employed in different sparse transform
domains, each with low weight (56), mitigating biased artifacts producing from one heavy
constraint. We have been able to achieve the highest acceleration rate reported in the
literature to date, 36× , with excellent image quality. This enabled near isotropic voxel
dimensions with wholebrain coverage for DCEMRI.
Despite the promise of PICS methods, these techniques are poorly validated. Most
validation has been done by retrospectively discarding data from fully sampled data sets or
using phantom simulation. Both approaches provided ground truth, but are imperfect due
to unrealistic data acquisition or anatomical features. Prospectively undersampled studies
29
demonstrate the potential of the method, but lack ground truth images and have not been
well validated. In this work, we present the first, to the best of our knowledge, prospective
clinical evaluation of accelerated DCEMRI using constrained reconstruction in brain
tumor patients.
Materials and Methods
3.2.1 HighResolution WholeBrain DCEMRI
The experimental DCEMRI scan was based on a conventional Cartesian T1weighted
3D SPoiled fast GRadient echo (SPGR) sequence. The flip angle was 15⁰, TE was 2 ms,
and TR was 6 ms. The phase encoding order was altered to follow a Cartesiangrid golden
angle radial scheme (53, 71, 72), which acquired kykz phase encodes following golden
angle rotating radial spokes. The frequency encoding direction kx was fully sampled.
Cartesian SPGR scans with flip angles of 2⁰, 5⁰, and 10⁰, were performed sequentially for
DESPOT1 (73) T1 mapping prior to both the conventional and experimental DCE scans.
The Cartesian goldenangle radial sampling scheme provides incoherent kspace
sampling, even at very high acceleration rates, where only a few spokes are presented
within one time frame (49, 52). This approach provides comparable image quality as
Poissonellipse sampling, which has been used in similar l1based reconstruction (72). All
samples fall on a 3DFT Cartesian grid, therefore the Fast Fourier Transform operator can
be directly applied. We have implemented this sampling scheme on a clinical scanner (3T
GE Signa Excite HDx scanner), where the phase encode order table can be generated in
response to operator input. The implementation is straightforward, since no additional
modification of the existing gradient waveform is needed.
The undersampled raw data was reconstructed using a sparse SENSE reconstruction
scheme that utilizes multiple l1norm constraints with very low weights, as described in
Ref (56). These sparsity penalties were chosen based on expected spatial and temporal
characteristics of DCEMRI images. Reconstruction involved solving the minimization
problem in Equation (3.1), where the final image, x, remains consistent with the acquired
data, y, yet is sparse in the temporal finite difference (V) domain, the spatial ‘db2’ wavelet
domain (Ψ), and the spatial total variation domain (TV). The image is related to the acquired
data using known coil sensitivities encoding (E) and the undersampling Fourier transform
UF. Coil sensitivity maps were generated by computing a density compensated average of
all kspace data acquired at all DCE time points. This resulted in highSNR timeaveraged
3D kspace dataset that was fullysampled. Individual coil maps were then computed in the
30
standard way, by dividing each individual coil image by the root sumofsquares image
(17).
2
2 1 1 2 1 3 1
argmin        
t
t t t t t
S
S UFES y VS TVS S
(3.1)
This optimization problem is then solved by an efficient augmentedLagrangian method,
Alternating Direction Method of Multipliers (ADMM), which performs variable splitting
twice (74). This algorithm is one of many stateoftheart algorithms to solve these l1
constrained minimization problems (19, 44). This particular algorithm was chosen because
it provides fast convergence (74–76).
All sparsity transforms can be denoted by a tall matrix T, that is
,,
T
T V TV .
Dummy variables u, v are used to split Equation (3.1) to Equation (3.2), where is a long
vector that has
1 2 3
, and in the corresponding transform location.
2
21
u,
min     . ,
tt
v
y UFu v s t v TS u ES
(3.2)
Then Lagrangian method is used to convert Equation (3.2) to an unconstrained problem
in Equation (3.3):
12
2 2 2
2 1 1 1 2 2 2 2
u, , , ,
min            
t
t
v S e e
y UFu v TS v e Ex u e
(3.3)
Where
12
, ee are Lagrange multiplier terms, and
12
, are penalties for the axillary
dummy variables. With the Lagrange terms (
12
, ee ), the penalty parameters (
12
, ) need
not tend to large values for the equivalence of Equation (3.1) and Equation (3.2) to hold;
the value of these parameters does not affect the final solution, just the rate of convergence.
We empirically chose them to be 0.05.
Denoting the objective function in Equation (3.3) as
,,
t
L S u v , the problem can be
decoupled to simpler welldefined subproblems, where every subproblem takes an
analytical form, and can be solved in a single step. The algorithm relies on iterating
between these subproblems until convergence to a guaranteed global minimum. These
steps are:
31
0 0 0 0 0
11
2 1 2 2 3 3
1 1 1 1
2 2 2
1
Select , , , and 0
argmin ( , , ) [ ' ' ] ( '( ) '( ))
argmin ( , , ) [ ' ' ] ( ' ( ))
argmin
t t t
n n n n n n n
tt
x
n n n n n
t u u u t
u
n
S u ES v TS n
S L S v u E E T T C u e T v e
u L S v u F F E E F y ES e
v
Initialization :
Repeat :
1 1 1
11
1 1 1
11
1 1 1
22
( , , ) ( ; / )
()
()
1
stopping criterion is met
n n n n
tt
v
n n n n
t
n n n n
t
L S v u shrink TS e
e e v TS
e e u CS
nn
Until
Note that E’E=I because of the properties of the sensitivity maps, the above steps were
simplified into 5 steps compared to 7 steps in Ref (74), which accelerated the reconstruction
and simplified the workflow.
In this study, regularization penalties were chosen empirically based on retrospective
studies. A fully sampled DCE data set was retrospectively undersampled at the same
acceleration rate as our prospective data, then repeatedly reconstructed with a range of
constraint penalties. Normalized root Mean squared error (nRMSE) was calculated
between the fullysampled and reconstructed data sets. Penalties were chosen to maintain
nearminimal nRMSE, yet provide traction during reconstruction  essentially locating the
corner point of the lcurve (77). We employed penalties of 0.01, 0.0001, 0.0001 for
temporal finite difference, spatial TV and spatial wavelet respectively for all subsequent
reconstructions. In retrospective studies, we found that reconstruction converged well
within 100 iterations. We therefore allowed a maximum of 100 iterations in this
prospective study, to provide control of the maximum reconstruction time.
3.2.2 Experimental Methods
Fifteen brain tumor patients were recruited from three of our affiliated sites. Informed
consent was obtained from patients prior to MRI scan. Our Institutional Review Board
approved this study and all procedures.
MRI scans were performed on a clinical 3T scanner (HDxt, GE Healthcare, Waukesha,
WI) with an eightchannel head coil. Two DCE scans were performed: a standard
(“conventional”) scan using the vendor supplied sequence and our highly accelerated
(“experimental”) scan. Both DCE acquisitions used a 3D SPGR sequence. Prior to each
32
scan, T1 maps were acquired using variable flip angle DESPOT1 method (73). The standard
clinical postcontrast T1 weighted scans (Coronal T1weighted FSPGR sequence, 1.0 mm
3
isotropic resolution, 22x22x20 cm
3
FOV) were used as the reference for lesion
identification (4, 7).
Table 3.1 shows the conventionalthenexperimental imaging protocol, which required
less than 50 minutes. The contrast agent, Gadobenate dimeglumine (MultiHance Bracco
Inc.) was administered with a dose of 0.05 mMol/kg, followed by a 20 ml saline flush in
the left arm by intravenous injection for each scan (this results in a total of 0.1 mMol/kg
which would be the standard dose for a DCEMRI with contrast). The two DCE scans were
separated by approximately 20 minutes, resulting in residual GBCA present in the second
scan. This residue results in underestimation of K
trans
values for the second scan.
Table 3.1 Comprehensive MRI protocol for a brain tumor patient in our study with a
conventionalthenexperimental DCEMRI order. Conventional and experimental DCE
scans are in bold, and separated by roughly 20 minutes. The conventional and experimental
scans were reversed for the last two patients.
Duration
(mm:ss)
Sequence
1:00 Localization and SENSE calibration
4:00 Precontrast Axial T1w FSPGR
5:00 Precontrast Axial T2w FSE
5:00 Precontrast Axial FLAIR
0:45 T1 mapping
4:08 Conventional DCE
15:00 Diffusion Tensor Imaging
0:45 T1 mapping
5:29 Experimental DCE
2:00 Postcontrast Axial T1w FSPGR
3:00 Postcontrast Coronal T1w FSPGR
<50:00 TOTAL
In this study, a Patlak model was used to estimate TK parameters for conventional and
experimental scans. The following equation shows the fitting of K
trans
and vp from contrast
concentration curves of plasma Cp(t) and tissue Ct(t):
0
( ) ( ) ( )
t
trans
t p p p
C t K C d v C t
(3.4)
33
Cp(t) is the arterial input function (AIF), where we used a populationbased analytic
form (36), and ΔCt(t) is the GBCA concentration change in the tissue, which is derived
from the signal intensity in the dynamic images. K
trans
is the volume transfer constant and
vp is the fractional plasma volume.
The Patlak model is used because it is robust to noise, in part due to its simplicity and
the dependence on linear (v.s. nonlinear) estimation(40, 78). The model is based on the
assumption that there is no backflux from the interstitium during a short scan (40). This
assumption is not satisfied for the second injection in which residual GBCA in the
interstitium from the first scan causes backflux that cannot be ignored. Use of this model
causes underestimation of K
trans
values for the second scan, as we illustrate below using
simulations.
The doubleinjection experiment is simulated using the more accurate twocompartment
exchange model (2CXM). Figure 3.1 (a) shows simulated Cp(t) and Ct(t) with 20 minutes
separation between injections. Figure 3.1 (b) shows the simulated ΔCt(t) from first and
second injections, where the arrows indicate the altered curve due to residual GBCA. This
residual, indicated with yellow circles in Figure 3.1 (a), is directly related to the
extracellular extravascular volume fraction and the GBCA concentration in the
extracellular extravascular space prior to the second injection. The end result is 16% to
50% underestimation of K
trans
for the second scan. This projected range is based on TK
parameter values reported in the literature for brain tumor (78). Figure 3.1 (c) shows
measured ΔCt(t) curves from a representative tumor voxel in one case in our study.
To further determine if this underestimation is caused by injection order, we performed
experiments in three brain tumor patients, where the conventional DCEMRI scan was
performed twice in a single session. In the three subjects, tumor K
trans
was underestimated
in the second scan by 15%, 0%, and 24%. We consider this data to be anecdotal, but
consistent with the simulations and the observed underestimation in the second scan of the
conventionalthenexperimental and experimentalthenconventional patient scans. Figure
3.2 (ad) contains zoomed anatomic images from the 24% case before and after contrast
injection for the two scans, where contrast residue can be seen in the second scan. Figure
3.2 (e) and (f) contain zoomed K
trans
maps, where the second scan produced lower estimated
K
trans
compared to the first scan. In this case the mean K
trans
values in tumor ROI from the
second scan was 24% lower than the first scan. Figure 3.2 (g) shows measured ΔCt(t)
curves in these two conventional scans, which matches the trend observed in simulations,
conventionalthenexperimental scans, and experimentalthenconventional scans.
34
Figure 3.1 (a) The 2CXM simulated contrast concentration curves Cp(t) (scaled by 0.2× )
and Ct(t) for double injections, separated by 20 minutes. (b) ΔCt(t) calculated from the 1st
and 2nd injections, where ΔCt(t) from the 2nd injection is lower primarily due to high
backflux from the tissue to the plasma. (c) Actual measured ΔCt(t) curves in one tumor
voxel from a representative patient data, which matched the trend of the simulation.
Figure 3.2 (a), (b): Pre and post contrast images of first conventional DCE scan, cropped
around an enhancing tumor. (c), (d): Pre and post contrast images from the second
conventional DCE scan. (e): K
trans
maps from first DCE scan, (f): K
trans
maps from second
DCE scan. (g): Measured ΔCt(t) curves in the tumor ROI for the two conventional scans.
Conventionalthenexperimental protocols were performed in the first 13 patients, for
the last two cases, we switched the order of the conventional and experimental scans and
created an experimentalthenconventional protocol. This was for the purpose of verifying
that K
trans
underestimation in the second scan was due to the scan order and residual GBCA,
and not due to the imaging methods.
Table 3.2 lists acquisition parameters for the two DCE scans. The experimental scan
achieved much smaller voxel dimensions and wholebrain coverage, while maintaining the
same temporal resolution as the conventional scan. A net acceleration factor of 30× was
35
achieved. All 15 brain scans had the same fieldofview, matrix size, voxel dimensions,
scan time, and phase encode order. (k,t)space was undersampled in the exact same fashion
for all subjects. The injection delay for the experimental method was 60 sec, compared to
20 sec for the conventional method, in order to allow time for fully sampling of 25% of
phase encodes prior to contrast arrival.
Table 3.2 Scan parameters for the standard conventional and experimental highresolution
wholebrain DCEMRI scans. The experimented scan slice thickness is less than one third
of that of conventional scan, and the slice number is 17 times greater than that of
conventional scan. This enables a wholebrain nearisotropic coverage of the experimental
scan while keeping the same temporal resolution.
Conventional Experimental
TR/TE 6ms/2ms 6ms/2ms
Flip Angle 15⁰
15⁰
Matrix size 256× 186× 6 256× 256× 100
FOV (cm
3
) 22× 22× 4.2 22× 22× 19
Voxel dimensions (mm
3
) 0.93× 1.3× 7 0.93× 0.93× 1.9
Temporal resolution 5 s 5 s
Injection delay 20 s 60 s
Total scan time 4min 08s 5min 29s
Time frames 50 61
Sampling Pattern
Cartesian 3DFT
linear order
Cartesian 3DFT
goldenangle radial order
Acceleration factor 2x 30x
3.2.3 Comparisons / Evaluation
The conventional and experimental DCEMRI scans were registered based on the peak
contrast images using MATLAB image registration toolbox. As the experimental scan had
wholebrain coverage, slices within the FOV of the conventional scan were located by
registration. For fair evaluation, three adjacent slices of the highresolution experimental
scan were averaged to match the slice thickness of the conventional scan. Then K
trans
maps
were computed using the Patlak analysis (22, 40, 79), and a populationaveraged analytic
arterial input function (AIF) (36).
Quantitative evaluation was performed using the K
trans
histogram within radiologist
defined ROI (66, 80–82). Under the guidance of an experienced neuroradiologist (20 years
36
of experience), tight ROIs were drawn on K
trans
maps to included the highest K
trans
values.
The maximum K
trans
values were calculated for both scans and compared.
Qualitative evaluation was performed using radiologists’ rating. Two experienced
neuroradiologists (10 and 20 years of experience respectively) from our institution
reviewed and scored the images. For each subject, three types of images were shown for
the conventional and experimental scans (a total of six image sets): 1) timeresolved images
of one slice through the tumor (three slices averaged in experimental scan), 2) post
contrast–enhanced images (no slice average for experimental scan) and 3) BBB
permeability K
trans
maps. Radiologists were blinded to the acquisition type (conventional
or experimental), and the presentation order for every scan was randomized. A 4point
Likert scale was used to quantify the general image quality, where 3 = good, 2 = average,
1 = poor, 0 = nondiagnostic.
In the second round of qualitative evaluation, the same Neuroradiologists were shown
the full resolution conventional and experimental DCEMRI results, and were asked to
evaluate three subcategories of image quality: (1) SNR, which incorporates the visual
appearance of noise, (2) apparent spatial resolution, which incorporates sharpness of the
images, and (3) conspicuity of tumor enhancement, which incorporates the detectability
and sensitivity of contrast enhancement in the tumor. For each category, readers were asked
to determine if the conventional scan was superior, if the two were equal, or if the
experimental scan was superior. Conventional and experimental scan images (postcontrast
images, timeresolved images, and K
trans
maps) were shown to the radiologists, who were
not blinded to the scan type.
Results
3.3.1 Experimental Results
Table 3.3 contains the demographic and clinical information for the 15 patients that
were included in the data analysis. All experimental data sets were reconstructed using the
same empirical constraint penalty values (the λ’s in Equation (3.1)). Reconstruction time
was roughly 8 hours per dataset.
37
Table 3.3 Patient demographic and clinical information of the fifteen brain tumor patients
participated in the study.
No. Sex Age Symptom
Tumor size
(cm)
*
001 F 46 Glioblastoma 2.4
002 M 71 Glioblastoma 2.3
003 F 76 Meningioma 1.9
004 F 53 Metastasis 3.4
005 M 26 Astrocytoma 1.6
006 M 77 Meningioma 0.8
007 M 72 Metastatic melanoma 1.0
008 M 65 Glioblastoma 6.0
009 M 71 Glioblastoma 1.5
010 F 65 Metastatic ovarian cancer 1.3
011 F 38 Glioblastoma 0.5
012 F 72 Meningioma 1.4
013 F 22 Glioblastoma 2.4
014 F 78 Metastatic Melanoma 1.6
015 F 59 Meningioma 3.9
* Tumor size was measured by the longest diameter in the postcontrast T1w images.
Figure 3.3 shows typical postcontrastenhanced DCEMRI images from two subjects,
one with a large 6 cm glioblastoma; the other with multiple metastatic melanoma tumors
scattered throughout the entire brain. The experimental approach was able to provide
detailed depiction of the entire tumor body and tumor boundary, and capture all the possible
small lesions (14 in total) in the whole coverage of the brain (see arrows). The tumor
boundary was clearly visualized in any scan plane, and small, scattered lesions were easily
identified in coronal/sagittal reformats. Conversely, the conventional scan provided limited
spatial coverage (only 4 lesions were captured), and the sagittal and coronal reformats had
extremely low resolution in the slice encoding direction. In all cases, the experimental scan
provided clearer and crisper depictions of all lesions that were presented.
38
Another benefit of wholebrain coverage is that it essentially eliminates inflow
enhancement artifacts. Figure 3.3 (c) shows bright signal in the sagittal sinus stemming
from inflow enhancement in a conventional scan where only a 6cm axial slab is imaged.
This inflow enhancement was consistently observed in every conventional scan and is
strongest in slices at the edge of the imaging slab. In contrast, Figure 3.3 (d) shows that
inflow enhancement was not present in the experimental scans.
Figure 3.4 demonstrates registered anatomic images and K
trans
maps from two other
representative patients, with conventionalthenexperimental and experimentalthen
conventional protocols respectively. Note that the experimental images were blurred in the
sliceencoding direction to match the conventional protocol. Anatomic images yielded
similar quality in the regions of interest. K
trans
measurements from the second scan were
consistently 42% to 66% of that from the first scan in the ROI of the tumor, regardless of
the scan order. The K
trans
maps, despite intensity differences (due to scan order), provided
comparable information and superior image quality based on the radiologists’ ratings,
described below.
39
Figure 3.3 Final DCEMRI time frames from two patients illustrating the volume coverage
of (a,c) experimental and (b,d) conventional scans. (a,b): 65/M patient with a 6 cm
glioblastoma tumor (subject #008). Experimental scan (a) shows significantly larger
coverage than the spatial coverage of the conventional scan (b), where on the sagittal and
40
coronal reformats only a thin slab of brain can be covered. (c,d): 78/F patient with
metastatic melanoma (subject #014), only 4 lesions were captured by the clinical scan (d),
all 14 lesions were captured by the experimental scan (c). Note that the conventional scans
show bright signal in the sagittal sinus (red arrows) and other blood vessels due to inflow
enhancement that are not present in the experimental scans because of the wholebrain
coverage. This was consistently observed in all subjects.
Figure 3.4 Registered anatomic images and K
trans
maps from two other representative
patients. On the left are the registered anatomic images; on the right are K
trans
maps. The
top row is from patient #013 with a glioblastoma (see arrows) using conventionalthen
experimental protocol. The bottom row is from patient #015 with a meningioma (see
arrows) using experimentalthen conventional protocol. Both sets of anatomical images
provide comparable features and image quality, and the K
trans
maps convey comparable
diagnostic information despite intensity differences.
3.3.2 Quantitative Assessment
Figure 3.5 shows a scatter plot of the maximum K
trans
value within manually segmented
tumor ROIs between two DCEMRI scans. Each pair of scans was performed using one of
three possible orderings: conventionalthenexperimental (blue circle), experimentalthen
conventional (green diamond), and double conventional (red star). Although the second
scan underestimates K
trans
, the two measurements were still highly correlated with
correlation coefficient r=0.513. The mean difference between the two scans was 0.036,
which corresponds to a consistent negative bias consistent with contrast residue from the
previous injection.
41
Figure 3.5 Scatter plot of the maximum K
trans
in tumor ROIs for the 1
st
and 2
nd
DCEMRI
scans. ConventionalthenExperimental (13 cases), ExperimentalthenConventional (2
cases) and Double Conventional (2 cases) are all shown. The correlation coefficient was
0.5132. The second scan experienced a consistent underestimation of K
trans
due to contrast
residue.
3.3.3 Qualitative Assessment
Table 3.4 lists the two radiologists’ ratings of the images. The Likert scale scores of
overall image quality were averaged across image types, and shown for all 15 patients. The
three subcategories, SNR, apparent resolution, and conspicuity of tumor enhancement,
were shown as experimentalbetter (denoted by ‘+’), scans are equal (denoted by ‘=’), or
conventionalbetter (denoted by ‘–‘). Both radiologists consistently rated the experimental
scans as higher or equal in quality to the conventional scans in terms of SNR, effective
resolution, and fine details. The conventional scan was not deemed superior to the
experimental scan in any of the cases by either radiologist. The radiologists also indicated
that they observed “better white/gray matter contrast,” “improved resolution and edge
sharpness,” “reduced phaseencoding artifacts,” “reduced noise level,” and “better
detection of tumor with large coverage of the brain,” from the experimental scan images.
42
Table 3.4 Two radiologists’ scores across the 15 patients for conventional and
experimental DCEMRI scans. A 4point Likert scale was used to score the overall image
quality (3 = good, 2 = average, 1 = poor, 0 = non diagnostic), and the average of this score
was taken for each patient across the three image types. Three subcategories of image
quality (SNR, apparent resolution, conspicuity of tumor enhancement) were scored as: +:
experimentalbetter, =: equal, : conventionalbetter. (I: Radiologist 1. II: Radiologist 2).
NO.
Average Likert score
SNR Resolution Conspicuity
Conventional Experimental
I II I II I II I II I II
1 1 1.67 2 1 + + + + + +
2 1 1 2.33 1.67 + + = + + +
3 1 1 2 1.67 + + + + = +
4 1.67 1.33 2.67 1.67 + + + + = +
5 1 1.33 2.33 1.67 + + + + = +
6 1.33 1 2.67 2 + + + + = +
7 1 2.33 3 2 = = = + + +
8 1 1.33 3 2.33 + + + + + +
9 1 1.67 2.33 2 + + = + + +
10 1 1.33 2.33 1.33 + + = + = +
11 1 1.33 3 2 + + + + = +
12 1 1.33 3 2 + + + + = +
13 1 1 3 2.67 + + + + = +
14 1 1 3 2 + + + + = +
15 1 0.67 2 2 + = + + = =
Histograms of the Likert scores for both conventional and experimental scans are shown
in Figure 3.6. The scores for the three image types were combined to show the overall
performance of conventional and experimental scans. Qualitative evaluation from the
experimental scans (mostly 2 and 3) clearly outperformed that of the conventional scans
(mostly 1).
43
Figure 3.6 Histogram of all conventional and experimental scores combing the three image
types. The statistics of conventional and experimental scans are 1.2± 0.6 and 2.2± 0.7,
respectively.
Discussion
We have implemented a prospective undersampling and constrained reconstruction
scheme for highresolution wholebrain DCEMRI. We have performed a pilot comparison
study in fifteen brain tumor patients that has demonstrated the strength of this technique
and its potential impact on clinical DCEMRI.
Specifically, we have shown that the DCEMRI with constrained reconstruction is able
to provide much higher spatiotemporal resolution and wholebrain spatial coverage
compared to current DCEMRI methods. This is extremely important when imaging large
tumors or patients with multiple metastatic lesions, which the conventional scans fail to
completely capture due to poor coverage and/or low spatial resolution. Within the same
ROI of the conventional scan, Radiologists’ reported improved quality, comparable or
better diagnostic information in both anatomic images and K
trans
maps of the experimental
scan.
The experimental images show better SNR, resolution and lesion conspicuity, as well
as overall image quality score, despite significant undersampling. This is due to novelty
in the way that raw data are reconstructed. In the conventional scan, each time frame is
reconstructed independent of every other time frame. In the experimental scan, each time
44
frame is undersampled in kspace, but all time frames are reconstructed in a single step
and the mutual information between time points is leveraged through the use of temporal
constraints. This is precisely the reason that compressed sensing and constrained
reconstruction methods can achieve high acceleration factors with equivalent or superior
image quality, despite undersampling each time frame. In DCEMRI (and similarly, time
resolved MR angiography) rapid temporal changes are limited to spatial positions
containing vessels, and temporal changes elsewhere are smooth. This makes the images
sparse after a temporal finite difference or high pass filter operator is applied. The
experimental method leverages information from several time points, which boosts the
effective signaltonoise ratio and image quality, even at high acceleration rates.
Regularization of this kind, with spatial and temporal constraints, is common in image and
video denoising applications not limited to MRI (83–85).
In a typical clinical protocol, conventional highresolution wholebrain 3D T1weighted
static image volumes were acquired before the GBCA injection and near the end of the
exam, following all GBCA administration. It is worth noting that the first and last time
frames of the experimental highresolution wholebrain DCE scan also constitute pre and
postcontrastenhanced T1weighted images. This could enable the experimental scan to
take place of additional static pre and postcontrastenhance image acquisitions, saving
scan resources and time.
There are several factors that preclude an optimal comparison with the conventional
protocol. First, both conventional and experimental scans were performed in the same
session with a relatively short interval. The rationale for this was practicality. Brain tumor
patients were unwilling to come back on another day for a second MRI scan that was purely
for research purposes. The consequence is that there is significant residual GBCA on board
prior to the second DCE MRI scan, which could lead to underestimation of K
trans
values in
the second scan. The order of conventional and experimental scans was reversed in the last
two cases; and the same K
trans
underestimation was observed in the second scan, consistent
with residual GBCA (not the imaging technique) being the cause. Additionally, due to the
short interval between the two scans, tumors that take up contrast slowly may show higher
conspicuity on the second scan due to temporal order (e.g. Figure (3.2), top row). A
limitation of this study is that only two cases were performed with the experimentalthen
conventional order. In both cases, the experimental (first) scan was rated as having higher
tumor conspicuity than the conventional (second) scan, but this should be considered
anecdotal because of the small sample size.
Another limiting factor for constrained reconstruction is the high computational
complexity. Using a powerful workstation and an efficient MATLAB implementation, data
reconstruction per dataset required 8 hours at the time the studies were performed. Since
45
then, we have been able to shorten the reconstruction time to 1.5 hours by incorporating
coil compression (86), reduced temporal segmentation, and optimization of the MATLAB
implementation. In general, reconstruction time is a limitation for iterative constrained
reconstruction techniques, making them problematic for realtime imaging or applications
that require immediate availability of image data. Parallelization, GPU based computation,
and efficient Cbased implementation are being explored by many groups, and these have
provided reconstruction speedup on the order of 3 to 200 (52, 87).
This study did not include a quantitative assessment of spatial and contrast resolution,
and this remains an important next step in the evaluation of constrained reconstruction
techniques in DCEMRI. Because the experimental method involves nonlinear
reconstruction, the definition of true spatial and contrast resolution is nontrivial and is an
open research question. Partial solutions such as characterization of local point spread
functions (88), and validation with digital or physical reference objects (89) exists but
require dedicated investigation and are highly objectdependent. There are efforts
underway to develop anatomically realistic brain DCEMRI digital reference objects (90,
91), which will facilitate such work.
Undersampling and constrained reconstruction for MRI is a rapidly developing area.
Several groups, including ours, are developing constraints and quality evaluation
techniques to maximize the reconstruction quality. The constraints used in this study have
been previously validated in a retrospective study (56), and the prospective study results
shown here further demonstrate the application and feasibility of constrained
reconstruction for clinical DCEMRI
Conclusion
We conclude that highresolution wholebrain DCEMRI using constrained
reconstruction is clinically feasible and provides superior image quality compared to the
current conventional DCEMRI technique. In our study, the experimental approach
provided superior image and pharmacokinetic map quality without compromising
diagnostic information compared to the current DCEMRI approach. The experimental
approach also provided complete characterization of all normal and abnormal tissues, and
allowed for arbitrary multiplanar reformatting of data. This was a significant advantage in
two of the fifteen cases, one with a large glioblastoma multiforme that exceeded the spatial
coverage of the conventional scan and had a narrow enhancing margin, the other with 14
metastatic lesions of which only 4 were characterized by the conventional scan. This study
represents, to the best of our knowledge, the first prospective evaluation of brain DCE
MRI with constrained reconstruction. Compared to current best practices, this new
46
approach has the potential to vastly improve visualization and characterization of brain
lesions with DCEMRI.
47
Chapter 4. Modelbased Direct Reconstruction for DCEMRI
Introduction
As described in previous chapters, current DCEMRI with Nyquist sampling is unable
to simultaneously provide high spatiotemporal resolution and adequate volume coverage.
Compressed sensing (44) and parallel imaging (17, 20) based schemes have been proposed
to accelerate acquisition process, primarily to achieve better spatial resolution and coverage
while maintaining the same temporal resolution. Notably, Lebel et al. (56) used a temporal
high pass filter and multiple spatial sparsity constraints to achieve an undersampling rate
(R) of 36x, and showed excellent quality of anatomic images in brain tumor cases. A recent
pilot study in brain tumor patients indicated that this approach performs superior to
conventional techniques with no apparent loss of diagnostic information (92). The works
of Feng et al. (49), Chandarana et al. (93), Rosenkrantz et al. (52), used a goldenangle
radial sampling pattern, compressed sensing, and parallel imaging to achieve a comparable
acceleration rate of 19.1 to 28.7. These studies showed improved resolution and reduced
motion sensitivity in breast, liver, and prostate DCEMRI, compared to either parallel
imaging alone or coilbycoil compressed sensing alone. We will refer to these techniques
as “indirect” methods, since the anatomic image series are reconstructed first, followed by
a separate step for TK parameter fitting on a voxelbyvoxel basis.
In this chapter, we propose a framework for “direct” estimation of TK parameter maps
from fullysampled or undersampled (k,t)space data. We employ a full forward model that
converts the TK maps to (k,t)space, and we pose the estimation of TK maps as an error
minimization problem. Our approach is motivated by two factors: 1) spatial TK parameter
maps have much lower dimensionality than those of dynamic image series (24 parameters,
compared to 50100 time points, per voxel), and 2) TK modelbased reconstruction directly
exploits what is known about contrast agent kinetics. These allow for robust parameter
estimation from an information theoretic perspective, and has the potential to provide the
most accurate restoration of TK parameter values, and allow for the highest acceleration.
Modelbased direct reconstruction has been previously explored in other applications
such as MRI relaxation parameter estimation (94–99), TK parameter estimation in PET
(100–102), and TK parameter estimation in MRI (103, 104). Notably, for MRI relaxation
parameter estimation, Sumpf et al. (97) used a modelbased nonlinear inverse
reconstruction to estimate T2 maps from highly undersampled spinecho MRI data; Zhao
et al. (98) estimated T1 parameters directly from undersampled MRI data with a sparsity
constraint on the parameter maps. For dynamic PET imaging, Kamasak et al. (101) directly
estimated TK parameter images from dynamic PET data using a kinetic modelbased
48
reconstruction optimization; Lin et al. (102) used a sparsity constrained mixture model to
estimate TK parameters from dynamic PET data, and evaluated with both simulated and
experimental PET data.
For DCE MRI parameter mapping, Felsted et al. (103) proposed to use a modelbased
reconstruction algorithm to solve TK parameters directly from undersampled MRI kspace
with a modified gradient descent algorithm. An undersampling factor of R=4 was
demonstrated on simulated data; Dikaios et al. (104) proposed a Bayesian inference
framework to directly estimate TK maps from undersampled MRI data, and achieved 8x
acceleration in phantom and invivo prostate cancer data.
While the prior studies demonstrate promise, the full potential of TK model based
reconstruction has lacked validation, both at higher undersampling rates, or with
prospectively undersampled data from patients. In this study, we explore the maximum
potential benefit of the direct TK estimation by testing very high undersampling rates. We
validate the approach using retrospective undersampling of fullysampled data and using
prospectively undersampled DCE data sets from brain tumor patients. Compared to prior
work, we are able to demonstrate much higher undersampling rates (up to 100× ) in the
retrospective study, using a more efficient gradientbased algorithm. We use quantitative
evaluation (root Mean Square Error (rMSE) in TK parameters) to provide a systematic
comparison against a stateoftheart compressed sensing method that uses spatial and
temporal sparsity constraints in 13 brain tumor patients. We also uniquely provide a
prospective invivo study showing that wholebrain coverage with high spatial resolution
can be achieved to capture complete pathological information. We demonstrate the
potential of direct reconstruction to enable “parameterfree” reconstruction, when no
sparsity constraints are added.
Theory
4.2.1 Direct TK Mapping
We propose to integrate TK modeling, specifically the Patlak model, into the image
reconstruction process. Figure 4.1 illustrates the forward model that relates TK parameter
maps to undersampled (k,t)space.
49
Figure 4.1 DCEMRI forward model flowchart illustrating the conversion from TK
parameter maps to undersampled (k,t)space. Patlak model is used to convert TK parameter
maps to contrast concentration over time, after which the T1weighted signal equation is
used to obtain the dynamic anatomic images. Fourier transform, sensitivity maps and
sampling pattern connect anatomic images to multicoil (k,t)space measurements.
We use the vector ( , , ) r x y z to represent image domain spatial coordinates,
( , , )
x y z
k k k k to represents kspace coordinates; t, c are the time and coil dimensions. The
variables beneath the arrows of each step are known or predetermined. The steps indicated
(above the arrows) in Figure 4.1 are explained below:
1). Contrast agent concentration over time ( , ) CA t r is assumed to follow the Patlak
model:
0
( , ) ( ) ( )d ( ) ( )
t
trans
p p p
CA t K C v C t
r r r
(4.1)
where ()
p
Ct is the arterial input function (AIF). In this work, we used a population
based AIF from (36). Notice that the AIF requires specifying a delay time. This is estimated
from the kspace origin, which is acquired in every time frame, and has shown to accurately
detect the time of contrast bolus arrival (105). We assume that the Patlak is appropriate for
all voxels in the imaging volume. We have observed that image regions outside of vessels
and tumor typically experience no enhancement during the DCEMRI acquisition which
results in a fit to vp=0, K
trans
=0.
2). Dynamic anatomic images ( , ) st r are related to ( , ) CA t r by the steady state spoiled
gradient echo (SPGR) signal equation:
50
1 1 1
1 1 1
[ ( ,0) ( , ) ] ( ,0)
00
[ ( ,0) ( , ) ] ( ,0)
( )sin (1 ) ( )sin (1 )
( , ) ( ,0)
1 cos 1 cos
TR R CA t r TR R
TR R CA t r TR R
M e M e
s t s
ee
r r r
r r r
rr
rr
(4.2)
where TR is the repetition time, α is the flip angle,
1
r is the contrast agent relaxivity,
1
( ,0) R r and
0
() M r are the precontrast R1 (reciprocal of T1) and the equilibrium
longitudinal magnetization that are estimated from a T1 mapping sequence. In this work,
we used DESPOT1 (29) immediately prior to the DCEMRI scan. ( ,0) s r is the precontrast
firstframe, which is fullysampled in this work. The bracketed [] term resolves differences
between the precontrast signal and the predicted precontrast signal based on the baseline
T1 and M0 maps (from DESPOT1 sequences) (31).
3). The undersampled raw (k,t)space data ( , , ) S t c k is related to ( , ) st r by the coil
sensitivities ( , ) Cc r , and undersampling Fourier transform (Fu):
( , , ) ( , ) ( , )
u
S t c F C c s t k r r (4.3)
In this work, ( , ) Cc r is estimated from time averaged data using the standard root sum
ofsquares method (17). The image phase information is assumed to be captured by the
complexvalued sensitivity maps.
Combining Equations 4.13, we reach a general function f to denote the relationship
between TK maps ( ), ( )
trans
p
Kv rr and undersampled (k,t)space ( , , ) S t c k .
1 0 1
( , , ) ( ( ), ( ); ( ), , , ( ,0), ( ), , ( , ))
trans
pp
S t c f K v C t TR R M r C c k r r r r r (4.4)
where
1 0 1
( ), , , ( ,0), ( ), , ( , )
p
C t TR R M r C c r r r are variables that are known or pre
determined as mentioned above. We solve for the unknown ( ), ( )
trans
p
Kv rr via least
square optimization, formulated as follows:
2
2
( ), ( )
( ( ), ( )) arg min  ( , , ) ( ( ), ( )) 
trans
p
trans trans
pp
Kv
K v S t c f K v
rr
r r k r r
(4.5)
This nonlinear optimization problem is solved by a quasiNewton limitedmemory
BroydenFletcherGoldfarbShannon (lBFGS) method (106), where ( ) and ( )
trans
p
Kv rr
are solved alternatingly. The details of the optimization algorithm and gradient calculation
are provided in Appendix 4A.
51
Direct reconstruction by itself is “parameterfree”. This is in contrast to compressed
sensing based algorithms that require tuning of one or more regularization parameters. It is
possible to incorporate additional spatial sparsity constraints on the TK maps themselves.
In this work, we test the potential value of adding a spatial ‘db2’ wavelet constraint (Ψ) to
the parameter maps. The optimization problem with sparsity constraint is formulated as
follows:
2
2 1 1 2 1
( ), ( )
( ( ), ( )) argmin  ( , , ) ( ( ), ( ))   ( )   ( ) 
trans
p
trans trans trans
p p p
Kv
K v S t c f K v K v
rr
r r k r r r r
(4.6)
4.2.2 Indirect TK Mapping
Current stateoftheart methods for highly accelerated DCE involve reconstructing
intermediate images prior to TK modeling. These indirect methods are the most relevant
alternatives to direct TK modeling and serve as a performance benchmark. A basic model
for indirect reconstruction solves the minimization problem in Equation (4.7), where the
final image, ( , ) st r , remains consistent with acquired (k,t)space data ( , , ) S t c k , yet is
sparse in the temporal finite differences (V) domain and spatial wavelet domain (Ψ).
2
2 1 1 2 1
( , )
( , ) argmin  ( , , ) ( ,c) ( , )   ( , )   ( , ) 
u
s
s t S t c F C s t Vs t s t
rt
r k r r r r
(4.7)
The image is related to the acquired data using known coil sensitivities ( , ) Cc r and the
undersampling Fourier transform Fu. TK modeling (e.g. using a Patlak model) is performed
in a last step to estimate the spatial TK maps (e.g. ( ), ( )
trans
p
Kv rr ) from s( , ) t r . This
optimization problem was solved by an efficient augmented Lagrangian method,
alternating direction methods of multipliers (ADMM) (74), to get the anatomic images.
Methods
4.3.1 Digital Phantom
We simulated realistic DCEMRI data using a digital phantom with known TK
parameter maps and using the Patlak TK model. We used a process identical to Ref. (90),
where the segmentation is extracted from patient data. Realistic sensitivity maps were used,
and noise was added to each channel according to noise covariance matrix estimated from
the patient data. A precontrast whitematter SNR level of 20 was chosen to mimic the
SNR level in actual DCE data sets.
52
We retrospectively undersampled (k,t)space with rates R of 20x to 100x. Ten noise
realizations were generated for each R. Undersampling was in the k
x
k
y
plane, simulating
the k
y
k
z
plane as in a prospectively undersampled 3D case, using a randomized golden
angle radial sampling pattern (53, 57). Detailed description and the videos of the
undersampling strategies can be found in and the supporting materials. Direct and indirect
methods were used to generate the TK parameters from both fullysampled and
undersampled data. TK map rMSE were computed over an ROI containing the entire tumor
boundary.
4.3.2 InVivo Retrospective Evaluation
We reviewed 110 fullysampled DCEMRI raw data sets from patients with known or
suspected brain tumor, receiving a routine brain MRI with contrast on a clinical 3T scanner
(HDxt, GE Healthcare, Waukesha, WI). The data sets were from patients receiving routine
brain MRI with contrast (including DCEMRI) at our Institution, and the demographics
reflect our local patient population. Our Institution follows standard exclusion criteria for
MRI with Gadoliniumbased contrast (107, 108) which includes: medically unstable, renal
impairment, cardiac pacemaker, internal ferromagnetic device that is contraindicated for
use in MRI, claustrophobia, and any other condition that would compromise the scan with
reasonable safety. The retrospective study protocol was approved by our Institutional
Review Board.
The sequence was based on a 3D Cartesian fast spoiled gradient echo sequence (SPGR)
with FOV: 22× 22× 4.2cm
3
, spatial resolution: 0.9× 1.3× 7.0mm
3
, temporal resolution: 5s, 50
time frames, and 8 receiver coils. The flip angle is 15° , and TE is 1.3ms, TR is 6ms.
DESPOT1 was performed before the DCE sequence, where three images with flip angle of
2° , 5° , 10° were acquired to estimate T1 and M0 maps before the contrast arrival. The
contrast agent, Gadobenate dimeglumine (MultiHance Bracco Inc., which has relaxivity
1
r
=4.39 s
1
mM
1
at 37⁰C at 3 Tesla (109)) was administered with a dose of 0.05 mMol/kg,
followed by a 20 ml saline flush in the left arm by intravenous injection.
Of the 110 cases, we found 18 that had visible tumor larger than 1cm by bi
directional assessment. TK parameter maps K
trans
and vp were calculated from the fully
sampled images, and TK model fitting error was computed by taking the l 2 norm between
the contrast concentration curves from fullysampled images, ( , ) CA t r , and the fitted
concentration curves generated from the TK parameter maps,
ˆ
( , ) CA t r . We then examined
the Patlak modeling error, defined as:
53
2
2
2
2
ˆ
 ( , ) ( , )
100%
 ( , ) 
CA t CA t
CA t
rr
r
(4.8)
Of the 18 cases with visible tumor larger than 1cm, 13 cases had Patlak modelling error
less than 1%, suggesting that the Patlak model with the population AIF was appropriate.
The analysis below was performed on the 13 cases which were fullysampled, had at least
one tumor larger than 1cm, and for whom the Patlak model fitted the fullysampled data
with less than 1% error.
For each selected case, three sets of TK maps were generated from: 1) standard Fourier
reconstruction of fully sampled data. These served as the gold standard reference maps. 2)
Direct reconstruction method using retrospective undersampling; 3) Indirect reconstruction
using retrospective undersampling. We examined R of 20× to 100× , with at increments of
20× . For each R, 10 realizations of the sampling pattern were generated using a different
initial angle in the randomized goldenangle radial scheme. This effectively creates
multiple noise realizations, since there is almost no overlap in the (k,t)space sampling
pattern (except for the one sample at the kspace origin, which is included in every sampling
scheme at every undersampling factor).
For indirect method in Equation (4.7), the regularization parameters were empirically
set as λ1=0.01 and λ2=0.0001. These values are motivated by empirical observations made
on retrospective undersampling studies on a number of datasets (around 15) based on a
criterion of achieving minimal rMSE between the reconstructed dynamic images from sub
sampled data and fullysampled data. Both the regularizations were used in Equation (4.7),
during experimental comparisons against the direct method with spatial sparsity constraint
Equation (4.6). For direct method, regularization parameters were also empirically set as
λ1=0.03 and λ2=0.00001. For fair comparison, when the direct method had no spatial
wavelet constraint (Equation (4.5)), λ2 in Equation (4.7) of the indirect method was set to
0.
The quantitative metric rMSE was computed on TK parameter maps, within a region
ofinterest (ROI) containing enhancing tumor. BlandAltman plots were generated using
the difference of the reconstructed K
trans
maps with respect to the fullysampled K
trans
maps
within the ROI to test for any systematic bias.
A twotailed paired STUDENT’s ttest was performed based on the rMSE of the two
methods in the 13 patients. The null hypothesis was equivalence of the two methods, with
the null value being zero. The significance criterion was P value less than 0.05. The
assumptions of normality were validated using ShapiroWilk test. Bonferroni correction
54
was applied to correct for multiple comparisons, that is, the significance level for each
individual test was set to 0.05/13.
For one data set we applied spatial wavelet sparsity constraints for both methods
(Equation (4.6) and Equation (4.7)) to demonstrate the feasibility and determine any
possible improvement.
4.3.3 InVivo Prospective Evaluation
Prospectively undersampled data were acquired in 4 brain tumor patients (65 M, 71 M,
46 F, 22 F, all Glioblastoma) with Cartesian goldenangle radial kspace sampling (57, 92).
3D T1w SPGR data was acquired continuously for 5 minutes. Wholebrain coverage was
achieved with a FOV of 22× 22× 20 cm
3
and spatial resolution of 0.9× 0.9× 1.9 mm
3
. The
prospective study protocol was approved by our Institutional Review Board. Written
informed consent was provided by all participants.
Five second temporal resolution was achieved by grouping raw (k,t)space data acquired
within consecutive 5 sec intervals. This prospective acquisition undersampled each 5 sec
temporal frame by 30× . Note that undersampling is not being used to shorten the scan time,
but rather to significantly increase the spatial coverage and spatial resolution. For
comparison, the standard clinical protocol at our institution that utilizes Nyquist sampling
and 5 sec temporal resolution achieves FOV 22× 22× 4.2cm
3
and spatial resolution
0.9× 1.3× 7.0mm
3
. For DCEMRI, the scan time and temporal resolution is kept the same to
capture dynamic changes during contrast arrival and washout.
Direct estimation of TK maps was performed using the proposed method. Threeplane
and panning volume K
trans
and vp maps for the 4 data sets are presented for visual
assessment. The first frame is necessary in the direct reconstruction process. A detailed
description of how to obtain this frame, utilizing the properties of the goldenangle radial
sampling pattern, can be found in the supporting materials. These prospectively
undersampled DCEMRI data were obtained 20 minutes after a standardofcare
conventional DCEMRI scan, therefore there was some residual contrast on board.
Results
Figure 4.2 shows the convergence performance for the direct method at R=20. Objective
function changes are plotted against iteration number, and final results where minimum
gradient values reached are shown with different initial conditions. In the experiments,
140–180 iterations were needed to reach the stopping criteria. Convergence was achieved
regardless of the initial condition. The reconstruction time for indirect and direct method
55
was 265s and 296s respectively on Linux workstation (24 core 2.5GHz, 192GB RAM).
K
trans
is reported in the units of min
1
, and vp is reported as percent fraction.
Figure 4.2 Objective function versus iteration number for three initial TK parameter
estimates (left), and cropped portions of these initial and final TK maps (shown is the
Ktrans map) at a undersampling rate of 20x (right). All initial conditions converged to the
same solution. No benefit was observed using a zeropadded kspace relative to a null
starting condition.
Figure 4.3 top shows the phantom results (cropped at the tumor part) of indirect and
direct reconstruction of the K
trans
and vp maps at R=1 (fullysampled data) and R=60 for
both direct and indirect method. At R=1, it is shown that the direct method does not
introduce additional error by enforcing the Patlak model into the reconstruction. At R=60,
the direct method performs better than indirect method for both K
trans
and vp maps,
overcoming the large errors and noise introduced by indirect method. Figure 4.3 bottom
shows the rMSE (calculated in tumor boundary ROI) performance across different R.
Across all R tested, the direct method outperformed indirect method at high undersampling
rates.
56
Figure 4.3 Retrospective evaluation of indirect and direct methods on phantom data. The
top row contains ground truth K
trans
and vp maps that are used to generate the phantom,
Patlak fitting results from fullysampled but noiseadded data, and R=1× and R=60×
reconstruction results for both direct and indirect methods. Realistic noise (SNR=20) were
added to the simulated kspace data. The bottom row contains rMSE across undersampling
rates for a region of interest containing the entire tumor boundary (761 voxels). The
proposed direct reconstruction produced lower mean rMSE at all sampling rates.
Figure 4.4 shows one representative example of the image results of direct and indirect
reconstruction. K
trans
and vp maps of a glioblastoma patient at three different undersampling
rates obtained from fully sampled data, the proposed direct reconstruction, and via indirect
reconstruction are shown. Across all undersampling rates, the direct method qualitatively,
and quantitatively depicted equal or more accurate restoration of TK parameter maps
compared to indirect method. At R of 20× or less, the direct and indirect methods had
equivalent performance. At higher undersampling rates, the indirect method failed to
57
capture critical tumor signals and tumor shapes in K
trans
maps and small vessel information
in vp maps, while the direct method was able to provide accurate restoration. The results
with spatial wavelet constraints (only 100× shown here) provide improved noise
performance and image quality. It is also worth noting that the indirect method tends to
underestimate K
trans
values, while the direct method overcomes this underestimation. This
is better observed in the BlandAltman plots in Figure 4.6.
Figure 4.4 Retrospective evaluation of direct and indirect reconstruction of K
trans
and vp
maps. Both reconstructions are shown without spatial wavelet sparsity constraints in the
first three columns, and with the sparsity constraint in the last column. By visual inspection,
direct reconstruction outperformed indirect results at all undersampling rates. The direct
method provided superior delineation of the tumor boundary and other high resolution
features in the K
trans
maps than did the indirect method. This was particularly true at the
highest undersampling rate (see arrows at 100× ). For vp maps, the direct method better
preserved small vessel signals (see arrows at 100× ) compared to the indirect method.
Spatial wavelet constraints, shown on the rightmost column with 100× , provide additional
noise suppression.
58
Figure 4.5 shows tumor ROI K
trans
rMSE plots for both methods for a range of
undersampling rates, and all 13 data sets that we used for retrospective evaluation. The
error bars show the variance introduced by varying the initial angle of the sampling patterns.
The direct method outperformed the indirect method for all the cases at high undersampling
factors (>80× ). The cases are ordered (left to right) by decreasing performance of direct
reconstruction compared to indirect reconstruction. More detailed analysis can be found
in .
Figure 4.5 Tumor ROI (fullysampled) and rMSE plot for direct and indirect
reconstruction results across different R (20× to 100× ) and 13 data sets for K
trans
s values.
The error bar indicate the mean and variance of the rMSE for each R where 10 different
realization of sampling patterns were used. The direct method outperformed the indirect
method in most cases, especially at R>60× . The variance for the two methods are
comparable.
Figure 4.6 shows the BlandAltman plots of direct and indirect methods for the K
trans
values combining the ROIs of all the 13 cases. The indirect method tended to suppress
59
K
trans
values smaller than 0.02 min
1
, and tended to underestimate values larger than 0.02
min
1
. This is similar to a softthresholding operation, and is illustrated by the green line
in Figure 4.6. This may be a sideeffect of the temporal finite difference constrained
reconstruction that suppresses small temporal changes of concentration. Such a trend was
not observed in the direct reconstruction, suggesting that tighter integration the TK model
is able to identify and restore low K
trans
values.
Figure 4.6 BlandAltman plots of the difference between estimated K
trans
and reference
K
trans
(from fullysampled reconstruction) for both direct and indirect method at R=100x.
Each dot corresponds to one voxel within the tumor ROI of one of the 13 cases. The
indirect reconstruction (right) demonstrated a pattern of underestimating K
trans
values,
especially for K
trans
< 0.02 min1 (see green line). This may be a sideeffect of the temporal
finite difference constraint suppressing small concentration changes. In contrast, the direct
reconstruction (left) did not demonstrate any considerable bias patterns, and had a lower
variance.
Table 4.1 lists the patient demographic information, and the performance of direct and
indirect methods at R=60× , evaluated by rMSE for both K
trans
and vp maps. The mean and
standard deviation of rMSE are listed. The order of presentation in this table matches the
order of Figure 6. The normality assumptions are met by ShapiroWilk test, and only 4 out
of 56 cases (13 patients, 5 undersampling rates) reject the null hypothesis of composite
normality assumption at significance level of 0.05. A twotailed paired STUDENT’s ttest
was performed between the two methods based on rMSE of K
trans
, and the last two columns
show the smallest R where the direct method started to statistically significantly outperform
60
the indirect method (p<0.0038 after Bonferroni correction), and the corresponding pvalue.
The direct method was consistently better than indirect method across multiple data sets at
high undersampling rates, and the difference was statistically significant at R>20× for 6
cases, R>40× for 3 cases, R>60× for 1 case, and R>80× for 3 cases.
Table 4.1 Patient demographic information and the rMSE performance (mean and standard
diviations) of K
trans
and vp for direct and indirect methods at R=60× . The patient order is
sorted as the direct reconstruction performance degraded (same as figure 6). At a
significance level of p<0.05 (for individual case, p<0.0038 after Bonferroni correction),
the direct method performed better than the indirect method for all the cases, with the cut
off undersampling rates varying between 20× to 80× . The pvalue for the cutoff
undersampling rate is shown in the last column.
No.
Age
/Sex
Diagnosis
Indirect 60x Direct 60x
Significant
different at
R>
pvalue
K
trans
rMSE
(× 10e4)
v p rMSE
(× 10e4)
K
trans
rMSE
(× 10e4)
v p rMSE
(× 10e4)
1 75M Meningioma 117± 2.1 237± 4.1 96± 2.0 219± 6.3 20x
3.72× 10

7
2 74M Glioblastoma 107± 2.0 250± 6.2 76± 1.2 188± 4.2 20x 6.3× 10
9
3 73M Glioma 120± 3.8 305± 9.5 99± 3.2 258± 9.4 20x
2.55× 10

8
4 69F Meningioma 178± 4.5 313± 8.3 103± 3.4 201± 6.2 20x
4.1× 10

10
5 77M Glioma 90± 4.7 191± 8.3 82± 5.2 178± 7.4 20x
6.86× 10

6
6 39F Meningioma 114± 5.0 192± 10.5 83± 4.7 191± 8.2 20x
1.75× 10

6
7 54F Glioma 167± 5.1 358± 14.5 156± 4.4 376± 17.3 40x
2.41× 10

6
8 44F Meningioma 155± 6.4 327± 13.1 123± 5.6 281± 14.3 40x
6.57× 10

6
9 60M Glioblastoma 130± 5.0 306± 10.2 113± 5.5 261± 10.1 40x
9.73× 10

7
10 38F Glioma 134± 5.4 457± 6.8 129± 5.3 424± 6.2 60x
7.42× 10

5
11 63M Meningioma 135± 5.4 361± 10.6 136± 6.9 352± 22.0 80x
1.25× 10

6
12 73F Glioma 171± 6.1 436± 18.2 171± 6.6 433± 16.5 80x
3.84× 10

4
13 79F Glioma 133± 9.5 318± 13.0 137± 4.7 316± 12.2 80x
1.31× 10

4
Figure 4.7 illustrates direct reconstruction of K
trans
and vp in two representative
prospectively undersampled wholebrain DCE data sets. Panning videos TK maps for all 4
cases are provided in supporting materials. The wholebrain highresolution TK maps
enable visualization of the tumor on any arbitrary reformatted plane, providing a complete
61
depiction of the pathological information, and evaluation of narrow enhancing margin and
small lesions. The reconstruction time was approximately 10 hours. This pilot study
demonstrates the feasibility of applying direct TK parameter reconstruction to wholebrain
DCEMRI.
Figure 4.7 Direct reconstruction of K
trans
((a), (c)) and vp ((b), (d)) maps from 2
representative prospectively undersampled data. Although lacking gold standard for the
prospective studies, the direct reconstruction provided reasonable K
trans
values and
complete depiction of the entire tumor region. Panningvolume videos for all 4 cases are
available in supplementary materials.
62
Discussion
We have presented a novel and potentially powerful TK parameter estimation scheme
for DCEMRI, where the TK parameter maps are directly reconstructed from
undersampled (k,t)space. By integrating the full forward model connecting the TK maps
to the (k,t)space data in the reconstruction, this method is able to provide excellent TK
map fidelity at undersampling rates up to 100× . Higher rates were not tested. The forward
model contained the analytic TK model along with specification of the AIF, coil sensitivity
maps, and precontrast T1, and M0 maps obtained from prescans. The optimization has the
flexibility to incorporate additional spatial sparsity constraints on the TK maps, as
demonstrated with a spatial wavelet transform. In the retrospective study, this method
outperformed an indirect reconstruction using parallel imaging and compressed sensing.
We also uniquely demonstrated the use of this method for prospectively undersampled
wholebrain DCE data, where wholebrain TK parameter maps can be produced with
excellent image quality.
The proposed method is a “parameterfree” reconstruction, when no spatial constraints
are applied to the TK maps. We demonstrated that by simply enforcing the TK model
during reconstruction, performance is improved relative to a stateoftheart compressed
sensing reconstruction, without the need to select a constraint or to tune associated
regularization parameters. It is straightforward to add sparsity constraints to the
optimization problem as shown in Equation (4.7). Such constraints improve the noise
performance, but at the expense of tuning regularization parameters. These constraints
were found to improve convergence at very high undersampling rates (>50× ) where the
TK map estimation problem becomes illposed.
The prospective study demonstrates that the proposed method can be used to achieve
substantially higher spatial resolution and broader spatial coverage DCEMRI, while
maintaining the same temporal resolution and overall scan time. Although not studied in
this work, this approach could potentially be used to improve the temporal resolution of
DCEMRI which is known to provide improvements in patientspecific AIF measurement
and TK parameter precision (110).
The proposed method estimates the TK maps directly, and rather than reconstructing
image time series as an intermediate step. This enables robust parameter estimation and
easyofuse in clinical application. Clinically, intermediate images (typically 50100
volumes) are not always viewed. The extracted TK maps are of primary interest as they
succinctly describe the behavior of the intermediate images. It is worth noting that the
proposed method can provide intermediate images using the full forward model described
in Figure 4.1. Figure 4.8 compares synthesized images from reconstructed TK maps against
63
fullysampled anatomic images. The synthesized images show close resemblance to the
fullysampled images.
Figure 4.8 Illustration of intermediate anatomic images from fullysampled (a),
synthesized direct reconstruction at R=30× (b). Images from pre, peak, post, and last
contrast arrival time points are shown from left to right. Anatomic images synthesized from
direct estimation of TK parameter maps show similar quality to the fullysampled and
reconstructed anatomic images.
The proposed direct reconstruction scheme requires apriori definition of the arterial
input function (AIF). In this work, we used a populationbased AIF, and the timeofarrival
was automatically detected as described in (105). However, other extensions which are
blind to the choice of AIF could be explored. For instance, Fluckiger et al. (32) proposed a
modelbased blind estimation of both AIF and TK parameters from DCE images for fully
sampled data. This approach may be combined with ours to provide joint reconstruction of
AIF and TK parameters directly from undersampled data.
This study has a few important limitations. First, we have thus far only demonstrated
effectiveness of this approach using the Patlak model. Patlak was chosen because it is
widely used and can be linearized and gradients can be readily computed. We also
restricted the retrospective study to datasets that fit the Patlak model. It will be important
to develop support for more sophisticated models, and utilize data that does not fit the
presumed model to fully characterize failure modes. Use of more sophisticated models (e.g.
extended Tofts model or 2compartment exchange model) may fit the data better. Their
64
inclusion will make the reconstruction problem nonlinear, and possibly nonconvex.
Gradient descent algorithms may not be applicable, as they require analytic solution for the
first derivative of CA(r,t) with respect to each TK parameter in the model (step 1 for Patlak
model in Figure 1). Dikaois et al. (111) recently demonstrated the use of a Bayesian
formulation of direct TK parameter estimation in DCE. The rationale was to use an optimal
model for different tissue types. The additional complexity of more sophisticated models
will necessitate longer reconstruction time, and convergence will require further
investigation, and remains as future work.
A second limitation is that the tight integration of TK modeling in our reconstruction
could be sensitive to data inconsistencies, such as patient motion. This is equally true for
indirect reconstruction. Prospective motion compensation could be added to the proposed
model but the complexities involved and efficacy have not been investigated here.
A third limitation is that the intermediate anatomic images are computed during the
reconstruction and thus require a similar amount of memory and computation time as
indirect methods. This approach is not currently solving the computational limitation of
constrained reconstructions but rather provides a framework for improved image quality.
A fourth limitation is that T1 and M0 maps are required in the forward model, and were
estimated using a separate multiple flipangle sequence (DESPOT1) performed
immediately prior to the DCE scan. Future work could include joint estimation of the pre
contrast T1 maps and the TK maps, as suggested by Dickie et al. (112).
We compared the direct method with a stateoftheart indirect method that utilized a
temporal finite difference constraint. Compressed sensing techniques are expected to
improve steadily as better constraints are identified (113–115); however, if used for TK
parameter estimation, a TK model will be applied to the data and model inconsistencies
introduced by the intermediate sparsity transforms are likely to propagate into the final
parameter maps. Although untested, we hypothesize that our direct estimation method is
likely to meet or exceed the image quality of any compressed sensing method.
Our proposed direct reconstruction scheme provides a method for highly accelerated
DCE. Extremely high acceleration rates have been demonstrated (up to 100× ), enabling
full brain DCE with high spatial and temporal resolution. Our method is able to provide a
parameterfree reconstruction and so avoids the empirical tuning required in other methods.
This technique is easily extendable to DCEMRI in other body parts such as the breasts,
prostate, etc. Future work will include exploration of these additional clinical applications,
optimization of data sampling schemes, and integration of more sophisticated TK models.
65
Conclusion
We have presented a novel and efficient reconstruction scheme to directly estimate TK
parameter maps from highly undersampled DCEMRI data. By comparison with a state
oftheart indirect compressed sensing method, we demonstrate that the proposed direct
approach provides improved TK map fidelity, and enables much higher acceleration. With
the prospective study, this method is shown to be clinically feasible and provide high
quality wholebrain TK maps.
Appendix 4A
The optimization problem in Equation (4.5) is solved alternatively by a quasiNewton
limitedmemory BroydenFletcherGoldfarbShanno (lBFGS) method. That is, solving
one while keeping the other fixed, as the pseudo codes indicated below:
Input initial guess
(0)
()
trans
K r
(0)
()
p
v r
k=0, while “stopping criteria not met” do {
( 1) ( ) 2
2
()
( ) arg min  ( , , ) ( ( ), ( ) ) 
trans
trans k rans k
p
K
K S t c y K v
r
r k r r
( 1) ( 1) 2
2
()
( ) argmin  ( , , ) ( ( ) , ( )) 
p
k trans k
pp
v
v S t c y K v
r
r k r r
k=k+1; }
The gradient of the cost function is evaluated analytically. This can be derived from the
model and signal equations. For notational simplicity, the coordinate notations r, k and c
are neglected (t is preserved to show the difference in dimension between TK parameter
maps and dynamic images). For example, is for ( ) and ( ) is for ( , , )
trans trans
K K S t S t c rk .
In Equation (4.5) we denote the cost function as:
2
2
( , )  ( ) ( , ) 
trans trans
pp
y K v S t f K v
(4.A1)
where for one iteration, vp and all other known variables are kept constant, and we focus
on deriving the gradient of y w.r.t. K
trans
. We use the derivative chain rule:
66
( ) ( )
( ) ( )
trans trans
y CA t s t y
K K CA t s t
(4.A2)
, where
[ ( ) ( , )]
()
H H trans
u u p
y
C F k t f K v
st
(4.A3)
1
()
1 1 1 1
1 0 1 2
1
( ) (1 cos ) (1 ) cos
sin ,
( ) (1 cos )
TR R t
s t TR E E E TR E
r M E e
CA t E
(4.A4)
0
0
()
()
L
t tt
p trans
t
CA t
Cd
K
(4.A5)
Sparsitybased constraints can be optionally applied to the TK maps as shown in
Equation (4.7). In this study, we demonstrate the use of wavelet transform, and we denote
the wavelet constrained part as y1=ψx1. For the evaluation of y1, the l1 norm is relaxed
as in Ref (44):
1 1
1
()
trans
H trans
trans
yK
WK
K
(4.A6)
And the i
th
diagonal element of W is calculated as:
*
( ( ) ( ) )
trans trans
i i i
W K K
(4.A7)
where μ is a small relaxation parameter.
The gradient w.r.t. vp is very similar:
( ) ( )
( ) ( )
pp
y CA t s t y
v v CA t s t
(4.A8)
where all other parts are the same as above except:
()
()
p
p
CA t
Ct
v
(4.A9)
67
Chapter 5. Joint Arterial Input Function and Tracker Kinetic
Parameter Estimation from UnderSampled DCEMRI using a
Model Consistency Constraint
Introduction
Dynamic contrast enhanced (DCE) MRI is a powerful technique for probing subvoxel
vascular properties of tissue including fractional plasma volume, fractional extracellular
extravascular volume, and clinically important transfer constants. DCEMRI involves
capturing a series of images before, during, and after administration of a T1shortening
contrast agent. Tracerkinetic (TK) parameter maps are then computed from the dynamic
images to provide information for diagnosis and monitoring treatment response (8, 22, 116).
DCEMRI is used throughout the body, however, in the brain, it has shown value in the
assessment of brain tumor, multiple sclerosis, and Alzheimer disease (10, 26, 61).
With conventional Nyquist sampling, DCEMRI is often unable to simultaneously
provide adequate spatiotemporal resolution and spatial coverage. A typical brain DCE
MRI uses 510 second temporal resolution, which is a minimum requirement for accurate
TK modeling (25, 35). Using Cartesian sampling at the Nyquist rate, only 510 slices are
achievable. This is typically inadequate in large Glioblastoma cases, and cases with
scattered metastatic disease that are often spread throughout the brain (92). Many
techniques involving undersampling and constrained reconstruction have been proposed
to overcome this limitation. Earlier, compressed sensing and parallel imaging techniques
have been used to reconstruct dynamic images from undersampled data (49, 50, 56), and
wholebrain or highresolution TK maps can be fitted based on the reconstructed images
(92, 93). Another approach is to utilize the TK modeling information to directly estimate
TK parameters from undersampled data. Such modelbased reconstruction approaches
have been explored in MRI relaxometry (97, 117), PET kinetic parameter estimation (100,
102), and recently, in DCEMRI TK parameter estimation (103, 104, 118, 119). Compared
to conventional compressed sensing techniques that reconstruct dynamic images first, the
modelbased approach for DCEMRI reconstruction provides superior results and allows
higher undersampling rates (118, 119).
In conventional DCEMRI, images are reconstructed for each time point. Patient
specific arterial input functions (AIF) can be identified from vessel pixels using either
manual region of interest (ROI) selection or automatic clusterbased ROI selection (33).
Some centers use a fixed populationaveraged AIF (36), and institutionallyderived AIF,
or a delay and dispersion corrected version of these (92). The use of patientspecific AIF
(patAIF) is known to provide more accurate TK mapping (120). However, the estimation
68
of patAIF from undersampled data is severely challenged by undersampling artifacts.
The current modelbased TK reconstruction approaches therefore rely on population
averaged AIF (popAIF), which greatly limits its full potential.
In this work, we propose a novel DCEMRI reconstruction approach, in which TK maps
and dynamic images are simultaneously reconstructed, and TK model consistency is
applied as a reconstruction constraint. Furthermore, this method allows for extraction of
patientspecific AIFs, a key advantage compared to previous methods. This approach is
inspired by recent works in accelerated quantitative MR relaxometry (121, 122), where
physical or physiological model consistency was applied as a reconstruction constraint.
This consistency constraint allowed for the data fit to deviate from the model, which made
the scheme robust to scenarios with model inconsistencies (eg. presence of motion).
The DCEMRI model consistency constraint requires a regularization parameter that
balances the tradeoff between data consistency and model consistency. We show that the
formulation allows for easy integration of different TK models and/or different TK solvers.
The method is validated using an established physiologicallyrealistic brain tumor digital
reference object and retrospectively undersampled clinical braintumor DCEMRI data.
We also demonstrate its application to prospectively undersampled highresolution whole
brain DCEMRI data.
Theory
5.2.1 Model Consistency Constraint
We jointly estimate contrast concentration vs time images (C) and TK parameter maps
( ) from the undersampled data (y) by solving a leastsquares problem:
22
0 2 2
,
( , ) arg min  ( )   ( ) 
C
C UFE C S y P C
(5.1)
where Ψ is contrast concentration (C) to image difference (ΔS) conversion, S0 is the
fullysampled precontrast image, U is the undersampling mask, F is the Fourier transform,
E is the sensitivity encoding, and P represents the forward TK modeling (linear for Patlak,
nonlinear for all other TK models). This formulation can be simplified to:
where A=UFEΨ, b=(y+UFES0) is the known data.
To solve the leastsquare optimization problem in Equation (5.2), we alternatively solve
for each variable while keeping others constant. For each iteration n:
22
22
,
( , ) arg min    ( ) 
C
C AC b P C
(5.2)
69
1 2 2
22
arg min    ( ) 
nn
C
C AC b P C
(5.3)
1 1 1
()
nn
PC
(5.4)
Note that Equation (5.3) is regularized SENSE reconstruction with an l2 norm constraint
that can be solved efficiently using conjugate gradients (CG) (30). Equation (5.4) is
backward TK modeling that can be solved using any DCEMRI modeling toolbox. Because
forward modeling (P) and backward modeling (P
1
) are used iteratively, the modeling
solver should not utilize linearization or other forms of approximation. For example,
Rocketship (123) and TOPPCAT(124) are two suitable solvers. Detailed substeps and
variants of Equation (5.3) and (5.4) can be found in Appendix I.
5.2.2 Joint AIF and TK Parameter Estimation
The proposed formulation allows for joint estimation of the patientspecific AIF.
Equation (5.2) can be modified to estimate C, θ and AIF from undersampled data by
solving the following leastsquares problem:
22
22
,,
( , , ) arg min    ( , ) 
C AIF
C AIF AC b P AIF C
(5.5)
Similar to the above, we solve each variable alternatively as follows (n
th
iteration):
1 2 2
22
arg min    ( , ) 
n n n
C
C AC b P AIF C
(5.6)
1 1 1
, ( )
n n n
AIF P C
(5.7)
Equation (5.7) is backward TK modeling from contrast concentration including patAIF
estimation. This can be performed by identifying an arterial ROI once, using the time
averaged image or postcontrast image. Within each iteration, it is then possible to: 1) apply
this ROI to C to estimate the AIF (averaging the pixels), and 2) use the updated AIF during
TK modeling. This is a common procedure in TK modeling for DCEMRI. The only
difference is identification of the arterial ROI prior to the reconstruction of the dynamic
images.
5.2.3 Theoretical Benefits
The proposed method formulates model consistency as a constraint with a penalty β,
and decouples it from data consistency. There are multiple benefits of this formulation: 1.
algorithm complexity is greatly reduced compared to recently proposed direct
reconstruction techniques that require complex gradient evaluations (103, 111, 118); 2.
70
different TK models can easily be included in this formulation, as described above; 3.
patientspecific AIFs can be estimated jointly with TK maps, as described above; and 4.
the penalty β can allow for TK model deviation, reducing errors that may be caused by
strict model enforcement (122). This work specifically demonstrates #2 and #3.
Methods
5.3.1 Data Sources
Digital Reference Object: Anatomicallyrealistic brain tumor DCEMRI digital
reference object (DRO) was generated based on the method and data provided by Bosca
and Jackson (125). The ExtendedTofts (eTofts) model was used to generate contrast
concentration curves with known TK parameter maps and popAIF (36). Coil sensitivity
maps measured on our MRI scanner (3T, 8channel head coil) were coregistered to the
DRO and used to generate realistic MRI kspace data.
Retrospective: Nine fullysampled brain tumor DCEMRI data sets were acquired from
patients receiving routine brain MRI with contrast (including DCEMRI) at our Institution.
The study protocol was approved by our Institutional Review Board. The acquisition was
based on a 3D Cartesian fast spoiled gradient echo sequence (SPGR) with FOV:
22× 22× 4.2cm
3
, spatial resolution: 0.9× 1.3× 7.0mm
3
, temporal resolution: 5s, 50 time
frames, and 8 receiver coils. The flip angle was 15° , and TE was 1.3ms, TR was 6ms.
DESPOT1 was performed prior DCEMRI, with flip angle of 2° , 5° , 10° to estimate pre
contrast T1 and M0 maps. The contrast agent, Gadobenate dimeglumine (MultiHance
Bracco Inc., relaxivity r1=4.39 s1mM1 at 37⁰C at 3 Tesla (109)) was administered with a
dose of 0.05 mMol/kg, followed by a 20 ml saline flush in the left arm by intravenous
injection.
Prospective: Prospectively undersampled data were acquired in one brain tumor
patient (65 M, Glioblastoma) with Cartesian goldenangle radial kspace sampling (57, 92).
3D SPGR data was acquired continuously for 5 minutes. Wholebrain coverage was
achieved with a FOV of 22× 22× 20 cm
3
and spatial resolution of 0.9× 0.9× 1.9 mm
3
. The
prospective study protocol was approved by our Institutional Review Board. Written
informed consent was provided by the participant.
71
5.3.2 Demonstration of TK Solver Flexibility
To demonstrate TK solver flexibility, DRO data was retrospectively undersampled
using a randomized goldenangle sampling pattern at R=20× (57). The proposed method
with eTofts modeling was used to reconstruct TK parameter maps at R=1× (fullysampled)
and R=20× . An inhouse gradientbased algorithm and an opensource TK modelling
toolbox, Rocketship (123), were used for the eTofts solver in the proposed algorithm
(Equation (5.4)). Tumor ROI rMSE for the K
trans
maps were calculated for evaluation.
5.3.3 Demonstration of TK Model Flexibility
The nine fullysampled patient data were fitted to both Patlak and eTofts model to
calculate model fitting error, and Ftest was performed in the tumor ROI to determine
whether Patlak or eTofts model is an appropriate fit (39). If more than 50% of the tumor
pixels were appropriately fitted for certain model, this model was selected for the data set.
We reconstructed the corresponding TK parameter maps for fullysampled data (used as
reference) and at undersampling rates of 20× , 60× and 100× for all 10 cases. A randomized
goldenangle sampling pattern (57) was used in the kxky plane, simulating kykz phase
encoding in a 3D wholebrain acquisition. Images were reconstructed using a popAIF (36)
with patientspecific delay corrected by the delay estimated from kspace center (105).
ROIbased K
trans
rMSE and K
trans
histograms were calculated based on the reference K
trans
maps. K
trans
histogram skewness and 90% percentile K
trans
values were measured for
evaluation (66, 126).
5.3.4 Demonstration of Joint AIF and TK estimation
The cases following Patlak model were reviewed with special attention to vessel signal.
Cases that showed significant precontrast inflow enhancement were identified and
subsequently excluded. With the remaining cases, we performed joint estimation of AIF
and Patlak parameter maps from undersampled data across sampling rates of 20× , 60× ,
and 100× . For each undersampling rate, 15 realizations were generated by varying the
initial angle of the goldenangle radial sampling pattern (57). Reconstructed patient
specific AIFs were compared to the fullysampled reference using rMSE and bolus peak
difference. ROIbased K
trans
relative rMSE (normalized by reference 90% percentile K
trans
)
were also calculated for evaluation.
72
5.3.5 Demonstration with Prospectively Undersampled Data
We demonstrate application of the proposed method for joint AIF and TK parameter
estimation on prospectively 30x undersampled highresolution wholebrain DCEMRI data.
Five second temporal resolution was achieved by grouping raw (k,t)space data acquired
within consecutive 5 sec intervals, effectively 30× undersampling compared to Nyquist
sampling (127). Patientspecific AIF and TK maps were jointly reconstructed using the
proposed model consistency constraint approach. PatAIF ROI was selected based on time
averaged images. Threeplane of K
trans
and vp maps and patAIF are presented for visual
assessment.
Results
Figure 5.1 shows the reconstruction results of the DRO at R=1× and R=20× . Etofts
model was used to generate the simulated DCEMRI data, and also for modelbased
reconstruction. Empirical β value and iteration number were chosen. Computation time
for the conversion from concentration vs time to TK maps, was 3.44s for the inhouse
gradientbased method, and 31.62s for Rocketship with parallel computing turned on (4
workers). At R=1× , by comparing to the true TK values that are used to generate the
phantom, both gradientbased method and Rocketship were able to restore accurate TK
maps. AT R=20× , both methods were able to restore goodquality K
trans
and vp maps, with
slightly degraded Kep maps because of the nonlinearity of the Kep parameter. The proposed
method provides better rMSE performance for K
trans
maps. These results shows that the
proposed method is compatible with a thirdparty solver.
73
Figure 5.1 The proposed method is compatible with thirdparty TK solvers. Shown are
results using an inhouse gradientbased solver and the Rocketship solver, both using the
model consistency constraint method on an anatomicallyrealistic braintumor DCEMRI
digital reference object. At R=1x (fullysampled), both methods are able to restore TK
maps accurately. At R=20x, both methods were able to restore goodquality K
trans
and vp
maps, with slightly degraded Kep maps. The gradientbased method provides better results
in terms of K
trans
rMSE in less computation time.
Figure 5.2 illustrates the impact of regularization parameter β for one representative in
vivo brain tumor dataset, using the Patlak model, at R=20x. The cost function values as a
function of iteration number, lcurve, and the final reconstructed TK maps are plotted for
different β values. A large β resulted in slow convergence, while a smaller β provided fast
convergence. This behavior is expected as illconditioning of the problem in Equation (5.3)
increases with β (128). TK maps obtained with a large β show poor fidelity as data
consistency is violated, while the maps with a small β is equivalent to a SENSE
reconstruction without constraints, and demonstrated gfactor related artifacts at R=20× .
The Lcurve shows the balance between the data consistency and model consistency, based
on which the β values in the range of 0.1 to 1 (green highlighted) show similar performance.
We then tune the β value in this range for different cases. We found the acceptable range
74
to be roughly 1 order of magnitude, and consistent among the 4 cases that we carefully
examined. This suggests that the parameter could be standardized.
Figure 5.2 Performance for different β values at R=20× . (a) The lcurve shows that β value
controls the balance between model and data consistency. (b, c, d) Convergence of the cost
function to within 1% of its final value required 116, 24, 10, 4, and 2 iterations for β values
of 10, 1, 0.1, 0.01, and 0.001, respectively. (e) Tumor ROI rMSE (K
trans
) are 9.1, 6.5, 6.4,
8.7, 9.3 (x10
3
) respectively for β value of 10, 1, 0.1, 0.01, 0.001. Reconstruction with small
β values converged quickly, and is closer to a SENSE reconstruction with associated g
factor losses and undersampling artifacts. Reconstruction with large β values shows slow
convergence, and provides inaccurate TK maps since the data consistency is violated.
Based on the tumor ROI Ftest, the Patlak model was appropriate for 6 cases, while the
eTofts model was appropriate for 3 cases. Figure 3 shows two representative cases of
Patlak and eTofts model at R=60× and R=100× . K
trans
and vp maps on the zoomedin tumor
region are shown (Kep for eTofts is not shown). Based on the K
trans
rMSE and histograms,
the proposed method is able to restore good quality TK parameter maps at undersampling
rates up to 100× , including for the nonlinear eTofts model.
75
Figure 5.3 Reconstruction of the TK maps of two representative case for Patlak and eTofts
model at R=60× and 100× . Tumor ROI (indicated in the reference images) histogram were
shown below the respective cases. Quantitative evaluation across all the cases were shown
in figure 4. Judging by the histogram and rMSE, the proposed method was able to restore
goodquality TK maps at a high undersampling rate of 100× for both Patlak and eTofts
model.
Figure 5.4 shows quantitative evaluation of the reconstruction results focusing on K
trans
values. For Patlak model reconstruction, the 90 percentile K
trans
values match well with the
reference values across all cases, the histogram skewness are also reasonably matched, and
tumor ROI rMSE is almost 0.023 for all cases and undersampling rates. For eTofts mode,
the 90 percentile K
trans
matched well with reference for one case, and have larger deviation
for other cases at R=100× . The maximum rMSE value were at most 0.073.
76
Figure 5.4 Quantitative evaluation of Patlak (top row) and eTofts (bottow row)
reconstruction. 90 percentile of the reconstructed K
trans
values for different cases were
plotted against the reference 90 percentile K
trans
. For Patlak model the values matched well
for all cases and undersampling rates (a). For eTofts model the values matched well for
R=20× and 60× , and have larger deviation for R=100× (d). The K
trans
histogram skewness
were also plotted against the reference histogram skewness (b), (e). The tumor ROI K
trans
rMSE were plotted against different R’s across different cases, and rMSE is at most 0.023
for Patlak reconstruction (c) and 0.073 for eTofts reconstruction (f).
Figure 5.5 shows the selection of AIF ROI from undersampled data, and the
comparison of popAIF and patAIF, and the resulting TK maps. It is known that popAIF
may not represent patient specific characteristics, and the difference in popAIF and pat
AIF can lead to significant errors in the resulting TK maps.
77
Figure 5.5 Left: Extraction of patAIF (b) from a manually selected ROI on the peak
contrast frame of fullysampled images (a). The popAIF show in (b) was delay and
amplitude corrected. In undersampling scenario, a time averaged image can be generated
(c), and even at R=100× (d), it is straightforward to select an artery ROI from this image
for the joint AIF and TK maps reconstruction. Right: Different AIFs can result in different
TK maps (e, f, g, h), and patAIF is preferred for more accurate TK modeling.
Figure 5.6 shows the reconstruction results of TK maps and patAIF (same case in figure
5) at different undersampling rates. Comparing to the AIF extracted from fullysampled
data, the proposed method is able to restore accurate AIF up to R=100× , with goodquality
TK maps restored at the same time.
78
Figure 5.6 Joint reconstruction of TK maps (cropped portion of the case in Figure 5) and
AIF at R=20× , 60× and 100× . Comparing to the fullysampled reference, the proposed
method is able to restore both AIF and TK maps at the same time, even at a high under
sampling rate of 100× .
Figure 5.7 show the quantitative evaluation of joint AIF and TK reconstruction across
the 4 data sets. Based on the relative rMSE of the TK maps, goodquality TK maps can be
restored even at R=100x, with a relative rMSE less than 20%. Multiple noise realizations
show that the method is robust noise, with an expected increase in variance at higher under
sampling rates. The shape of the AIF can be estimated at up to R=100× , with AIF rMSE
below 0.08 mMol. The peak of the AIF shows larger variance for different noise realization,
since the peak is only one point. However, the proposed method is still able to restore the
AIF peak up to R=60× with the error at most 0.25mMol across all cases.
79
Figure 5.7 Quantitative evaluation of the joint AIF and TK reconstruction for the 4 cases
across R=20× , 60× and 100× . (a) K
trans
relative rMSE was calculated as the spatial rMSE
across all tumor pixels, divided by the 90 percentile of the reference tumor K
trans
value. (b)
AIF rMSE was calculated as the temporal rMSE with respect to the reference AIF. (c) AIF
peak error was calculated as the reference peak minus the estimated peak. Across different
ceases, the rMSE mean and variance all increased with undersampling rate, as expected.
Figure 5.8 shows reconstruction of patAIF and TK maps from prospectively under
sampled data. This demonstrates that wholebrain TK maps can be reconstructed jointly
with patientspecific AIF, with no obvious undersampling artifacts in the final TK maps.
The clinicallymeaningful benefits of undersampling can be best demonstrated in
prospective study, where arbitrary reformats of the 3D TK maps are made possible thanks
to the ability to achieve high spatial resolution and wholebrain coverage.
Figure 5.8 Joint reconstruction of PatAIF and TK maps from prospective undersampled
data. Wholebrain highresolution TK maps can be provided together with patientspecific
AIF using the proposed modelbased reconstruction approach.
80
Discussion
We have described, demonstrated, and evaluated a novel modelbased reconstruction
approach for DCEMRI, where the TK model is posed as a penalized consistency constraint.
By this formulation, we decoupled the TK model consistency from the (k,t)space data
consistency. The two subproblems can be solved using existing techniques, namely TK
modeling (including AIF estimation) and regularized SENSE reconstruction. The proposed
approach allows for easy inclusion of different TK solvers, including thirdparty solvers,
and also allows for joint estimation of the patientspecific AIF. We have demonstrated the
robustness of the proposed method in one anatomicallyrealistic braintumor DRO, and
retrospective study across nine brain tumor DCEMRI datasets. We also demonstrated the
application of the proposed method to prospectively undersampled data.
Limitations of the Study: The proposed method also has a few important limitations.
First, the alternating algorithm proposed is a twoloop iteration, where an iterative solver
is needed for each subproblem. Comparing to a gradientbased direct reconstruction (118),
this formulation takes longer computing time. This issue can be addressed by using
powerful computers, implementing in C, or using GPU acceleration.
Second, although we demonstrate that the proposed method is compatible with a third
party solver, it requires that the solver does not use any linearization or approximation for
the modeling. For complicated TK models, a few linearized approaches have been
proposed for fast computation (129, 130). Unfortunately, those methods cannot be
integrated into the proposed framework.
Third, although we have shown that this method can include different TK solver, it may
be difficult to use a nested model that selects between several different local model based
on local fitting errors (39). This type of approach has been shown in recent literature to be
advantageous. The quality of intermediate anatomic images in the proposed method,
especially in the first few iterations, may make it challenging to generate a modeling mask
needed for nested models
Conclusion
We have demonstrated a novel modelbased reconstruction approach for accelerated
DCEMRI. Posing the TK model as a model consistency constraint, this formulation
provides flexible use of different TK solvers, joint estimation of patientspecific AIF, and
straightforward implementation. Evaluated in digital reference object, retrospective and
81
prospective study of actual brain tumor patient data, this method provides good quality TK
maps and patAIF up to undersampling rates of 100x for both Patlak and eTofts model.
Appendix 5A
For the proposed method, we used an alternative approach to solve C and θ from under
sampled kspace, and this appendix explains the detailed steps for subproblems in
Equation (5.3) and Equation (5.4).
In Equation (5.3), we solve the contrast concentration from the measured data using the
following equation:
1 2 2
22
argmin    ( ) 
nn
C
C AC b P C
(5.A1)
Here A=UFEΨ, and we first solve the image difference (ΔS) from b (since we moved
S0 to b) by solving the following leastsquare problems using CG (or other iterative
algorithm for leastsquare problems). We can use result from previous iteration as initial
guess for faster convergence.
22
22
arg min  ( )   ( ) 
n
S
S UFE S b S P
(5.A2)
The first term is SENSE, and the second term is an identity constraint to ΨP(
n
) that is
constant in this step. P is the forward modeling from TK maps to contrast concentration C,
and Ψ is the conversion from contrast concentration C to signal difference ΔS following
the steadystate SPGR signal equation:
0 1 0
0 1 0
()
00
()
sin (1 ) sin (1 )
()
1 cos 1 cos
TR R C r TR R
TR R C r TR R
M e M e
SC
ee
(5.A3)
Where TR is the repetition time, α is the flip angle, r1 is the contrast agent relaxivity. R0
and M0 are the precontrast R1 (reciprocal of T1) and the equilibrium longitudinal
magnetization that are estimated from a T1 mapping sequence. In this work, we used
DESPOT1 (29) prior to the DCEMRI scan.
Note that Ψ is a onetoone mapping for each voxel, and its inversion (C= Ψ
1
(ΔS)) is:
82
Equation (5.A5) is used to compute C after solving for ΔS using Equation (5.A2), this
completes the detailed algorithm for solving Equation (5.3).
After C (C(t) used in this equation to avoid confusion) is estimated, Equation (5.4)
represents backward TK modelling. For Patlak model, it is:
0
( ) ( ) ( , ) ( )d ( )
t
trans trans
p p p p
C t P P K v K C v C t
(5.A6)
Where Cp(t) is the arterial input function (AIF). Patlak model is linear, and a pseudo
inverse can be used to solve =P
1
(C).
And for extendedTofts (eTofts) model, it is:
()
0
( ) ( ) ( , , ) ( ) d ( )
ep
t
Kt
trans trans
p ep p p p
C t P P K v K K C e v C t
(5.A7)
where an extra TK parameter Kep is modeled for better fitting. eTofts is nonlinear and
iterative algorithm can be used to solve this model fitting:
2
2
arg min  ( )  PC
(5.A8)
We use a gradientbased lBFGS algorithm to solve Equation (5.A8), where we derive
the gradient for each TK parameter. Opensource toolbox, Rocketship, was also used for
comparison.
0
0
0
0
0
0
01
1
1
sin 1 cos
1
ln
1
1 cos
sin 1 cos
(  ) /
TR R
TR R
t
TR R
TR R
t
Se
Me
R
TR
Se
Me
C R R r
(5.A5)
83
Chapter 6. Conclusion
In this thesis, we proposed several novel methods targeted to improve brain DCEMRI
towards highresolution, wholebrain coverage, and highfidelity TK maps. DCEMRI is a
valuable dynamic imaging technique that is able to provide powerful biomarker for
evaluation of brain tumor, helping diagnosis and therapy response. However, DCEMRI is
not currently the standardofcare protocol in many centers conducting oncology, due to its
low resolution, limited spatial coverage, and low reproducibility for kinetic parameter
mapping.
We have proposed a specially tailored constrained reconstruction technique for DCE
MRI in chapter 3. The sampling scheme, sparsifying transform constraints and the
reconstruction algorithm are all specially tailored for dynamic imaging to provide truthful
dynamic information from vastly undersampled data. We have conducted the first (to the
best of our knowledge) evaluation of prospectively undersampling reconstruction of DCE
MRI in brain tumor patients, and compared to conventional technique that is sampled at
Nyquist rate. With 30× undersampling rate, we were able to achieve a wholebrain high
resolution DCEMRI. In radiologists’ rating, the denoising properties of the sparsifying
constraints also provided improved dynamic image quality in terms of fine details and
signaltonoise ratio comparing to conventional imaging technique. The highresolution,
wholebrain coverage enabled us 3D imaging in any arbitrary plane, providing a complete
pathological depiction of some large Glioblastoma tumor, or scattered lesions that spread
the wholebrain like metastatic brain tumor.
In chapter 4, We proposed an innovative and rigorous technical development approach
in which DCEMRI acquisition and reconstruction are tailored from an estimation
theoretic point of view to create the most reproducible tracker kinetic (TK) parameter maps.
Unlike conventional approaches that optimize the quality of intermediate dynamic images,
our proposed methods fully integrated TK models with DCEMRI acquisition and
reconstruction. The TK maps were directly estimated from undersampled kspace data,
and the TK model was implicitly enforced during the reconstruction process. Comparing
to earlier indirect method using stateoftheart compressed sensing techniques, this direct
reconstruction technique is able to provide even higher undersampling rates and better
quality TK parameter maps. We have shown up to 100× undersampling rate in
retrospective study, a rate that is never reported in DCEMRI literature. This ultrahigh rate
can enable us to further improve the spatiotemporal resolution of DCEMRI, achieving
submillimeter isotropic resolution with close to 1 second temporal resolution. We have
84
also shown the feasibility of the direct reconstruction in prospective undersampled data,
where wholebrain highresolution TK maps can be reconstructed directly from 30× under
sampled kspace.
In chapter 5, we developed and evaluated a modelbased DCEMRI reconstruction
framework for joint arterial input function (AIF) and kinetic parameter estimation from
undersampled data. This method poses the TK model as a model consistency constraint,
enabling the inclusion of different TK solvers and the joint estimation of the arterial input
function (AIF) from highly undersampled data. Comparing to previously proposed
indirect constrained reconstruction, or direct modelbased reconstruction, this model
consistency constrained reconstruction combined the benefits of both methods, while
overcoming the limitations of both methods. By posing the TK model as a model
consistency constraint instead of directly forcing the model. The framework can allow for
model deviation and provide better quality TK maps. It decouples the direct reconstruction
problem into two welldefined subproblems, thus greatly reducing the complexity of the
algorithm. Benefiting from the intermediate images similar to indirect reconstruction, the
proposed method also allows for easy inclusion of different TK solver, and joint estimation
of AIF from undersampled data. In digital reference object study, the proposed method is
shown to be compatible with a thirdparty TK solver. In the retrospective study, the
proposed method is demonstrated to produce TK maps and patientspecific AIF with
minimal error at up to R=100× for both Patlak and eTofts model. In prospective study, the
proposed method was able to provide highresolution wholebrain TK maps and patient
specific AIF at R=30× .
The proposed methods in this dissertation are still far from providing a complete and
universal DCEMRI approaches, especially in clinical translation practice. Several
limitations and possible improvement exist to further facilitate the clinical translation of
improved DCEMRI, some of them were discussed in previous chapters. For example, a
rigorous reproducibility study is still lacking to prove the proposed method is able to
provide reproducible TK mapping in the same patients. For higher undersampling rate up
to 100× , prospective study was not performed to validate if higher spatial or temporal
resolution can truly be achieved. Also, the TK model selection is not perfect in the proposed
study, and a spatial information based model selection is still a challenging in the proposed
methods.
85
References
1. Centers for Disease Control and Prevention [http://www.cdc.gov/nchs/fastats/leading
causesofdeath.htm]
2. Brain Tumor: Diagnosis [http://www.cancer.net/cancertypes/braintumor/diagnosis]
3. Ricard D, Idbaih A, Ducray F, Lahutte M, HoangXuan K, Delattre JY: Primary brain
tumours in adults. Lancet 2012; 379:1984–1996.
4. Ellingson BM, Bendszus M, Boxerman J, et al.: Consensus recommendations for a
standardized Brain Tumor Imaging Protocol in clinical trials. Neuro Oncol 2015;
17(August):1188–1198.
5. Essig M, Nguyen TB, Shiroishi MS, et al.: Perfusion MRI: the five most frequently
asked clinical questions. AJR Am J Roentgenol 2013; 201:W495510.
6. Wen PY, Macdonald DR, Reardon D a., et al.: Updated response assessment criteria for
highgrade gliomas: Response assessment in neurooncology working group. J Clin
Oncol 2010; 28:1963–1972.
7. Lin NU, Lee EQ, Aoyama H, et al.: Response assessment criteria for brain metastases:
proposal from the RANO group. Lancet Oncol 2015; 16:e270–e278.
8. Heye AK, Culling RD, Herná ndez CV, Thrippleton MJ, Wardlaw JM: Assessment of
blood – brain barrier disruption using dynamic contrastenhanced MRI . A systematic
review. NeuroImage:Clinical 2014; 6:262–274.
9. Tofts PS, Brix G, Buckley DL, et al.: Estimating kinetic parameters from dynamic
contrastenhanced T(1)weighted MRI of a diffusable tracer: standardized quantities and
symbols. J Magn Reson Imaging 1999; 10:223–32.
10. Law M, Yang S, Babb JS, et al.: Comparison of cerebral blood volume and vascular
permeability from dynamic susceptibility contrastenhanced perfusion MR imaging
with glioma grade. Am J Neuroradiol 2004; 25:746–755.
11. Larsson C, Kleppestø M, Grothe I, Vardal J, Bjø rnerud A: T1 in highgrade glioma and
the influence of different measurement strategies on parameter estimations in DCEMRI.
J Magn Reson Imaging 2014; 0.
12. DCE MRI Technical Committee. DCE MRI Quantification Profile, Quantitative
Imaging Biomarkers Alliance. Version 1.0. Reviewed Draft. QIBA, July 1, 2012.
[http://rsna.org/QIBA_.aspx]
13. Weishaupt D, Kö chli VD, Marincek B: How Does MRI Work? 2008.
14. Stanisz GJ, Odrobina EE, Pun J, et al.: T1, T2 relaxation and magnetization transfer in
tissue at 3T. Magn Reson Med 2005; 54:507–12.
15. Wansapura JP, Holland SK, Dunn RS, Ball WS: NMR relaxation times in the human
brain at 3.0 tesla. J Magn Reson Imaging 1999; 9:531–538.
16. Chen L, Bernstein M: Measurements of T1 relaxation times at 3.0 T: implications for
clinical MRA. Proc 9th … 2001; 9:2001.
17. Pruessmann KP, Weiger M, Scheidegger MB, Boesiger P: SENSE: sensitivity
encoding for fast MRI. Magn Reson Med 1999; 42:952–962.
18. Griswold MA, Jakob PM, Heidemann RM, et al.: Generalized autocalibrating partially
parallel acquisitions (GRAPPA). Magn Reson Med 2002; 47:1202–1210.
19. Lustig M, Pauly JM: SPIRiT: Iterative selfconsistent parallel imaging reconstruction
from arbitrary kspace. Magn Reson Med 2010; 64:457–71.
86
20. Uecker M, Lai P, Murphy MJ, et al.: ESPIRiTan eigenvalue approach to
autocalibrating parallel MRI: Where SENSE meets GRAPPA. Magn Reson Med 2014;
71:990–1001.
21. Haldar JP: Lowrank modeling of local kspace neighborhoods (LORAKS) for
constrained MRI. IEEE Trans Med Imaging 2014; 33:668–81.
22. Tofts PS, Kermode AG: Measurement of the bloodbrain barrier permeability and
leakage space using dynamic MR imaging. 1. Fundamental concepts. Magn Reson Med
1991; 17:357–67.
23. Cramer SP, Larsson HBW: Accurate determination of bloodbrain barrier permeability
using dynamic contrastenhanced T1weighted MRI: a simulation and in vivo study on
healthy subjects and multiple sclerosis patients. J Cereb blood flow Metab 2014;
34:1655–1665.
24. Yang S, Law M, Zagzag D, et al.: Dynamic contrastenhanced perfusion MR imaging
measurements of endothelial permeability: differentiation between atypical and typical
meningiomas. AJNR Am J Neuroradiol 2003; 24:1554–9.
25. Cramer SP, Simonsen H, Frederiksen JL, Rostrup E, Larsson HBW: Abnormal blood
brain barrier permeability in normal appearing white matter in multiple sclerosis
investigated by MRI. NeuroImage Clin 2014; 4:182–189.
26. Montagne A, Barnes SR, Law M, et al.: BloodBrain Barrier Breakdown in the Aging
Human Report BloodBrain Barrier Breakdown in the Aging Human Hippocampus.
Neuron 2015; 85:296–302.
27. Cunningham CH, Pauly JM, Nayak KS: Saturated doubleangle method for rapid B1+
mapping. Magn Reson Med 2006; 55:1326–1333.
28. Sacolick LI, Wiesinger F, Hancu I, Vogel MW: B 1 mapping by Bloch
Siegert shift. Magn Reson Med 2010; 63:1315–1322.
29. Deoni SCL, Peters TM, Rutt BK: Highresolution T1 and T2 mapping of the brain in a
clinically acceptable time with DESPOT1 and DESPOT2. Magn Reson Med 2005;
53:237–41.
30. Pruessmann KP, Weiger M, Bö rnert P, Boesiger P: Advances in sensitivity encoding
with arbitrary kspace trajectories. Magn Reson Med 2001; 46:638–651.
31. Li KL, Zhu XP, Waterton J, Jackson a: Improved 3D quantitative mapping of blood
volume and endothelial permeability in brain tumors. J Magn Reson Imaging 2000;
12:347–57.
32. Fluckiger JU, Schabel MC, DiBella EVR: Modelbased blind estimation of kinetic
parameters in Dynamic Contrast Enhanced (DCE)MRI. Magn Reson Med 2009;
62:1477–1486.
33. Shi L, Wang D, Liu W, et al.: Automatic detection of arterial input function in dynamic
contrast enhanced MRI based on affinity propagation clustering. J Magn Reson Imaging
2014:1327–37.
34. Van Osch MJP, Vonken EJP a, Bakker CJG, Viergever M a.: Correcting partial volume
artifacts of the arterial input function in quantitative cerebral perfusion MRI. Magn
Reson Med 2001; 45:477–485.
35. Henderson E, Rutt BK, Lee TY: Temporal sampling requirements for the tracer kinetics
modeling of breast disease. Magn Reson Imaging 1998; 16:1057–1073.
36. Parker GJM, Roberts C, Macdonald A, et al.: Experimentallyderived functional form
for a populationaveraged hightemporalresolution arterial input function for dynamic
87
contrastenhanced MRI. Magn Reson Med 2006; 56:993–1000.
37. Calamante F, Gadian DG, Connelly A: Delay and dispersion effects in dynamic
susceptibility contrast MRI: simulations using singular value decomposition. Magn
Reson Med 2000; 44:466–73.
38. Shiroishi MS, Habibi M, Rajderkar D, et al.: Perfusion and permeability MR imaging
of gliomas. Technol Cancer Res Treat 2011; 10:59–71.
39. Bagherebadian H, Jain R, Nejaddavarani SP, et al.: Model Selection for DCET1
Studies in Glioblastoma. Magn Reson Med 2012:241–251.
40. Sourbron SP, Buckley DL: Classic models for dynamic contrastenhanced MRI. NMR
Biomed 2013; 26:1004–1027.
41. Sourbron S, Ingrisch M, Siefert A, Reiser M, Herrmann K: Quantification of cerebral
blood flow, cerebral blood volume, and bloodbrainbarrier leakage with DCEMRI.
Magn Reson Med 2009; 62:205–217.
42. Sourbron SP, Buckley DL: On the scope and interpretation of the Tofts models for
DCEMRI. Magn Reson Med 2011; 66:735–745.
43. Sodickson DK, Manning WJ: Simultaneous acquisition of spatial harmonics (SMASH):
fast imaging with radiofrequency coil arrays. Magn Reson Med 1997; 38:591–603.
44. Lustig M, Donoho D, Pauly JM: Sparse MRI: The application of compressed sensing
for rapid MR imaging. Magn Reson Med 2007; 58:1182–1195.
45. Awate SP, DiBella EVR: Spatiotemporal dictionary learning for undersampled
dynamic MRI reconstruction via joint framebased and dictionarybased sparsity. In
Biomed Imaging; 2012:318–321.
46. Caballero J, Price AN, Rueckert D, Hajnal J V: Dictionary learning and time sparsity
for dynamic MR data reconstruction. IEEE Trans Med Imaging 2014; 33:979–94.
47. Gamper U, Boesiger P, Kozerke S: Compressed sensing in dynamic MRI. Magn Reson
Med 2008; 59:365–373.
48. Jung H, Sung K, Nayak KS, Kim EY, Ye JC: kt FOCUSS: a general compressed
sensing framework for high resolution dynamic MRI. Magn Reson Med 2009; 61:103–
16.
49. Feng L, Grimm R, Block KT, et al.: GoldenAngle Radial Sparse Parallel MRI :
Combination of Compressed Sensing , Parallel Imaging , and GoldenAngle Radial
Sampling for Fast and Flexible Dynamic Volumetric MRI. Magn Reson Med 2014;
72:707–717.
50. Zhang T, Cheng JY, Potnick AG, et al.: Fast pediatric 3D freebreathing abdominal
dynamic contrast enhanced MRI with high spatiotemporal resolution. J Magn Reson
Imaging 2015; 41:460–473.
51. Wang H, Miao Y, Zhou K, et al.: Feasibility of high temporal resolution breast DCE
MRI using compressed sensing theory. Med Phys 2010; 37:4971.
52. Rosenkrantz AB, Geppert C, Grimm R, et al.: Dynamic contrastenhanced MRI of the
prostate with high spatiotemporal resolution using compressed sensing, parallel imaging,
and continuous goldenangle radial sampling: Preliminary experience. J Magn Reson
Imaging 2015; 41:1365–1373.
53. Winkelmann S, Schaeffter T, Koehler T, Eggers H, Doessel O: An optimal radial
profile order based on the Golden Ratio for timeresolved MRI. IEEE Trans Med
Imaging 2007; 26:68–76.
54. Yeh EN, Stuber M, McKenzie C a, et al.: Inherently selfcalibrating nonCartesian
88
parallel imaging. Magn Reson Med 2005; 54:1–8.
55. Lustig M, Donoho D, Pauly JM: Sparse MRI: The application of compressed sensing
for rapid MR imaging. Magn Reson Med 2007; 58:1182–95.
56. Lebel RM, Jones J, Ferre JC, Law M, Nayak KS: Highly accelerated dynamic contrast
enhanced imaging. Magn Reson Med 2014; 71:635–644.
57. Zhu Y, Guo Y, Lingala SG, Lebel RM, Law M, Nayak KS: GOCART: GOldenangle
CArtesian randomized timeresolved 3D MRI. Magn Reson Imaging 2016; 34:940–950.
58. Lin FH, Wang FN, Ahlfors SP, Hä mä lä inen MS, Belliveau JW: Parallel MRI
reconstruction using variance partitioning regularization. Magn Reson Med 2007;
58:735–744.
59. Lustig M, Donoho DL, Santos JM, Pauly JM: Compressed Sensing MRI. (March
2008):72–82.
60. O’Connor JPB, Jackson A, Parker GJM, Roberts C, Jayson GC: Dynamic contrast
enhanced MRI in clinical trials of antivascular therapies. Nat Rev Clin Oncol 2012;
9:167–77.
61. Larsson HB, Stubgaard M, Frederiksen JL, Jensen M, Henriksen O, Paulson OB:
Quantitation of bloodbrain barrier defect by magnetic resonance imaging and
gadoliniumDTPA in patients with multiple sclerosis and brain tumors. Magn Reson
Med 1990; 16:117–131.
62. Khouli RH El, Macura KJ, Jacobs MA, et al.: Dynamic ContrastEnhanced MRI of the
Breast: Quantitative Method for Kinetic Curve Type Assessment. Am J Roentgenol 2009;
193:295–300.
63. Rosenkrantz AB, Lim RP, Haghighi M, Somberg MB, Babb JS, Taneja SS:
Comparison of interreader reproducibility of the prostate imaging reporting and data
system and likert scales for evaluation of multiparametric prostate MRI. AJR Am J
Roentgenol 2013; 201:W6128.
64. Paldino MJ, Barboriak DP: Fundamentals of Quantitative Dynamic ContrastEnhanced
MR Imaging. Magn Reson Imaging Clin N Am 2009:277–289.
65. Miller JC, Pien HH, Sahani D, Sorensen AG, Thrall JH: Imaging angiogenesis:
applications and potential for drug development. J Natl Cancer Inst 2005; 97:172–187.
66. Thomas AA, ArevaloPerez J, Kaley T, et al.: Dynamic contrast enhanced T1 MRI
perfusion differentiates pseudoprogression from recurrent glioblastoma. J Neurooncol
2015; 125:183–190.
67. Barger A V, Block WF, Toropov Y, Grist TM, Mistretta CA: Timeresolved contrast
enhanced imaging with isotropic resolution and broad coverage using an undersampled
3D projection trajectory. Magn Reson Med 2002; 48:297–305.
68. Haider CR, Hu HH, Campeau NG, Huston J, Riederer SJ: 3D high temporal and spatial
resolution contrastenhanced MR angiography of the whole brain. Magn Reson Med
2008; 60:749–760.
69. Trzasko JD, Haider CR, Borisch EA, et al.: SparseCAPR: highly accelerated 4D CE
MRA with parallel imaging and nonconvex compressive sensing. Magn Reson Med
2011; 66:1019–32.
70. Lee GR, Seiberlich N, Sunshine JL, Carroll TJ, Griswold MA: Rapid timeresolved
magnetic resonance angiography via a multiecho radial trajectory and GraDeS
reconstruction. Magn Reson Med 2013; 69:346–359.
71. Doneva M, Stehning C, Nehrke K, Bö rnert P: Improving Scan Efficiency of
89
Respiratory Gated Imaging Using Compressed Sensing with 3D Cartesian Golden
Angle Sampling. In ISMRM. Volume 19; 2011:641.
72. Zhu Y, Guo Y, Lebel RM, Law M, Nayak KS: Randomized Golden Ratio Sampling
For Highly Accelerated Dynamic Imaging. In ISMRM; 2014:4365.
73. Deoni SCL: Highresolution T1 mapping of the brain at 3T with driven equilibrium
single pulse observation of T1 with highspeed incorporation of RF field
inhomogeneities (DESPOT1HIFI). J Magn Reson Imaging 2007; 26:1106–11.
74. Ramani S, Fessler JA: Parallel MR image reconstruction using augmented Lagrangian
methods. IEEE Trans Med Imaging 2011; 30:694–706.
75. Afonso M V, BioucasDias JM, Figueiredo MAT: Fast image recovery using variable
splitting and constrained optimization. IEEE Trans Image Process 2010; 19:2345–56.
76. Afonso M V, BioucasDias JM, Figueiredo MAT: An augmented Lagrangian approach
to the constrained optimization formulation of imaging inverse problems. IEEE Trans
Image Process 2011; 20:681–95.
77. Hansen PC, Dianne Prost O’Leary: The use of the lcurve in the regularization of
discrete illposed problems. Soc Ind Appl Math 1993; 14:1487–1503.
78. Sourbron SP, Buckley DL: On the scope and interpretation of the Tofts models for
DCEMRI. Magn Reson Med 2011; 66:735–45.
79. Ewing JR, Knight RA, Nagaraja TN, et al.: Patlak plots of GdDTPA MRI data yield
bloodbrain transfer constants concordant with those of 14Csucrose in areas of blood
brain opening. Magn Reson Med 2003; 50:283–292.
80. Yun TJ, Park CK, Kim TM, et al.: Glioblastoma treated with concurrent radiation
therapy and temozolomide chemotherapy: differentiation of true progression from
pseudoprogression with quantitative dynamic contrastenhanced MR imaging.
Radiology 2015; 274:830–840.
81. Kickingereder P, Wiestler B, Graf M, et al.: Evaluation of dynamic contrastenhanced
MRI derived microvascular permeability in recurrent glioblastoma treated with
bevacizumab. J Neurooncol 2015; 121:373–380.
82. Bonekamp D, Deike K, Wiestler B, et al.: Association of overall survival in patients
with newly diagnosed glioblastoma with contrastenhanced perfusion MRI: Comparison
of intraindividually matched T1  and T2 (*) based bolus techniques. J Magn Reson
Imaging 2014; 42:87–96.
83. Chang SG, Yu B, Vetterli M: Adaptive wavelet thresholding for image denoising and
compression. IEEE Trans Image Process 2000; 9:1532–1546.
84. Elad M, Aharon M: Image denoising via sparse and redundant representations over
learned dictionaries. IEEE Trans Image Process 2006; 15:3736–45.
85. Chan SH, Khoshabeh R, Gibson KB, Gill PE, Nguyen TQ: An Augmented Lagrangian
Method for Total Variation Video Restoration. IEEE Trans Image Process 2011;
20:3097–3111.
86. Zhang T, Pauly JM, Vasanawala SS, Lustig M: Coil compression for accelerated
imaging with Cartesian sampling. Magn Reson Med 2013; 69:571–82.
87. Gai J, Obeid N, Holtrop JL, et al.: More IMPATIENT: A GriddingAccelerated
Toeplitzbased Strategy for NonCartesian HighResolution 3D MRI on GPUs. J
Parallel Distrib Comput 2013; 73:686–697.
88. Frahm J, Schä tz S, Untenberger M, et al.: On the Temporal Fidelity of Nonlinear
Inverse Reconstructions for Real Time MRI – The Motion Challenge. Open Med
90
Imaging J 2014; 8:1–7.
89. Wech T, Stab D, Budich JC, et al.: Resolution evaluation of MR images reconstructed
by iterative thresholding algorithms for compressed sensing. Med Phys 2012; 39:4328–
38.
90. Zhu Y, Guo Y, Lingala SG, et al.: Evaluation of DCEMRI data sampling,
reconstruction and model fitting using digital brain phantom. In Proc Int Soc Magn
Reson Med; 2015:3070.
91. Bosca RJ, Jackson EF: An extensible methodology for creating realistic
anthropomorphic digital phantoms for quantitative imaging algorithm comparisons and
validation. In ISMRM; 2015:797.
92. Guo Y, Lebel RM, Zhu Y, et al.: Highresolution wholebrain DCEMRI using
constrained reconstruction: Prospective clinical evaluation in brain tumor patients. Med
Phys 2016; 43:2013–2023.
93. Chandarana H, Feng L, Ream J, et al.: Respiratory MotionResolved Compressed
Sensing Reconstruction of FreeBreathing Radial Acquisition for Dynamic Liver
Magnetic Resonance Imaging. Invest Radiol 2015; 50:749–756.
94. Haldar JP, Hernando D, Liang Z: Superresolution Reconstruction of MR Image
Sequences with Contrast Modeling. IEEE Int Symp Biomed Imaging From Nano to
Macro 2009:266–269.
95. Ma D, Gulani V, Seiberlich N, et al.: Magnetic resonance fingerprinting. Nature 2013;
495:187–92.
96. Welsh CL, Dibella EVR, Adluru G, Hsu EW: Modelbased reconstruction of
undersampled diffusion tensor kspace data. Magn Reson Med 2013; 70:429–440.
97. Sumpf TJ, Uecker M, Boretius S, Frahm J: Modelbased nonlinear inverse
reconstruction for T2 mapping using highly undersampled spinecho MRI. J Magn
Reson Imaging 2011; 34:420–8.
98. Zhao B, Lam F, Liang ZP: Modelbased MR parameter mapping with sparsity
constraints: Parameter estimation and performance bounds. IEEE Trans Med Imaging
2014; 33:1832–1844.
99. Peng X, Liu X, Zheng H, Liang D: Exploiting parameter sparsity in modelbased
reconstruction to accelerate proton density and T2 mapping. Med Eng Phys 2014;
36:1428–1435.
100. Wang G, Qi J: Direct estimation of kinetic parametric images for dynamic PET.
Theranostics 2013:802–815.
101. Kamasak ME, Bouman C a, Morris ED, Sauer K: Direct reconstruction of kinetic
parameter images from dynamic PET data. IEEE Trans Med Imaging 2005; 24:636–50.
102. Lin Y, Haldar J, Li Q, Conti P, Leahy R: Sparsity Constrained Mixture Modeling for
the Estimation of Kinetic Parameters in Dynamic PET. IEEE Trans Med Imaging 2013;
33:173–185.
103. Felsted BK, Whitaker RT, Schabel M, DiBella EVR: Modelbased reconstruction for
undersampled dynamic contrastenhanced MRI. Proc SPIE 2009; 7262:1–10.
104. Dikaios N, Arridge S, Hamy V, Punwani S, Atkinson D: Direct parametric
reconstruction from undersampled (k, t)space data in dynamic contrast enhanced MRI.
Med Image Anal 2014; 18:989–1001.
105. Lebel RM, Guo Y, Zhu Y, et al.: The comprehensive contrastenhanced neuro exam.
In Proc Int Soc Magn Reson Med; 2015:3705.
91
106. Liu DC, Nocedal J: On the limited memory BFGS method for large scale optimization.
Math Program 1989; 45:503–528.
107. Shellock FG, Kanal E: Guidelines and recommendations for MR imaging safety and
patient management. III. Questionnaire for screening patients before MR procedures. J
Magn Reson Imaging 1994; 4:749–751.
108. Kanal E, Barkovich AJ, Bell C, et al.: ACR guidance document for safe MR practices:
2007. Am J Roentgenol 2007; 188:1447–1474.
109. Stanisz GJ, Henkelman RM: GdDTPA relaxivity depends on macromolecular
content. Magn Reson Med 2000; 44:665–667.
110. De Naeyer D, De Deene Y, Ceelen WP, Segers P, Verdonck P: Precision analysis of
kinetic modelling estimates in dynamic contrast enhanced MRI. MAGMA 2011; 24:51–
66.
111. Dikaios N, Punwani S, Atkinson D: Direct parametric reconstruction from (k, t)space
data in dynamic contrast enhanced MRI. In Proc Int Soc Magn Reson Med; 2015:3706.
112. Dickie BR, Banerji A, Kershaw LE, et al.: Improved accuracy and precision of tracer
kinetic parameters by joint fitting to variable flip angle and dynamic contrast enhanced
MRI data. Magn Reson Med 2015; 76:1270–1281.
113. Lingala S, Hu Y, DiBella E, Jacob M: Accelerated Dynamic MRI Exploiting Sparsity
and LowRank Structure: kt SLR. IEEE Trans Med Imaging 2011; 30:1042–1054.
114. Hu Y, Lingala SG, Jacob M: A fast majorizeminimize algorithm for the recovery of
sparse and lowrank matrices. IEEE Trans Image Process 2012; 21:742–53.
115. Miao X, Lingala SG, Guo Y, Jao T, Nayak KS: Accelerated Cardiac Cine Using
Locally Low Rank and Total Variation Constraints. In Proc Int Soc Magn Reson Med;
2015:571.
116. O’Connor JPB, Jackson A, Parker GJM, Roberts C, Jayson GC: Dynamic contrast
enhanced MRI in clinical trials of antivascular therapies. Nat Rev Clin Oncol 2012;
9:167–177.
117. Velikina J V, Alexander AL, Samsonov A: Accelerating MR parameter mapping
using sparsitypromoting regularization in parametric dimension. Magn Reson Med
2013; 70:1263–73.
118. Guo Y, Lingala SG, Zhu Y, Lebel RM, Nayak KS: Direct Estimation of Tracer
Kinetic Parameter Maps From Highly Undersampled Brain Dynamic Contrast
Enhanced MRI. Magn Reson Med 2016.
119. Lingala SG, Guo Y, Zhu Y, Barnes S, Lebel RM, Nayak KS: Accelerated DCE MRI
using constrained reconstruction based on pharmacokinetic model dictionaries. In
ISMRM; 2015:196.
120. Port RE, Knopp M V., Brix G: Dynamic contrastenhanced MRI using GdDTPA:
Interindividual variability of the arterial input function and consequences for the
assessment of kinetics in tumors. Magn Reson Med 2001; 45:1030–1038.
121. Samsonov A: A Novel Reconstruction Approach Using Model Consistency Condition
for Accelerated Quantitative MRI (MOCCA). In ISMRM; 2012:358.
122. Velikina J V., Samsonov AA: Reconstruction of dynamic image series from
undersampled MRI data using datadriven model consistency condition (MOCCO).
Magn Reson Med 2015; 74:1279–1290.
123. Barnes SR, Ng TSC, SantaMaria N, Montagne A, Zlokovic B V., Jacobs RE:
ROCKETSHIP: a flexible and modular software tool for the planning, processing and
92
analysis of dynamic MRI studies. BMC Med Imaging 2015; 15:19.
124. Barboriak, D. P., J. R. MacFall, A. O. Padua, G. E. York, B. L. Viglianti and MWD":
Standardized software for calculation of Ktrans and vp from dynamic T1weighted MR
images. In Int Soc Magn Reson Med Work MR Drug Dev From Discov to Clin Ther
Trials; 2004.
125. Bosca RJ, Jackson EF: Creating an anthropomorphic digital MR phantom—an
extensible tool for comparing and evaluating quantitative imaging algorithms. Phys Med
Biol 2016; 61:974–982.
126. Jung SC, Yeom JA, Kim J, et al.: Glioma : Application of Histogram Analysis of
Pharmacokinetic Parameters from T1Weighted Dynamic ContrastEnhanced MR
Imaging to Tumor Grading. 2014:1103–1110.
127. Guo Y, Lebel RM, Zhu Y, et al.: HighResolution WholeBrain DCEMRI Using
Constrained Reconstruction: Prospective Clinical Evaluation in Brain Tumor Patients.
Med Phys 2016.
128. Bertsekas DP: Multiplier methods: A survey. Automatica 1976; 12:133–145.
129. Murase K: Efficient Method for Calculating Kinetic Parameters Using T 1Weighted
Dynamic ContrastEnhanced Magnetic Resonance Imaging. Magn Reson Med 2004;
51:858–862.
130. Flouri D, Lesnic D, Sourbron SP: Fitting the twocompartment model in DCEMRI
by linear inversion. Magn Reson Med 2016; 76:998–1006.
Abstract (if available)
Linked assets
University of Southern California Dissertations and Theses
Conceptually similar
PDF
Fast flexible dynamic threedimensional magnetic resonance imaging
PDF
Shiftinvariant autoregressive reconstruction for MRI
PDF
Correction, coregistration and connectivity analysis of multicontrast brain MRI
PDF
Functional realtime MRI of the upper airway
PDF
Model based PET image reconstruction and kinetic parameter estimation
PDF
New theory and methods for accelerated MRI reconstruction
PDF
Methods for improving reliability and consistency in diffusion MRI analysis
PDF
Highdimensional magnetic resonance imaging of microstructure
PDF
Seeing sleep: realtime MRI methods for the evaluation of sleep apnea
PDF
Optimization methods and algorithms for constrained magnetic resonance imaging
PDF
Fast upper airway MRI of speech
PDF
New methods for carotid MRI
PDF
Improving sensitivity and spatial coverage of myocardial arterial spin labeling
PDF
Human appearance analysis and synthesis using deep learning
PDF
Radiofrequency pulse performance for myocardial ASL
PDF
Measuring functional connectivity of the brain
PDF
Susceptibilityweighted MRI for the evaluation of brain oxygenation and brain iron in sickle cell disease
PDF
Novel theoretical characterization and optimization of experimental efficiency for diffusion MRI
PDF
Graphbased models and transforms for signal/data processing with applications to video coding
PDF
Multiscale dynamic capture for high quality digital humans
Asset Metadata
Creator
Guo, Yi
(author)
Core Title
Improved brain dynamic contrast enhanced MRI using modelbased reconstruction
School
Viterbi School of Engineering
Degree
Doctor of Philosophy
Degree Program
Electrical Engineering
Publication Date
07/17/2017
Defense Date
03/08/2017
Publisher
University of Southern California
(original),
University of Southern California. Libraries
(digital)
Tag
compressed sensing,DCEMRI,image processing and reconstruction,modelbased reconstruction,OAIPMH Harvest
Language
English
Contributor
Electronically uploaded by the author
(provenance)
Advisor
Nayak, Krishna S (
committee chair
), Haldar, Justin P (
committee member
), Li, Hao (
committee member
)
Creator Email
eagle13gy@gmail.com,yiguo@usc.edu
Permanent Link (DOI)
https://doi.org/10.25549/uscthesesc40401213
Unique identifier
UC11264385
Identifier
etdGuoYi5529.pdf (filename),uscthesesc40401213 (legacy record id)
Legacy Identifier
etdGuoYi5529.pdf
Dmrecord
401213
Document Type
Dissertation
Rights
Guo, Yi
Type
texts
Source
University of Southern California
(contributing entity),
University of Southern California Dissertations and Theses
(collection)
Access Conditions
The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...
Repository Name
University of Southern California Digital Library
Repository Location
USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 900892810, USA
Tags
compressed sensing
DCEMRI
image processing and reconstruction
modelbased reconstruction