University of Southern California Dissertations and Theses

Technology for improved 3D dynamic MRI

(USC Thesis Other)

Technology for improved 3D dynamic MRI

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content Technology for Improved 3D Dynamic MRI
by
Ziwei Zhao
A Dissertation Presented to the
FACULTY OF THE GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFORNIA
In Partial Fulfillment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
(ELECTRICAL AND COMPUTER ENGINEERING)
December 2024
Copyright 2024 Ziwei Zhao

In memory of my grandfather.
ii

Acknowledgments
”Reason is, and ought only to be, the slave of passions.” Looking back to my journey, I have
mostly emerged myself into the MRI knowledge and experiments etc. Now it is my best time
to express my appreciation to the people who have surrounded me and accompany with me
during these difficult but fruitful time. I am truly grateful for the great people for their
tremendous help and support, to together create the passionate, positive atmosphere that
fulfilled my PhD.
I would like to express my deepest gratitude to my PhD advisor, Krishna S. Nayak. Your
generous support, guidance and care are precious gifts that I have. I am grateful that being
given the freedom to choose the topics that I am interested and want to explore, and yet
being guided at certain important stages. Your positive energy, and the way you think,
share and deliver has shaped my work and life philosophy. You have offered the full trust
throughout my entire PhD especially during the very difficult periods. The more I grow, the
more that I constantly feel being blessed under your supervision during my PhD, and later
for the rest of my life.
I am also very grateful to Professors Justin Haldar, Richard Leahy, C.-C. Jay Kuo, and
John Wood for being my qualifying and defense committee members. The insights that
you have provided have been helpful for doing several interdisciplinary studies. Thank you
for your time, effort, and valuable feedback, and for your influences on my training and
research. I am also grateful for Professors Shri Narayanan and Dani Byrd for providing
insightful advice for the speech MRI collaboration, from perspectives of speech scientists.
Thank you Roberta Kato for providing insights for pulmonary physiology and help with
the experiments. I also want to thank Professor Peder Larson, for offering great advice on
different projects of lung MRI.
I’m extremely grateful to all lab members and alumni at MREL. Yongwan, my first
iii

mentor, thank you for guiding me my first PhD project, which has now become one whole
chapter in the thesis. Nam, thank you for being my late night scan boddy. You have always
encouraged me to pursue with the important research topics and has provided tremendous
help from a broad range of the topics. To Ye, thank you for always giving insightful feedbacks
and always be helpful in the lab. To Zhibo, thank you for accompanying me at the airport
waiting for the early morning flight, and together undertook the TA duties for Digital Signal
Processing class. To Kubra, Ecrin, and Sarina, thank you for setting up the ’Girl’s coffee
time’ within MREL, with smiling faces and support. To Bochao, thank you for lots of
interesting chats in the office, making EEB 416 lots of fun. I also thank Prakash for his
generous help to take care of the server and help me set up the processing scripts during the
rush hours before ISMRM deadline. Thank you, Bilal, for collaborating on the Pilot Tone
project. Your simple and elegant code has inspired me. To Sophia, thank you for bringing
lots of support from Siemens. To Mary, Diane, and Gloria, thank you for keeping track of
the logistics, and make my life easier. I also thank the past MREL members, Drs Brian
Ziyue Wu, Kyung Song, Wesley Zun, Weiyi Chen, Xin Miao, Yannick Bliesener for giving
advice on my PhD topics selection career and life advices. MREL has always been my second
home.
I would like to also thank my friends outside of lab which has fulfilled my life during
PhD. To Philip K. Lee, thank you for teaching me the diffusion-prep pulses and the logic
to debug sequences, it was enjoyable to work and chat about small little things in life with
you. To Xingfeng Shao and Yovnne Xiao, and Grant Yang, thank you for being friends
and gathering together during festivals. To Zalan Fabian, Hesameddin Mohammadi, Sen
Ma, Jiayang Wang and Chenyang Zhao, thank you for the great conversations we have, and
provide helps for varieties of the projects. I would also thank my friends Fei Tan, Yuke
Zhang, Ge Cui, Yijun Liu, Yu Cao, Xiaoyu Liu, Erum Mastaq, and Chenwei Tang for all
the friendship and happiness we have created.
My completion of my PhD would not have been possible without the unwavering support
iv

and encouragement of my grandma, dad, and especially my mom. From the countless phone
calls to the 18-hour flights just to visit, I am deeply grateful for all the efforts and sacrifices
you’ve made to give me more. Thank you from the bottom of my heart.
Lastly, I want to give my thank and love to my fianc´e Bowen. Your kindness, generous
sharing, your positive energy, the constant enthusiasm of life have always been inspiring me.
You have provided the closest care, happiness and love that nurture me to grow further.
Having you had been unexpected but the best things that happened to me.
v

Table of Contents
Dedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii
Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii
List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix
List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . x
Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiii
Chapter1:
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.1 3D Dynamic MRI . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2 Imaging of the vocal tract . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.2.1 Articulation of speech . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.2.2 Imaging modalities . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.3 Imaging of the lung . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3.1 Lung anatomy and functions . . . . . . . . . . . . . . . . . . . . . . . 6
1.3.2 Imaging modalities . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
Chapter2:
Magnetic Resonance Imaging . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
2.1 MRI fundamentals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
2.1.1 Nuclear magnetic resonance . . . . . . . . . . . . . . . . . . . . . . . 11
2.1.2 Excitation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.1.3 Acquisition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.1.4 Reconstruction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.2 Data sampling strategies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2.1 Non-Cartesian sampling . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.2.2 Reconstruction with undersampled Data . . . . . . . . . . . . . . . . 25
2.3 Field strength dependence . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.3.1 Equilibrium polarization and SNR . . . . . . . . . . . . . . . . . . . . 27
2.3.2 Susceptibility effects, T
∗
2
and intravoxel dephasing . . . . . . . . . . . 28
2.3.3 Concomitant field effects . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.4 MRI of Speech Production . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
2.4.1 Imaging requirements . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
vi

2.4.2 Imaging considerations . . . . . . . . . . . . . . . . . . . . . . . . . . 32
2.4.3 Unmet needs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
2.5 MRI of Lung Ventilation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
2.5.1 Imaging considerations . . . . . . . . . . . . . . . . . . . . . . . . . . 36
2.5.2 Unmet needs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
Chapter3:
Improved 3D Real-Time MRI of Speech Production . . . . . . . . . . . . . 40
3.1 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
3.1.1 Experimental Methods . . . . . . . . . . . . . . . . . . . . . . . . . . 40
3.1.2 Data Sampling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.1.3 Image Reconstruction . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.1.4 Regularization Parameter Selection . . . . . . . . . . . . . . . . . . . 44
3.1.5 Experimental Optimization . . . . . . . . . . . . . . . . . . . . . . . 44
3.1.6 Speech Production Experiments . . . . . . . . . . . . . . . . . . . . . 45
3.1.7 Speech Stimuli . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
3.1.8 Evaluation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
3.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.2.1 Regularization Parameter Selection . . . . . . . . . . . . . . . . . . . 49
3.2.2 Experimental Optimization . . . . . . . . . . . . . . . . . . . . . . . 50
3.2.3 Application to Speech Production . . . . . . . . . . . . . . . . . . . . 53
3.3 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
3.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
Chapter4:
Free-breathing 3D pulmonary ventilation mapping at 0.55T using Stackof-Spiral Out-in bSSFP . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
4.2 Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
4.2.1 Experimental methods . . . . . . . . . . . . . . . . . . . . . . . . . . 65
4.2.2 Pulse sequences . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
4.2.3 Respiratory navigator and data binning . . . . . . . . . . . . . . . . . 66
4.2.4 Image reconstruction . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
4.2.5 Ventilation analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
4.2.6 In-vivo experimental study design . . . . . . . . . . . . . . . . . . . . 70
4.2.7 Image quality analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 72
4.2.8 Statistical analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
4.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
4.3.1 Sequence and reconstruction parameter optimizations . . . . . . . . . 72
4.3.2 Structure and ventilation maps . . . . . . . . . . . . . . . . . . . . . 77
4.3.3 Repeatability and validation . . . . . . . . . . . . . . . . . . . . . . . 78
4.3.4 Gravity effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
4.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
4.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83
vii

Chapter5:
Multidimensional RF Pulse Design with Consideration of Concomitant
Field Effects . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
5.2 Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86
5.2.1 Concomitant fields . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86
5.2.2 Approximation of concomitant fields as a Bloch-Siegert shift . . . . . 89
5.2.3 Solving Bloch equations with concomitant fields in a rotating frame . 91
5.3 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93
5.3.1 Validation of Bloch-Siegert Approximation . . . . . . . . . . . . . . . 93
5.3.2 Evaluation of Bloch-Siegert Approximation - Multi-channel RF design 93
5.3.3 Evaluation using single-channel RF design and phantom experiment
at 0.55T . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94
5.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96
5.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103
5.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104
Chapter6:
Concluding remarks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
6.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
6.2 Future directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111
viii

List of Tables
1.1 Comparisons of vocal tract imaging modelities . . . . . . . . . . . . . . . . . 6
1.2 Comparisons of lung imaging modalities . . . . . . . . . . . . . . . . . . . . 10
3.1 Stimuli design I . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
3.2 Stimuli design II . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.1 In-vivo experimental design . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
ix

List of Figures
1.1 Examples of dynamic MRI applications . . . . . . . . . . . . . . . . . . . . . 2
1.2 Vocal tract imaging modalities . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.3 Lung anatomy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
1.4 Lung imaging modalites . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
2.1 Polarization of magnetic spins . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2.2 Forced precession of magnetization . . . . . . . . . . . . . . . . . . . . . . . 14
2.3 Relaxation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
2.4 Components of MR scanner . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.5 Selective excitation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.6 Relationship between k-space and object space . . . . . . . . . . . . . . . . . 21
2.7 Non-Cartesian trajectories . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
2.8 2D Spiral bSSFP diagram . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.9 3D Spiral bSSFP diagram . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.10 Spatiotemporal resolution requirements for various speech task . . . . . . . . 33
2.11 MR Speech images . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.12 MR Lung ventilation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
2.13 MRI of lung parenchyma at 1.5T and 0.55T . . . . . . . . . . . . . . . . . . 38
x

3.1 Data sampling patterns . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42
3.2 Results of parameter sweep . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.3 Pairwise Comparisons . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
3.4 Retrospective selection of temporal resolution . . . . . . . . . . . . . . . . . 53
3.5 Demonstration of speech application of the proposed method in stimuli set 1 55
3.6 Demonstration of speech application of the proposed method (Case V) in
stimuli set 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
3.7 Illustration of sublingual cavity during a rhotic flap and the tongue side channeling during a laterally release tap . . . . . . . . . . . . . . . . . . . . . . . 59
4.1 Free-breathing SOS out-in pipeline . . . . . . . . . . . . . . . . . . . . . . . 67
4.2 Flip angle optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
4.3 Structure and ventilation results for different choices of regularization parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
4.4 Impacts of motion-states on structure and ventilation . . . . . . . . . . . . . 76
4.5 Example structural images and ventilation maps from one healthy volunteer 77
4.6 Repeatability studies of ventilation . . . . . . . . . . . . . . . . . . . . . . . 79
4.7 Ventilation comparison with PREFUL . . . . . . . . . . . . . . . . . . . . . 80
4.8 Posture-related ventilation differences . . . . . . . . . . . . . . . . . . . . . . 81
5.1 Derivation of concomitant gradient terms from concomitant fields using the
Bloch-Siegert shift . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
5.2 Evaluation of the accuracy of the approximation of concomitant fields as concomitant gradient terms using the BS shift at various field strengths . . . . . 97
5.3 Off-isocenter dependence of concomitant field effects at 0.55T without and
with consideration and correction of the concomitant fields . . . . . . . . . . 98
xi

5.4 Simulated B0 dependence of concomitant fields obtained without and with
modeling of concomitant fields effects . . . . . . . . . . . . . . . . . . . . . . 100
5.5 Single-channel simulation and phantom validation results at 0.55T . . . . . . 102
5.6 Representative optimized waveforms obtained with the proposed method at
0.55T with parallel transmission . . . . . . . . . . . . . . . . . . . . . . . . . 105
5.7 Pulse duration dependence of T2 effects in conjunction with concomitant field
effects at 0.55T . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106
xii

Abstract
Magnetic resonance imaging (MRI), without ionizing radiation, is a non-invasive medical
imaging modality that can visualize structural anatomy as well as functional responses.
three-dimensional (3D) dynamic MRI provides an efficient solution to visualize dynamic
movements and to reveal functional responses volumetrically. This dissertation is focused
on 3D dynamic MRI methods on vocal tract during real-time, and 3D pulmonary ventilation mapping during free-breathing. It targets on improving three limitations faced by
these applications: limited spatial-temporal resolution, low parenchyma signals, and system
imperfections.
First, I describe a 3D real-time MRI of speech production at 1.5T with improved spatiotemporal sharpness using randomized, variable-density, stack-of-spiral (SOS) sampling combined with a 3D spatiotemporally constrained reconstruction. Real-time MRI has emerged
as one of the most powerful tools for studying human speech production. In this work,
I evaluated five candidate (k, t) sampling strategies using a previously proposed gradientecho stack-of-spiral sequence and a 3D constrained reconstruction with spatial and temporal
penalties. The strategy yielding highest image quality was chosen as the proposed method.
The proposed method provides the best qualities with the improved spatiotemporal resolution for better visualization of fast lip and tongue movements. The significant improvements
are demonstrated in tongue boundary sharpness (p<0.001) during normal and 1.5× speeded
speech. This allows researchers to demonstrate both functional and morphological aspects
of the vocal tract during speech and to understand the complex spatiotemporal coordination
of upper airway structures in motion.
Next, I describe an efficient free-breathing 3D functional lung imaging using SOS Outin balanced steady-state free-precession (bSSFP). MRI, as an alternative for computed
tomography, also provides functional response in screening and diagnostic assessment for
xiii

many lung diseases, such as cancer, emphysema, and chronic obstructive pulmonary disease, without ionizing radiation. The proposed method, for the first time, demonstrates 3D
regional pulmonary ventilation mapping at 0.55T within 5min scan using 3D bSSFP acquisition, constrained reconstruction and non-rigid registration. It overcome the challenges of
low parenchyma signals of lung MRI in mid- and high-field strengths (≥ 1.5T). In healthy
volunteers, the SOS out-in lung images provided sufficient vessel-parenchyma contrast and
boundary sharpness to support accurate ventilation estimation. Regional ventilation measurements from 3D SOS out-in demonstrated good repeatability (relative differences¡10%).
Ventilation maps from 3D SOS out-in strongly correlated with Phase-resolved Functional
Lung (PREFUL) MRI on a slice-matched basis, as well as with global tidal volume (R2 >0.7,
p<0.001). Therefore, ventilation measurements are sensitive, consistent, and in good agreement with PREFUL and spirometry.
Lastly, I describe my efforts in developing a small-tip multidimensional radiofrequency
(RF) pulse design procedure that incorporates linear time-invariant gradient imperfections
and concomitant field effects. I developed an extension of the small-tip excitation k-space
formalism, where concomitant fields were approximated as a Bloch-Siegert shift in the rotating frame. This was evaluated using realistic simulations of 2D selective excitation at
various field strengths (0.2T, 0.55T, 1.5T,3T, and 7T) with single and parallel transmit.
The extended formalism provides improved 2D excitation profiles in all scenarios simulated,
compared against the original formalism. The proposed method corrects the concomitant
field effects on 2D selective excitations for B0 ≥ 0.2T when the magnitude of the B0 is far
larger than that of nonrotating concomitant fields. The proposed method providing sharper
and more accurate excitation profiles at off-isocenter distances up to 15 cm. The impact of
the proposed method is greatest in scenarios where concomitant fields are substantial, such
as low field strengths and off-isocenter. This could be particularly important for contemporary low-field MRI systems with high-performance gradients, with potential applications in
off-isocenter navigators during real-time volumetric imaging.
xiv

A classic description of the principles of MRI, as well as the summary and future remarks
about this dissertation are also presented.
xv

Chapter 1
Introduction
Over the years, several imaging modalities have been developed. These modalities are used
in clinical diagnosis and treatment. They have various applications. Some examples are
computed tomography (CT), MRI, ultrasounds, optical coherence tomography (OCT),
single-photon emission computed tomography (SPECT), and positron emission tomography
(PET). Among these, MRI is unique because it does not use ionizing radiation. This sets it
apart from CT and X-ray. MRI provides superior soft tissue contrast. This makes it better
for imaging the brain and spinal cord. MRI also offers higher resolution than ultrasound.
This allows for better imaging of organs like the liver, kidneys, or heart. Additionally, MRI
can penetrate deeper than OCT. This provides detailed images of deeper tissues and organs.
MRI also produces detailed three-dimensional images.
1.1 3D Dynamic MRI
Dynamic MR imaging is a process that involves capturing and visualizing changing images
over time, commonly used in medical fields to observe real-time physiological processes.
This technique improves diagnostic accuracy by providing motion-based data on organs or
systems, crucial for understanding the timing and function of biological activities. Several
clinical applications are benefiting from dynamic imaging techniques, including understanding of human vocal tract shaping, cardiac structure and function, and dynamic movement
1

Figure 1.1: Examples of dynamic MRI applications. a) human vocal tract shaping [Adapted
from Ref. [1]]; b) cardiac structure and function [Adapted from Ref. [1]]; c) dynamic
movement of wrist [Adapted from Ref. [2]]; d) gastric motility [Adapted from Ref. [2]]; e)
characterization of fetal brain, heart development and growth [Adapted from Ref. [1]] and
f) capturing of dynamic lung imaging during breathing [Adapted from Ref. [1]]. (Figure
courtesy of Dr. Ye Tian, Los Angeles, CA, USA.)
of the wrist; characterizing fetal brain and heart development and growth; and capturing
dynamic lung imaging during breathing; and gastroenterology, which is used to evaluate
swallowing disorders and monitor gastrointestinal motility. Figure 1.1 shows examples of
two-dimensional (2D) dynamic MRI applications.
Dynamic imaging provides additional information about changes over time. Volumetric imaging offers a complete and realistic representation of medical scans. Traditional 2D
imaging visualizes organs in two dimensions. It is limited in height and width on a single
plane. In contrast, 3D imaging provides depth. Together, 2D and 3D imaging offer a com2

prehensive visualization in an x-y-z coordinate volume. 2D imaging typically has shorter
scan times. However, it offers limited visualization and provides less information. 3D volumetric imaging, on the other hand, gives clinicians, scientists, and healthcare professionals a
comprehensive perspective. It helps them understand the complex geometry of organs and
their related functional responses.
In this thesis, I focus on demonstrating the clinical value and technical requirements with
a focus on lung and speech, which are two main applications of 3D dynamic MRI.
1.2 Imaging of the vocal tract
1.2.1 Articulation of speech
Speech production is human’s unique process by which we generate and articulate speech,
turning our thoughts into spoken words. It involves a complex sequence of cognitive, neurological, and physical actions, allowing us to communicate verbally.
The articulation of the vocal tract is an important step in speech production. This
includes the coordination of the muscles required for producing speech sounds. For example,
the tongue, lips, jaw, vocal cords, and diaphragm are all involved in this motor planning.
During articulation, the vocal tract (including the lips, tongue, hard palate, soft palate,
larynx, and epiglottis etc.) creates the sounds that form words. Air from the lungs is pushed
through the vocal cords, which vibrate to create sound. The sound is then modified by
the movements of the articulators (e.g., tongue, lips, soft palate) to produce specific speech
sounds. The sound is then created when vocal cords (larynx) vibrate as air passes through
them.
The vocal tract consists of several soft tissues and structures, each with specific properties
that affect how sound is produced and modulated. When imaging the vocal tract for speech
production studies, several key parameters and features are measured to capture the tissue
properties of vocal tract, to understand the anatomy, physiology, and dynamics of speech.
3

1) measuring vocal tract anatomy and shape. 2) capturing tougue, lips shapes, positions
and movements. For example, researchers want to track how the tongue moves for different
consonants (e.g., /t/, /k/) and vowels (e.g., /i/, /a/). Lip shapes affect the articulation of
sounds such as /p/, /b/, and /m/. 3) measuring velum (soft palate) movement, which is
essential for controlling airflow between the nasal and oral cavities, and seeing pharyngeal
movements and laryngeal functions. 4) measuring the airflow and pressure. 5) Acoustic
analysis that allows for an understanding of the mapping between articulation (movement)
and acoustics (sound). Understanding different speech properties help researchers, linguistics
have a better understanding of human natural speech.
1.2.2 Imaging modalities
Tissue contrasts of vocal tract is important in imaging. Soft tissues such as the tongue, vocal
folds, and soft palate differ significantly from the bony structures of the skull and mandible.
These contrasts are crucial for visualizing the vocal tract’s soft tissue dynamics (e.g., tongue
movements). MRI and CT can well capture these specific contrasts. Muscle (e.g., tongue)
and fat (e.g., around the cheeks) can also be captured in MRI. In ultrasound, the contrast
is based on the differing densities and stiffness of tissues (e.g., soft tissues like the tongue
vs. more solid structures like the mandible). These differences allow ultrasound to visualize
tongue movements during speech.
There are several imaging modalities that have been used in speech imaging. Figure 1.2
provides a brief summarization. that includes electromagnetic articulography (EMA) [8],
endoscopy [9], computed tomography (CT) [10], ultrasound [11, 12], and MRI [13, 14, 15].
Table 1.1 lists and compares these modalities. EMA tracks the exact position of specific
articulator (e.g., the tongue tip) using attached sensor coils [8]. Endoscopy involves using
a small camera inserted into the throat (rigid or flexible) to visualize internal structures
like the vocal folds and larynx [9]. However, these modalities are invasive and provide
information on soft tissue surfaces. CT creates detailed cross-sectional images of the vocal
4

Figure 1.2: Vocal tract imaging modalities. a) EMA, electromagnetic articulography
[Adapted from Ref. [3]]; b) endoscopy [Adapted from Ref. [4]]; c) ultrasound [Adapted
from Ref. [5]]; d) CT, computed tomography [Adapted from Ref. [6]]; e) MRI, magnetic
resonance imaging: sagittal and coronal views [Adapted from Ref. [7]]. MRI is the most
promising imaging modality for imaging vocal tract airway anatomy and function, as it
provides detailed tomography, adequate spatiotemporal resolution, excellent soft-tissue contrast, and involves no ionizing radiation or invasive procedure.
5

Table 1.1: Comparisons of vocal tract imaging modelities. Modailty acronyms: EMA =
electromagnetic articulography; CT = computed tomography; MRI = magnetic resonance
imaging.
tract [10], particularly useful for visualizing both soft tissues and bony structures. However,
CT involves ionizing radiations, and provides low soft tissue contrasts compared to MRI.
Ultrasound is ideal for real-time visualization of tongue movements but is limited in depth
and soft-tissue resolution [11, 12]. Among these modalities, MRI is preferred for its highresolution, non-invasive imaging of soft tissue and vocal tract dynamics, without ionizing
radiation, and provides the excellent soft tissue contrasts.
1.3 Imaging of the lung
1.3.1 Lung anatomy and functions
Breathing is a vital and essential function for our daily life. The lungs are the primary organs
of the respiratory system and are crucial for breathing, as they facilitate the gas exchanges
between the blood and the air.
The lungs are located in the thoracic cavity, one on each side of the heart, and are
6

Figure 1.3: Lung anatomy. a) Respiratory system [Adapted from Ref. [16]] and b) gas
exchange at alveolar walls [Adapted from Ref. [17]]. The O2 travels from the alveoli into the
capillaries, while CO2 travels from the blood into the alveoli. The waste CO2 in the alveoli
is then expired, which means coming out of the mouth.
protected by the rib cage.They are divided into sections called lobes, and their primary role
is to facilitate the exchange of oxygen and carbon dioxide during respiration. Figure 1.3
shows a detailed lung structure. The lungs are surrounded by a double-layered membrane
called the pleura. The airways in the lungs are organized in a branching pattern. The trachea
divides into two main stem bronchi for the right and left lungs. Each main stem bronchus
then divides into smaller and smaller branches like a tree. The tiniest branch is called a
bronchiole. At the end of each bronchiole is a cluster of tiny air sacs called alveoli. These
alveolus has a lot of holes that allow air to pass through. Looking closer to the alveolus,
it has a mesh of tiny blood capillaries. The very thin walls of the alveoli and bronchioles
provide an extremely thin and large surface for gas exchange to occur. The exchange of
oxygen and carbon dioxide occurs at the alveolar-capillary membrane. Oxygen from the
inhaled air diffuses across the alveolar walls into the blood in the capillaries, while carbon
dioxide diffuses in the opposite direction to be exhaled.
The primary functions of the lungs are ventilation and gas exchanges. The term to
7

describe the gas exchange is called ventilation. It is a process of moving air into and out
of the lungs. During inspiration, also known as breathing in, the air flows into the lung,
causing the total volume of the lung to increase. During exhalation, the air flows out of
the lung. The recoil force brings the lungs back to the resting position. The exchange of
oxygen and carbon dioxide occurs at the alveolar-capillary membrane. Oxygen from the
inhaled air diffuses across the alveolar walls into the blood in the capillaries, while carbon
dioxide diffuses in the opposite direction to be exhaled.Though this breathing process, the
lung volume plays a major role in the gas exchange. This changing lung volume allows us
to identify the normal and abnormal lung conditions. Pulmonary circulation also happens
in the lungs as the right side of the heart pumps deoxygenated blood into the lungs for
oxygenation.
1.3.2 Imaging modalities
There are several imaging modalities that help pulmonary dysfunction diagnosis in the clinic.
Figure 1.4 displays the image qualities for each modality, and Table 1.2 compares these
modalities.
Pulmonary function test (PFT) is a standard way to detect obstructive and restrictive
lung diseases [22]. It is non-invasive and can measure the lung volume flow rate as well as the
amount of gas exchange, showing in a 1D curve reflecting the volume change. CT imaging,
which provides lung parenchyma density information, can be used for disease progression
tracking and therapy response monitoring [23, 24, 25]. Yet, the risks of ionizing radiation
remains a concern for patients, especially for young children and newborns [26]. Moreover,
the rapid growth and development of the lung in neonates requires frequent imaging, and
lung is highly radiation-sensitive, which would cause a significant radiation accumulation
if relying on CTs [27] to monitor lung disease. Ventilation-perfusion (V/Q) single-photon
emission computed tomography (SPECT) is a nuclear medicine test that evaluates the
ventilation (airflow) and perfusion (blood flow) in the lungs by using radioactive tracers [28].
8

Figure 1.4: Lung imaging modalites. a) PFT, pulmonary function test [Adapted from Ref.
[18]]; b) V/Q SPECT, Ventilation-perfusion single-photon emission computed tomography
[Adapted from Ref. [19]]; c) CT, computed tomography [Adapted from Ref. [20]]; d) MRI,
magnetic resonance imaging: coronal and sagittal views [Adapted from Ref. [20]]; e) More
possibilities of functional lung MRI [Adapted from Ref. [21]]. MRI is one emerging imaging
modality for imaging lung structure and function, as it provides detailed tomography, tissue
characterization, captures functional responses, and involves no ionizing radiation.
9

Table 1.2: Comparisons of lung imaging modalities. Modality acronyms: PFT = pulmonary
function test; SPECT = single-photon emission computed tomography; CT = computed
tomography; MRI = magnetic resonance imaging.
Among these modalities, MRI, provides superior soft-tissue contrasts, tissue characterization
and functional response in screening, diagnostic assessment for many lung diseases, such as
cancer, emphysema, and chronic obstructive pulmonary disease, without ionizing radiation.
10

Chapter 2
Magnetic Resonance Imaging
2.1 MRI fundamentals
In this chapter, a basic overview of MRI principles is introduced. The phenomenon of
nuclear magnetic resonance is briefly described. The fundamental procedure to generate an
MR image is illustrated. Several acceleration techniques are also introduced. The physics
of MRI can be generally understood and analyzed in a classical way. It is also described in
many popular textbooks [29, 30, 31, 32] with greater detail.
2.1.1 Nuclear magnetic resonance
The phenomenon of nuclear magnetic resonance (NMR) was first discovery of I.I. Rabi
[33] in 1937. Then, it was experimentally demonstrated by Felix Bloch [34] and Edward
Purcell [35] in 1946, independently. Since then, Richard Ernst developed magnetic resonance
spectroscopy (MRS) that uses NMR phenomenon to identify chemical compositions based
on different atomic resonance frequencies in 1966 [36]. Later, Kurt W¨uthrich determined the
three-dimensional structure of biological macromolecules in 1982 using MRS [37]. In 1973,
Paul Lauterbur generated the first MR image that uses gradient magnetic fields to spatially
encode the objects [38]. Then, Peter Mansfield introduced the mathematical formalism and
also developed echo-planar imaging (EPI), which opens a door for functional MRI [39]. So
11

Figure 2.1: Polarization of magnetic spins. a) The directions of spins are random when no
external magnetic field exists. b) When an external magnetic field is applied, all the spins
align to the applied field, This results in a net magnetization that points towards B0.
far, all aforementioned researchers and scientists have won several Nobel Prizes in physics,
chemistry and physiology or medicine.
Atoms with an odd number of protons and/or odd number of neutrons process a nuclear
spin angular momentum. The spin angular momentum can be visualized as a spinning
charged spheres that leads to a small magnetic moment. These MR-relevant nuclei will
often be called ‘spins’. The angular momentum arises the consequences of both spin and
mass. Hydrogen imaging is the most studied and imaged in MR imaging, since it is mostly
dominant because of the large portion of H2O in the body. Therefore, we assume 1H imaging
in this thesis.
Polarization. When interacting with static magnetic field B0, a magnetic moment M
is produced in the direction of B0, for example, the z direction. This bulk magnetization, M
is an ensemble behavior the sum of billions of magnetic dipoles: M =
Pµ. It is worthy to
notice that without external magnetic field is present, the directions of these charged spins
are random, therefore, the net magnetization, M, is zero. (See Figure 2.1) With an external
static magnetic field is applied, all the spins align either parallel (n−) or anti-parallel (n+)
12

to the applied field with a slightly greater number in the parallel direction. Their ratio is
given by Boltzmann distribution:
n−
n+
= e
−
γℏB0
kT (2.1)
with γ to be gyromagnetic ratio. For 1H, γ/2π = 42.58 MHz/T. ℏ represents Planck’s
constant, k is Boltzmann’s constant, and T is the absolute temperature. The equilibrium
nuclear magnetization M0 can be calculated as
M0 =
N γ2ℏ
2
Iz(Iz + 1)B0
3kT (2.2)
with N is the number of spins per unit volume, Iz = 1/2 for 1H. The direction of magnetization aligns with external magnetic field B0.
Precession. If the magnetization M0 is perturbed such that it doesn’t align with the
B0 direction, it experiences torque. The dynamics of an individual magnetic dipole under
externally applied B field can be described as
dµ
dt = µ × γB. (2.3)
If we consider the summation of µ over a unit volume, the equation can be derived as:
dM
dt = M × γB. (2.4)
The solution to this equation is known as the Larmor frequency, which is defined as
ω = γB. (2.5)
The resonant magnetization processes about B at certain Larmor frequency for different
nucleus, which is often noted as free procession.
13

Figure 2.2: Forced precession of magnetization. a) After certain RF pulse (B1) is applied, M
is perturbed from its equlibrium state. Meanwhile, it also precesses about the B0 axis at its
Larmor frequency ω. b) The same process is usually presented in the rotating frame whose
coordinates rotate about B0 at the spins’ Larmor frequency. As a result, B1 appears to be
static and M rotates about B1.
When certain external radiofrequency (RF) pulse is applied with a carrier frequency
that matches with the spin Larmor frequency, the spins will be excited at the applied angle
related to the RF pulse and perturbed from its equilibrium state. It will also rotate about
the axis that is align with its Larmor frequency. The applied RF pulse also rotates about
the B0 direction at the spin’s Larmor frequency, can it be often noted as B1. This forced
procession can be viewed in both Laboratory frame and rotating frame in Figure 2.2.
Relaxation. After the B1 field pulse is played, the transverse component of the magnetization decays away while the longitudinal component returns to its thermal equilibrium
state. These are considered as longitudinal and transverse relaxation, accordingly. As shown
in Figure 2.3, both of them follow exponential decay, and are written as:
Mz(t) = M0(1 − e
−t/T1
) (2.6)
Mxy(t) = M0e
−t/T2
(2.7)
Here we introduce two important parameters T1 and T2, which are used in describing the
relaxation time of tissues [40]. The two relaxation processes are also known for T1 recovery
14

Figure 2.3: Relaxation. a) Longitudinal relaxation. b) Transverse relaxation.
and T2 decay, respectively. T1 is dependent on field strength, while T2 is largely independent
of field strength. T2 also counts for the microscopic dephasing of the transverse component
and is usually shorter than T1. In MR literature, different tissue has different T1 and T2
values, engineers often use T1 and T2 values as references to generate the desired tissue
contrasts. In the clinic, T1 and T2 values are often treated as important indicators that help
the diagnosis.
Considering tissue relaxation in the Bloch Equation, the dynamics of the magnetization
can be described as:
dM
dt = M × γB −
Mxi + Myj
T2
−
(Mz − M0)k
T1
(2.8)
where i, j, and k are unit vectors.
With a static and homogeneous field that can be presented as B(t) = B0k, the solution
can be expanded as:
15

M(t) =


e
−t/T2 0 0
0 e
−t/T2 0
0 0 e
−t/T1


Rz(ω0t)M0 +


0
0
M0(1 − e
−t/T1 )


. (2.9)
Here, Rz(ω0t) represents the rotation matrix where M precessing about z-axis. M0 is an
initial magnetization vector. The receiver coils in MR system pick up the transverse magnetization, which often can be written as:
Mxy ≜ Mx + iMy (2.10)
with an initial condition of
M0 ≜ Mx,0 + iMy,0 (2.11)
Therefore, the transverse magnetization, i.e., detected signals can be written as
Mxy = M0e
−t/T2 e
−iω0t
e
−i∆ωt (2.12)
We introduce the last exponential term that considers the local field inhomogeneities which
can be potentially caused by nonuniform B0 and susceptibility.
The MR system has three main components as shown in Figure 2.4 : 1) main field
B0, 2) radiofrequency field B1, and 3) liner gradient fields G. The following chapters will
mainly describe how the spins are interact with these different fields to excite, encode and
reconstruct MR signals.
2.1.2 Excitation
With the presence of the main field B0, the B1 field is turned on and is designed to make the
magnetization rotate away from the z-axis to produce transverse components. This is often
referred to excitation. The general approach is applying an oscillating RF magnetic field B1
16

Figure 2.4: MR scanner and components. a) A whole body 0.55T system (prototype MAGNETOM Aera; Siemens Healthineers, Erlangen, Germany) equipped with high-performance
shielded gradients (45 mT/m amplitude and 200 T/m/s slew rate), which was used to conduct part of experiments in this dissertation. b) Key components of a typical MR scanner
[Adapted from Ref. [41]].
in the transverse direction, modulated by spins’ Larmor frequency. There are two types of
excitations. In the most basic case, RF is turned on with only the B0 exist, therefore excites
all spins in the volume. This is commonly referred to non-selective excitation. In a more
general case, with additional gradient field turned on plus RF fields, the spins that locate
only in certain slices of the volume get excited. This is often called selective excitation. The
designed RF pulses can have variants of shapes, depending on the excitation requirements.
In general, the phase of RF pulse decides which axis it is applied. For example, phase of 0
degree is about x-axis, 90 is about y-axis. The flip angle (θ) is defined as
θ =
Z τ
0
ω1(s)ds (2.13)
for general time-varying B1(t). Often the case, we denote 90◦
x pulse that specifies a rotation
of 90◦ about the x-axis in the rotating frame.
Selective excitation. To excite a selectable plane in the volume, the basic approach
17

Figure 2.5: Selective excitation. The gradient coils create a linear gradient field on top of B0
across the whole object. The gradient amplitude can be tuned so that the RF’s frequency
band matches the spins’ Larmor frequency range within the region of interest. Only the
region satisfies resonance condition will be excited.
is to apply an RF pulse with the presence of a static gradient, for example, Gz. With the
applied static gradient on z-axis, a spatial distribution of the frequencies is created. With a
band-limited RF pulse played at the same time, only the spins whose resonance frequencies
fall within the range of RF pulse’s bandwidth will be excited. As shown in Figure 2.5, the
excitation thickness can be adjusted by changing the gradient amplitude.
Small tip approximation. Solving the Bloch Equation with the presence of an applied on-resonance B1 field, a small tip approximation is usually assumed (Mz ≈ M0, and
dMz/dt ≈ 0) to simplify the solution to be:
Mr(t, z) = iM0e
iω(z)t
Z t
0
e
iω(z)sω1(s)ds (2.14)
where ω represents the Larmor frequency (on-resonance), ω1 = γB1 represents the excitation frequency of B1 field. This equation holds with the initial condition to be Mrot(0) =
[0, 0, M0]
T
. Intuitively, the slice profile can be therefore approximated by the Fourier trans18

form of the RF envelop. This provides the insights in designing RF pulses to achieve the
desired slice profile. For example, the Sinc shaped RF pulse is commonly used in the sliceselective excitation to create a near-rectangular slice profile.
2.1.3 Acquisition
In MR imaging system, the receiver coil is designed to detect the flux changes in the transverse direction, and it is ideally sentience over the entire volume-of-interest. Therefore, this
Free Induction Decay (FID) signal sr(t) is derived from the contributions of all transverse
magnetization, which can be written as:
sr(t) = Z Z Z M0(x, y, z)e
−t/T2(r)
e
−iω0t
dxdydz, (2.15)
which is the integral over all three dimensions of the volume. It cannot be resolved spatially
by FID. In order to form an image, each voxel needs to be spatially encoded. We consider a
2D function that we are interested in imaging which takes the integral over the slice width.
We define
m(x, y) ≡
Z zo+∆z/2
zo−∆z/2
M0(x, y, z)dz (2.16)
Therefore, we applying time-varying gradient fields G(t), which can create a linear phase
across the subject in an arbitrary direction:
sr(t) = Z Z m(x, y)e
−iτ R t
0 G(τ)rdτdxdy (2.17)
The term e
−t/T2(r)
is dropped if ignoring the relaxation. The term e
−iω0t
is dropped since a
phase sensitive detection is often dropped so that sr(t) is demodulated by ω0.
19

If considering an 1D static object, the signal becomes
s(t) = Z
x
m(x)e
−iγ R t
0 Gx(τ)xdτdx
=
Z
x
m(x)e
−i2πkx(t)x
dx
= s(kx)
(2.18)
with kx(t) = γ/2π
R t
0 Gx(τ )xdτ . This above equation follows a Fourier transform pair
with s(kx) in the frequency domain, and corresponding 1-D object m(x) in object domain.
By sampling s(kx) that often meets the Nyquist rate, m(x) can be calculated using the
inverse Fourier transform.
If considering a 2D object on x and y domains, the additional Gradient field that is played
on the y axis is introduced. Similarly, the signal can be written as
s(t) = Z
x
Z
y
m(x, y)e
−i2π[kx(t)x+ky(t)y]
dxdy
= s(kx, ky)
(2.19)
with ky(t) = γ/2π
R t
0 Gy(τ )ydτ , and kx(t) = γ/2π
R t
0 Gx(τ )xdτ . Here we introduce kx(t)
and ky(t) that represent the time integrals of the gradient waveforms along two dimensions.
The signal s(kx, ky) can be treated as the 2D Fourier transform of m(x, y).
Here we introduce an kx −ky spatial frequency domain, which is called k-space. In object
domain, x axis often is called Frequency-encoding (FE) or readout (RO) direction, while
y axis is called Phase-encoding (PE) direction. In designing of the gradient waveforms,
usually the Gy is turned on before each readout line is acquired, so that to create initial
phase differences between each phase encoding direction.
2.1.4 Reconstruction
The reconstruction process in MR imaging refers to how the object can be generated from the
encoded MR signals. From Equation (2.19), a 2D Fourier Transform relationship has been
20

Figure 2.6: Relationship between object space and k-space.
established between s(kx, ky) and m(x, y) on a Cartesian grid. K-space must be sufficiently
sampled to avoid aliasing and achieve the desired resolution:
F OVx = 1/∆kx; F OVy = 1/∆ky; (2.20)
∆x = 1/(2 × kmax,x); ∆y = 1/(2 × kmax,y). (2.21)
See (Figure 2.6) that shows a detailed relationship between k-space and object domain.
If considering the relaxation, the signals that are acquired on the Cartesian line will have
an additional T2 (and off-resonance) weightings over the time. Depending on the ordering
to sample the 2D Cartesian k-space, and different types of sequence, different contrasts
(depends the different tissue properties T1 and T2) will be generated additional to proton
density weighted images. If the designed readout time is comparable to T2, the signals will
experience an non-negligible decay, which doesn’t follow the Fourier transform relationship.
In these cases, more advanced sampling strategies have been developed over times, with
targets on resolving related artifacts [42].
21

Figure 2.7: Examples of non-Cartesian trajectories.
2.2 Data sampling strategies
2.2.1 Non-Cartesian sampling
Although 2D Cartesian sampling is the most used in the clinic and convenient overall, there
are advantages in non-Cartesian samplings for certain applications. Among these, spiral
trajectory and Echo-planar imaging (EPI) are often used in 2D imaging that travels the
k-space faster with more samples given the same amount of time compared to Cartesian.
Therefore, the sampling efficiency is higher. However, because of the longer readouts per
excitation, they are more sensitive to T2 decay and off-resonance will accumulate over the
readout time. Radial trajectory is slower but robust to motion and flow. When non-Cartesian
sampling is used, a process known as gridding [43] is required before discrete 2D Discrete
Fourier Transform can be applied. See Figure 2.7 for demonstrated several typical nonCartesian trajectories. In this thesis, spiral trajectories have been utilized in variant of
applications, e.g. speech and lung.
2D Spiral. Spiral sampling typically begins by acquiring data from the center of k-space
and gradually moves outward along a spiral path, allowing it to cover the entire k-space with
one or a few passes to achieve Nyquist sampling. It provides higher efficiency compared to
Cartesian,as a large fraction of k-space is acquired during each TR. This technique helps
22

reduce motion artifacts, as it naturally oversamples the central region of k-space. Spiral
sampling is particularly well-suited for advanced reconstruction methods like compressed
sensing, especially when combined with strategies such as undersampling or a golden angle
scheme, similar to radial sampling. Because of these reasons, spiral sampling is commonly
used in real-time MRI, where the ability to capture fast-moving objects is essential, such
as in cardiac or speech production imaging. However, longer spiral readouts are sensitive
to off-resonance effects, which can cause signal loss and blurring in reconstructed images,
especially at high field strengths. At low field strength, because of the reduced susceptibility
differences between the air-tissue boundaries, the spiral design becomes more flexible, and is
well-suited for dynamic imaging in vocal tract and lung. Figure 2.8 shows simple 2D spiral
design for bSSFP imaging with selective excitation. Same spiral readouts can be adaptive
to other types of fast gradient echo sequences.
3D Stack-of-Spiral. An natural extension from 2D imaging to 3D imaging is possible
by adding an additional kz encoding realized by Gz gradient. Figure 2.9 illustrates representative examples of 3D stack-of-spiral (SOSP) trajectories extended from 2D spiral, along
with the sequence diagram. Here an non-selective RF pulse is used, therefore eliminating
the previous Gz gradients along with excitation in 2D spiral design. In an SOSP acquisition,
each plane is sampled using one or multiple spiral interleaves and the through-plane direction
is sampled with multiple phase encoding steps.
There are variant ways in designing 3D SOSP. Rotated spiral interleaves [44, 45] and
randomized sampling [46, 47] in k-t space with variable density [48, 49, 50] have previously
proved beneficial in volumetric imaging for different applications. For example, Deng et al.
[45] utilized rotated stack-of-spiral sampling to mitigate aliasing artifacts. Nayak et al. [47]
utilized random sampling in k-space to reduce coherence of undersampling artifacts. Liao et
al. [49] utilized increased sampling density in spiral trajectories to mitigate motion artifacts
in Cine MRI. Lee et al. [50] investigated applying variable density stack-of-spiral sampling
in limb perfusion to significantly reduce the scan time. Overall, this sequence is generally
23

Figure 2.8: Illustration of 2D Spiral bSSFP. Selective excitation is used with slice-selection
Gz gradient turned on paired with gradient prewinder and rewinder. Gx and Gy represents
spiral gradient waveforms as readouts. ADC represents the period when the signals are
received.
24

Figure 2.9: Illustration of 3D Spiral bSSFP. (left) Stack-of-spiral sampling [Adapted from
Ref. [51]] and (right) sequence diagram. Non-selective excitation is used with a hard RF
pulse. A pair of Gz gradients are inserted to encode the kz direction. Gx and Gy represents
spiral gradient waveforms as readouts. ADC represents the period when the signals are
received.
useful for volumetric imaging applications.
2.2.2 Reconstruction with undersampled Data
2.2.2.1 Parallel imaging
Parallel imaging is a widely-used technique for acceleration of scan time with multiple receiver
coils and their known placement and sensitivities [52]. With the help of an array of coils with
each of them having an unique spatial sensitivity, it is possible to recover the missing data
using an accelerated acquisition process, such as skipping some phase-encoding steps, instead
of having aliasing artifacts. There are two types of parallel imaging techniques. First type
utilizes the sensitivity directly for fast MRI in the image domain, such as Sensitivity encoding
25

(SENSE) [53]. It often requires an estimation of the receiver coil sensitivities to solve an
inverse problem. Another type is indirectly using the spatial encoding by coil sensitivity
in a system model, with the assumption that the missing k-space data can be recovered
by interpolation. Methods of this kind include generalized autocalibrating partially parallel
acquisitions (GRAPPA) [54] and more [55, 56].
2.2.2.2 Constrained Reconstruction
The constrained reconstruction is another efficient way to reduce the minimum acquired data
so that to reduce the scan time, with the assumption on the reconstructed image properties.
”Constrained”, means a prior information, or a prior models or bounds to place on the
images to be constructed. These constraints often take the form of regularization terms in
the objective function used for reconstruction. The reconstruction process mostly involves
solving the following optimization problem:
S = arg min
S
||UF ES − y||2
2 + λR(S) (2.22)
The first term is called data consistency, and it enforces the image to be reconstructed
(S) remains consistent with the acquired k-space data (y). U representing the sampling
operater, resulting in undersampled matrix, F is the Fourier transform, and E is the coil
sensitivity encoding. The second term is the additional constraints based on certain types
of the prior knowledge of the images. Several types of the constraints are often used. λ is
the regularization parameters that weights the regularizer R(S).
TV regularization. The TV regularization prior assumes sparsity in the gradient of the
object (i.e. that image is piece-wise constant). This prior information will apply to various
degrees to different images and signals, e.g. it will be more applicable to MR angiography
with contrast than to standard anatomical imaging.
Wavelet sparsity. The wavelet sparsity constraint assumes that images are sparse in the
wavelet domain [57]. Under this assumption, constrained reconstruction would remove noise
26

and artifacts from the images which are not sparse.
Temporal finite difference In dynamic imaging, where a series of images are reconstructed
over time, the temporal finite difference constrains are often added to the temporal domain
with the assumption of the redundancy along the temporal dimension. In the case of speech
imaging for example, dramatic motion occurs mostly in the tongue, lips and the airway while
the other parts of the head and neck remain mostly still.
In practice, constraints are often combined. This prior information can be incorporated
into the reconstruction as a sparsity constraint, so that less sampling is required to recover
the original image without undersampling artifacts.
2.3 Field strength dependence
In the current clinically-used MR scanners, mid-to-high field strength scanners (1.5T and
3.0T) have been dominated the market, providing diagnostic qualities which is capable for
MR spectroscopy, and other forms of functional MRI, high speed imaging, and high resolution
imaging. However, cost remains a significant limitation to the wider dissemination of high
field MRI. There are definite cost advantages to the use of lower field MRI. In this section,
physical principles of the field strength dependence of MRI will be reviewed in relation to
image quality.
2.3.1 Equilibrium polarization and SNR
Signal-to-noise-ratio (SNR) in MR imaging is one of the fundamental measures of the image
qualities. Physical parameters, such as the main field B0 is one the basic dependencies of
SNR. A general expression for SNR, ignoring the resistance from the coupling of the coils,
and pre-amplifier noise can be expressed as:
SNR =
ω0B1,xyM0(∆V )
p
4kT(∆f)(Rc + Rm)
(2.23)
27

where ω0 is the resonant frequency, M0 is the amount of nuclear polarization, ∆V is the voxel
volume, k is Boltzmann’s constant, T is the absolute temperature, ∆f is the bandwidth. Rc
and Rm are the resistance of the receiver coils, and of imaged body realized by the receiver
coil.
Based on Equation (2.2), both M0 and ω0 has dependence of the B0, the noise power for
the receiver coil, and body is proportional to ω
2
0
, and ω
1/2
0
. Hence,
SNR ∝
B2
q
0
αB1/2
0 + βB2
0
(2.24)
with α and β being constants. If body is the primary source of noise, SNR ∝ B0; if
coil resistance dominates, then SNR ∝ B
7/4
0
. In field strengths that are typically used in
the clinic (1.5T, 3T), or higher fields strength and sample size scanned by whole-body MR
systems, body noise is dominant. Whereas in low-field imaging and small-volume imaging,
noise from the receiver coil dominants. Overall, the SNR dependence on main field strength
B0 doesn’t follow the simple linear relationship, but varies with the scanner system design
[58].
2.3.2 Susceptibility effects, T
∗
2 and intravoxel dephasing
There are nonidealities and constrains in the MR imaging that need to be considered, such
as off-resonance conditions and T2 relaxation. In the off-resonance conditions, three main
sources are 1) main field inhomogeneities, 2) susceptibility-induced field variations, and 3)
chemical shift. I will discuss the field strength dependence of them and how the signals will
be effected.
Susceptibility effects. Susceptibility of a material is a measure of its ability to become
magnetized in an external magnetic field, the susceptibility differences cause inhomogeneities
in the imaging gradient fields. Resonance frequency differences may exist because of sampledinduced B0 variations even in a perfect magnet, due to the differences in bulk magnetic
28

susceptibility χ within the sample. The field inhomogeneities depends on the susceptibility
difference and the geometry, usually represented in terms of parts per million (ppm). The
resultant inhomogeneity is most severe near the tissue boundaries between two materials
with different susceptibilities such as air and tissue. Imaging body regions like abdomen
(intestinal gases), lungs, and vocal tract, relatively sharp local gradient inhomogeneities
may result.
Susceptibility effects is also proportional to field strength B0. For example, the magnetic
susceptibility difference (∆χ) at the air-tissue interface is around 8 ppm. At 3T, this leads
to 1022Hz of inhomogeneities. While at 0.55T, this results in 187.4 Hz differences. The field
inhomogeneities caused by susceptibility is linearly proportional to the main field strength.
T
∗
2
. Because the field inhomogeneities exist and introduce space-dependent resonance
offsets that lead to image distortions and artifacts, the received signals can be expressed as
s(t) = Z
vol
m(r)e
−t/T2(r)
e
−i∆ω(r)t
dV
=
Z
vol
m(r)e
−t/T ∗
2
(r)
dV.
(2.25)
The term e
−i∆ω(r)t
considers the local field inhomogeneities which can be potentially
caused by nonuniform B0 and susceptibility, which leads to a distribution of frequencies over
the object that leads to the loss of phase coherence of the spins. This results in a more
obvious amplitude decay than T2 processes and results in phase errors in signals. Thus, a
T
∗
2 decay is introduced to describe this phenomenon, where T
∗
2
is the effective reduced time
constant. T
∗
2
is always shorter than T2. In MR sequences, the spin-echo sequence is designed
to overcome the phase dispersion due to the T
∗
2
, thereby robust to off-resonance.
With the increased field strengths, the shorter of the T
∗
2
. For example, lung parenchyma
T
∗
2
is increased from less than 2ms at 1.5T [59] to roughly 10ms at 0.55T, due to the reduced
susceptibility effects.
Intravoxel dephasing. Field inhomogeneities give rise to both amplitude and phase
29

aberrations that affect image qualities in a complex manner that dependent on the specific
gradient waveforms and designed k-space trajectories. In MR sequence, when the readout
time is large than T
∗
2
, this results in a significant signal decay over the integral of the voxel to
be imaged. We introduce a phenomenon - ”Intravoxel dephasing”. Within a single voxel (a
small volume element in the image), different magnetic spins lose their phase coherence due
to variations in the magnetic field strength within that voxel, leading to a loss of signal and
potential image degradation, particularly near tissue boundaries or areas with significant
magnetic susceptibility differences.
When the inhomogeneity-induced frequency causes slowly varying over phase, signal loss
from dephasing will be minimal if voxel size is small. Thus, the signal strength may actually
increase by imaging at finer spatial resolution. Since the T
∗
2
is longer with decreased field
strength, the intravoxel dephasing may be less pronounced, resulting in the recovered of
certain signals. At lower field strength, intravoxel inhomogeneity will be less pronounced in
regions such as lungs, where there exists relatively large susceptibility-induced variations over
a small spatial range due to the many air/tissue interfaces. As a result, the lung parenchyma
signals will be recovered at low field strengths, compared to high field strength.
2.3.3 Concomitant field effects
Concomitant fields, also known as Maxwell fields, exist whenever the linear gradients are
active, with non-linear spatial dependence result. This is a consequence of Maxwell’s equations for the divergence and curl of the magnetic field. According to [60], the concomitant
field to the lowest order can be expressed as:
Bc(x, y, z, t) = 1
2B0
{G
2
x
z
2 + G
2
y
z
2 + G
2
z
x
2 + y
2
4
− GxGzxz − GyGzyz}. (2.26)
This assumes the symmetry of x and y gradient coils, with a 90◦
rotation between them.
From the expression, we observe that the concomitant field terms are inversely proportional
30

to the main field strength B0. The concomitant fields are increasingly important at low field
strengths. It results in unwanted phase accumulation during the acquisition, and excitation
(zero if isocenter). At lower field strengths, with the high performance gradient systems, if
the imaging is acquired using long readout time (>10ms), the concomitant field effects will
lead to image artifacts like blurring or phase shifts due to their spatial variations.
The effects of concomitant fields [61] have been addressed for image acquisition and
partially addressed for RF excitation [62, 63]. The impact of concomitant fields on image
acquisition at low and mid field strengths has been demonstrated by several groups [64, 65,
66, 67, 68], and includes spatial blurring in spiral imaging and spatially-varying shifts in
echo-planar imaging (EPI) [60, 69]. Since concomitant fields behave as another source of
off-resonance, their effects on the excitation profile during RF pulse design will be similar
to their effects on images during image acquisition (e.g., spatially-varying local blurring
for spiral trajectories). This is an important consideration especially at low field strengths,
strong gradient fields, off-isocenter excitation, and for long excitation k-space trajectories. To
partially address this, Nielsen et al. [62] modeled the concomitant fields as additional phase
terms during RF pulse design. This provided a solution for 3D RF pulse applications with a
large field of excitation and long pulse durations. Since this work does not provide a detailed
derivation of additional phase terms caused by concomitant fields and lacks experimental
validations, our work is motivated to provide a more complete solution.
2.4 MRI of Speech Production
Real-time MRI (RT-MRI) has emerged as one of the most powerful tools for studying
human speech production. It has allowed researchers to demonstrate both functional and
morphological aspects of the vocal tract during speech [70] and to understand the complex
spatiotemporal coordination of upper airway structures in motion [13, 71]. Speech production
RT-MRI has benefited from early adoption of technology for rapid and SNR-efficient imaging,
31

such as custom RF coils [72, 73], parallel imaging [74, 75], advanced data sampling [76,
77], constrained reconstruction [76, 78], and automated post-processing [79, 80, 81]. These
technologies have continued to provide new insights to linguistics research.
2.4.1 Imaging requirements
To image different dynamic movements of the speech production, it needs different spatial
and temporal resolution for various speech tasks. Figure 2.10 illustrates the detailed imaging requirements from the 2014 Speech MRI Summit. Here, each rectangular zone with an
approximate boundary represents a specific speech task in terms of the spatial and temporal
resolutions. The speech MRI generally requires high spatial and temporal resolution, although specific imaging parameters would be dictated by the vocal tract regions and speech
tasks of interest. It is recommended in the speech MRI community that an in-plane spatial
resolution of no more than 3.5 mm2 and time resolution of below 70 ms are required to
study very fast articulatory movements such as those during consonant constrictions and
co-articulation events that are major tasks of interest.
2.4.2 Imaging considerations
Trade-offs. Trade-offs exist such as choosing an appropriate trade-off between spatial and
temporal resolution or between signal-to-noise ratio (SNR) and sampling artifacts. Several
advanced sampling and reconstruction approaches have been explored to improve on this
tradeoff [82, 83, 84]. Fu et al. [82] used sparse sampling and low-rank modeling to achieve
166 frames per sec. This approached used navigators and leveraged repetition of consonantvowel utterances during each data collection. Burdumy et al. [83] utilized golden-angle stackof-radial sampling and constrained reconstruction to achieve 1.3s temporal resolution. This
approach was validated in static phantoms and during sustained and sung phonation. Lim et
al. [84] used 3D stack-of-spiral sampling with spatio-temporally constrained reconstruction
to achieve full vocal tract imaging with 2.4 × 2.4 × 5.8 mm3
spatial resolution and 61 ms
32

Figure 2.10: Spatiotemporal resolution requirements for various speech task [Adapted from
Ref. [71]].
temporal resolution. Lim et al. also demonstrated the value of 3D volume information to
visualize tongue grooving and doming during consonants /s/ and /l/. This did not require
any repetition of the utterances.
Pulse sequences. For speech imaging, there are multiple MR sequences that can be
used to capturing different characteristics of vocal tract. Figure 2.11 shows some examples.
(a) 2D real-time MRI of the vocal tract of the mid-sagittal slice are shown, with frame rate
83 frames per second during production of a comprehensive set of scripted and spontaneous
speech material, synchronized audio recording is shown at bottom; (b) T2-weighted volumetric images are acquired at rest position to capture the details of anatomical characteristics
of the vocal tract. (c) shows examples of 3D real-time volumetric images of the vocal tract,
showing in both mid-sagittal and coronal views.
There are two rapid gradient echo sequences that are widely used for real-time speech
imaging. Spoiled gradient recalled echo (GRE) and Balanced steady-state free precession
(bSSFP) sequences. Both sequences can be designed with very short readout times with
33

Figure 2.11: Examples of MR speech images. a) 2D RT-MRI, real-time magnetic resonance
imaging [Adapted from Ref. [72]]; b) T2-weighted static images; c) 3D RT-MRI, real-time
magnetic resonance imaging [Adapted from Ref. [84]].
efficient sampling, i.e., spiral, together enabling the temporal resolution to be high, and are
well-suited for dynamic imaging in general. GRE sequences provides T1-weighted contrast
with a short TR. bSSFP provides higher SNR efficiency than the spoiled sequences and
T2/T1 contrast, which is advantageous in applications such as in cardiac imaging because of
the excellent blood–myocardium contrast, although its sensitivity to off-resonance manifests
as banding artifacts, limiting its usage at higher field strength.
Field strengths. Imaging at higher field strength can provide high SNR. Prior work
has used RT-MRI at 1.5 T and 3 T, typically using spiral GRE at 1.5 T [85, 74, 86, 72]
and 3T, [87] or radial GRE at 3T [76, 88]. The main challenge is susceptibility at air–tissue
interfaces, which increases linearly with field strength and can cause blurring and signal loss
at speech articulator boundaries. To compensate for this, current speech RT-MRI studies are
most often conducted using very short readouts (2.5 ms) on 1.5T commercial MRI scanners.
MRI at lower field strength, equipped with high performance gradients, offers the advantage
34

of reduced off-resonance, enabling longer spiral readouts [89] with less blurring artifacts.
Speech imaging at low field strength can also benefit from balanced (bSSFP) feasible, which
at higher field strengths is avoided due to banding artifacts in the regions of interest. The
latter is especially important because bSSFP provides superior SNR efficiency compared
with traditional GRE. Recently, Lim et al. [90] has demonstrated speech-production realtime MRI using a contemporary 0.55T system, and to identify opportunities for improved
performance compared with conventional field strengths.
2.4.3 Unmet needs
2D mid-sagittal RT-MRI has had a particularly profound impact on our understanding
of speech production [72, 76, 77, 91, 92]. It is now a standard and mature method. A
typical setup can achieve 2.4 × 2.4 mm2
spatial resolution and 12 ms temporal resolution
[72]. However, there are some limitations. Visualization is limited to a single imaging
plane (typically mid-sagittal), which is inadequate for some articulations [74, 93], such as
/l/ that have critical maneuvers off the mid-sagittal plane. Asymmetric behavior of the
tongue body, which is only recognized from axial views, and volumetric properties critical
for generating the acoustic resonance structure of speech have also generated researcher
interest. These shapes carry important functionality for understanding speech acoustics.
Interleaved multiplanar imaging [91] has been achieved with 18-36 ms temporal resolution,
and this provides additional but still incomplete information. For these reasons, 3D RTMRI is desirable. However, finer spatial and temporal resolution are desirable, and remain
an unmet need in linguistics. This will be beneficial for understanding the dynamics of vocal
tract shaping, and the associated acoustics.
35

2.5 MRI of Lung Ventilation
Lung ventilation measures the gas exchange in and out of the lungs. Lung MRI can provide
structural and functional information for the screening, diagnosis, and longitudinal assessment of lung diseases without ionizing radiation [94, 95, 96, 97, 98]. MRI lung ventilation
has been explored for years [99, 100]. Hyperpolarized gas [101, 102] like 129Xe and 3He have
been used as contrast agents with single or multiple breath holds to acquire airway and ventilation images. However, this requires special equipment with high cost and scarcity, and
difficult to perform in patients who have pulmonary diseases. Oxygen enhanced MRI [103,
104] has been applied to quantify regional perfusion and ventilation. Free breathing 1H MRI
techniques offer a contrast-agent free methodology of assessing pulmonary ventilation by
measuring local changes in proton density following cardio-respiratory motion. Researchers
have been explored contrast agent-free 1H MRI for functional lung imaging at 1.5T and 3.0T
[105, 106]. These methods together improve the diagnosis and evaluation of patients with
pulmonary dysfunction.
2.5.1 Imaging considerations
Sequence design. The requirements for developing lung MR imaging requires high resolution, free-breathing with high SNR efficiency, resulting in motion-resolved images that
captures the dynamic movements of the lung (parenchyma, vessels) details. Structural lung
imaging at 0.55T has been recently developed and applied to different pulmonary diseases
by several research groups. Stack-of-spiral UTE has achieved 1.75mm isotropic resolution
in a 15.5-min acquisition with sufficient parenchyma SNR and sharp features [107]. This
technique has been improved using iterative concomitant field correction [108]. Balanced
steady-state free precession half radial dual-echo (bSTAR) [20, 109] lung imaging leverages the high SNR efficiency of bSSFP at 0.55T with very short TR (≤2.14ms). bSTAR
has been used to achieve 0.9mm isotropic resolution within a 13-min free-breathing scan,
36

Figure 2.12: Examples of MRI of lung ventilation. (a-c) Hyperpolarized 3He 129 Xe images of one healthy volunteer, chronic obstructive pulmonary disease (COPD) and asthma
[Adapted from Ref. [114]]. (d,g) Signal enhancement maps using Oxygen-enhanced MRI
at 0.55T in one healthy volunteer and patient with lymphangioleiomyomatosis (LAM)
[Adapted from Ref. [115]]. (e-f) Fractional ventilation maps using Phase-resolved Functional Lung (PREFUL) MRI at 0.55T in one healthy volunteer and patient with long
COVID [Adapted from Ref. [116]].
providing superior parenchyma details. The oversampling near the center of k-space and
relative inefficiency of the 3D radial trajectory is the reason for the long acquisition, which
can be impractical in some clinical settings. Recently, a stack-of-spiral out-in (SOS out-in)
bSSFP variant has been explored in the context of breath-hold, and provides superior vesselto-parenchyma contrasts [110]. In general, the spiral-based trajectory offers more efficient
sampling and flexible TR selection compared to the radial trajectory. A premise of this work
is that free-breathing SOS out-in is suitable for the evaluation of pulmonary function and
potentially more practical [111, 112, 113].
Field strengths. Proton-based lung MRI suffer from low signal-to-noise ratio (SNR)
at conventional field strengths (≥1.5T) [59]. High-resolution ultrashort echo time (UTE)
sequences, often combined with 3D non-Cartesian trajectories, e.g., radial sampling [117,
37

Figure 2.13: Example T2-weighted magnetic resonance imaging (MRI) images obtained
with a conventional MRI system (1.5-T MAGNETOM Aera; Siemens Healthcare) and the
new high-performance low-field MRI system (prototype 0.55-T MAGNETOM Aera; Siemens
Healthcare) from a patient with lymphangioleiomyomatosis (LAM). The visibility of lung
tissue using high-performance low-field MRI enables assessments of regional function and
tissue characterization [Adapted from Ref. [21]].
111, 118], have been utilized at 1.5T and 3T to overcome the SNR loss, and provide robustness to motion. These technical developments opened the door for MRI-based quantitative
evaluation of lung function [119, 120] and pulmonary diseases [121, 122, 123]. Contemporary
0.55T MR systems provide numerous advantages for pulmonary MRI (see Figure 2.13), and
are well-suited for structural and functional imaging [89, 124, 125]. The 0.55T field strength
provides higher parenchyma signals compared to conventional field strengths (≥1.5T) due
to the reduced susceptibility differences between air (alveoli) and tissue (parenchyma and
blood) [126]. Lung parenchyma T
∗
2
is roughly 10ms at 0.55T, compared to ≤2ms at 1.5T
[59, 127], resulting in substantially higher SNR. Moreover, contemporary gradient hardware
enables developing fast imaging based on spiral [90, 2, 107] and echo-planar imaging (EPI)
[128] that can improve acquisition efficiency.
38

2.5.2 Unmet needs
Pulmonary function testing remains the primary clinical routine for lung disease evaluation
[22], but it only provides global ventilation. The dependence on patient effort also limits
its use in young children and infants [129]. For years, Fourier decomposition (35), pencil
matrix decomposition [130], and phase-resolved functional lung (PREFUL) [131, 132, 133,
134] have been developed to resolve ventilation/perfusion (V/Q) weighted lung images simultaneously. These intensity-based quantification methods rely heavily on sufficient lung
parenchyma SNR, which are of limited use at conventional MRI (≥1.5T) due to the short
T
∗
2
and low proton density. PREFUL generates the full cycle of cardiac and respiratory
phase-resolved images based on the sorted signal curves extracted from registered image series. This method, combining with free-breathing 0.55T MRI, has been shown pulmonary
dysfunction in patients with persistent symptoms after COVID-19 and in adults [135] and
pediatric patients [113]. However, PREFUL is limited in 2D, offering a relative percentage
of the V/Q measurements. Therefore, an efficient, repeatable, and accurate assessment of
free-breathing 3D lung ventilation mapping method is needed, benefiting from the improved
parenchyma signals at low field strengths.
39

Chapter 3
Improved 3D Real-Time MRI of Speech
Production
In this study, we explore improvements in achievable spatio-temporal resolution for 3D RTMRI, based on the method demonstrated by Lim et al. [84]. This work utilized spiral
sampling along kx − ky (sagittal) plane, and linear phase encode order along kz (left-right).
Our hypothesis is that the spatiotemporal resolution can be further improved by altering
the data-sampling approach.
3.1 Methods
3.1.1 Experimental Methods
Experiments were performed on a commercial 1.5T scanner (Signa Excite HD, GE Healthcare, Waukesha, WI) with gradients of 40 mT/m strength and 150 mT/m/ms maximum slew
rate per axis. Experiments used a body coil for RF transmission and a custom eight-channel
upper airway coil for signal reception. Acquisition and reconstruction were implemented
using a real-time interactive imaging platform (RT Hawk, Heart Vista Inc, Los Altos, CA)
[136]. We used a spoiled gradient-echo stack-of-spiral (SOSP) trajectory with imaging parameters: spatial resolution = 2.4 × 2.4 × 5.8 mm3
, FOV = 200 × 200 × 70 mm3
, TR/TE
40

= 5.05/ 0.68 ms, and readout duration = 2.52 ms. Spiral readout was used due to its scan
efficiency and prior successful applications to speech RT-MRI [71, 72, 77, 137]. Shimming
is critical, and our protocol uses auto-calibration, followed by manual adjustment of the
center frequency by the scan operator to minimize visible blurring of the tongue boundary.
The imaging protocol was approved by USC’s Institutional Review Board, and all subjects
provided written informed consent. Subjects, while being imaged, read stimuli presented by
a mirror projector setup [72].
3.1.2 Data Sampling
Figure 3.1 illustrates sampling patterns. We consider five different (k,t) data sampling
patterns (Cases I - V), all of which are based on the SOSP sequence. Case I is identical to
Lim et al. (19)[84], and is taken as the baseline method. The selection of sampling patterns
enables pair-wise comparison that elucidate benefits of each specific data sampling change.
Case I utilizes a linear temporal order along kz with a constant spiral angle (CA) in the
kx − ky plane. After kz is fully sampled, the spiral angle is incremented by the golden angle,
θGA = 2π×2/(
√
5+1) (4,14,15) [72, 91, 92]. The increment angle is reset after N interleaves.
We use N = 34 in this work. Case II uses the same linear temporal order along kz as in Case
I, but the spiral angle is incremented by GA with the kz. Case III introduces the bit-reversed
temporal order (29)[138] along kz while utilizing the same CA strategy as in Case I. Case IV
adopts the bit-reversed temporal order along kz but with the GA increment in kx −ky. Case
V utilizes GA increment in kx − ky and introduces variable density (VD) random temporal
order along kz. Bit-reversed temporal order (29) [138] maximizes the time between acquiring
adjacent indices to mitigate motion artifacts. Although this is a coherent order, it provides
similar benefits to random temporal order. VD produces more frequent sampling of low
spatial kz spatial frequencies; specifically, the density along kz is implemented as:
p(kz) = 1
|αkz|
(3.1)
41

Figure 3.1: Data sampling patterns. Five sampling strategies are compared. I-V show kz
vs. time plots. Each dot represents one spiral arm, where the color represents its initial
angle in the kx − ky plane (see color disc). Case I. kz is in a linear order, and the golden
angle increment is applied for every 12 - TRs (full kz encoding). Case II. The linear order
along fully sampled kz is maintained while the golden angle (GA) in the kx − ky plane is
increased for each TR. Case III. Constant-angle remained in the kx − ky plane for 12-TRs;
a bit-reversed interleave temporal ordering is obtained in the kz. Case IV. Golden angle
increment is applied for each TR, while the temporal order of the kz is bit-reversed. Case V.
kz is acquired in a randomized and variable density fashion, and the golden angle increment
is applied every 1-TR.
42

where kz is normalized from -0.5 to 0.5. The scale factor α =
PNz
2 −1
n=0
4(Nz−1)
1+2n where Nz is
the number of kz steps (slices). We use Nz = 12 for all experiments.
3.1.3 Image Reconstruction
Image reconstruction was performed by solving the following constrained optimization:
arg min
f(r,t)
||A(f) − b||2
2 + λs||TVs(f)||1 + λt
||FDt(f)||1 (3.2)
where f(r, t) represents the dynamic images to be reconstructed, A is the encoding function
(including coil sensitivity and non-Cartesian Fourier transform), b is multi-coil k-t space
measurement data, the vector r ∈ (x, y, z) represents spatial coordinates, and t denotes
time. T Vs represents isotropic 3D spatial total variation, and F Dt represents first-order
temporal finite difference. λs and λt are corresponding spatial and temporal regularization
parameters, respectively.
This optimization problem was solved using the alternating direction method of multipliers (ADMM) algorithm [139], as implemented in the Berkeley Advanced Recon-struction
Toolbox (BART) [140]. Coil sensitivity maps were generated using Efficient iTerative Selfconsistent Parallel Imaging Reconstruction (ESPIRiT) [141] and were assumed to be timeinvariant. Reconstruction is fully 3D in this work, which allows for greater flexibility in
data sampling but carries substantially higher memory requirements. We divided the reconstruction into several time segments (100 time frames per segment) to manage the memory
requirement. Note that prior work [84] applied inverse Fourier transform along kz and exploited nonlinear reconstruction for each x − y slice separately, which was possible because
kz was always fully sampled within each temporal frame. All sampling patterns were reconstructed with a temporal resolution of 61ms per frame (12 TRs). We additionally explored
retrospective selection of temporal resolution for Case V by reconstructing images with temporal windows of 30.5 ms per frame (6 TRs) and 15.25 ms per frame (3 TRs).
43

3.1.4 Regularization Parameter Selection
We analyzed the sensitivity of regularization parameters λs and λt
in the constrained reconstruction using data acquired from Speaker 1 (31-year-old, Male, native Chinese speaker,
English as second language). This speaker was scanned while reading the English stimuli “/loo/-/lee/-/la/-/za/-/na/-/za/”, repeated twice at a natural rate. Regularization parameters λs and λt were chosen visually based on image quality in mid-sagittal views and
time-intensity plots by three academic linguists with 5 decades of experience using speech
RT-MRI. We performed parameter sweeps in two stages. We first performed a coarse sweep
spanning 6 orders of magnitude in log scale (λs = 0, 0.001, 0.1. λt= 0.0001, 0.01, 1.) and
then tested a finer region (λs= 0.001, 0.002, 0.004, 0.006,0.008, 0.010. λt= 0.005, 0.01, 0.02,
0.03, 0.04).
Table 3.1: Stimuli design I. Alternating low and high vowels with no intervening lingual
consonants. These stimuli were recorded at both a normal and speeded (roughly 1.5×)
speaking rate.
3.1.5 Experimental Optimization
We compared all of the candidate sampling approaches in a prospective experiment with
highly repeatable speech tasks produced by one volunteer (Speaker 2: 18-year-old, Male,
44

native speaker of American English). The subject was scanned with all five aforementioned
3D RT-MRI data sampling schemes, as well as interleaved three-slice 2D RT-MRI [91]. The
imaging parameters can be found in [91]. Table 3.1 lists the speech tasks. Each sentence
was spoken twice, once at a “natural” rate and once at a speaking rate of approximately
1.5× the initial rate, denoted “speeded.” This stimuli set was designed to test image quality
when we expect large signal fluctuations.
Results from the five candidate sampling patterns (Figure 3.1) were compared using
intensity vs. time plots of sagittal, axial and coronal views. The best sampling scheme was
chosen based on boundary sharpness and visual clarity of the dynamic articulators. This
(Case V) is denoted the proposed method for the subsequent experiments.
3.1.6 Speech Production Experiments
We applied the best performing sampling and reconstruction approach to a broad set of
speech stimuli that were chosen for assessing image improvements. The stimuli were selected
to elicit sweeping movements of the tongue body used to produce relatively wide vocalic
airway and rapid movement of the tongue tip coming into contact with the alveolar palate
to produce a consonant constriction.
Two subjects participated in the study of the evaluation for speech application. We
evaluated the results from stimuli set 1 (Speaker 2, Table 3.1) using the original method,
the proposed method, and the reference 2D interleaved method. These stimuli were recorded
at a normal speaking rate and at the speeded speaking rate. The stimuli set in Table 3.2
was recorded for Speaker 3 (20-year-old, Male, native speaker of American English) who
was scanned with the original and proposed method, as well as the reference interleaved
three-slice 2D RT-MRI [91] at normal speech rates.
45

3.1.7 Speech Stimuli
The stimuli set in Table 3.1 contains four sentences that were designed to elicit an alternating sequence of phonetically low (/æ/, /a/) and high, front vowels (/i/, /I/), with no
intervening or bordering lingual consonants. These vowel sounds are made with tongue postures that are for low vowels retracted and depressed in the mouth, and for high vowels,
bunched and domed in the oral palatal cavity. Thus, moving from one type to the other
(hi to low or low to hi) requires the speaker to produce sweeping lingual movements exhibiting large displacements. This allows for the examination of large lingual movements in
real, natural speech over a significant spatiotemporal span within the functional vocal tract
space. Specifically, in this dataset’s sagittal views, the /a/ sound is characterized by tongue
retraction producing a narrow constriction of the tongue root in the pharynx. And /æ/ is
articulated similarly but with a somewhat flatter tongue surface profile and slightly lesser
pharyngeal retraction. For high vowels, /i/ is characterized by a relatively narrow constriction in the oral cavity in close proximity to the hard palate, and /I/ behaves similarly with
an only slightly less constricted posture. These sentence stimuli also contain bilabial stop
(closure) consonants (/b/, /p/) intervening between the vowels, during which the lips are
approximated, then compressed, and then released, and /h/s, which are articulated without
any supralaryngeal constriction. These labial and glottal flanking consonants were chosen so
as to minimize any coarticulation with the target vowels by ensuring that the neighboring
consonants are not articulated with the tongue. Results were compared in mid-sagittal and
axial views at the time when specific vowels were articulated.
The stimuli set in Table 3.2 was designed to examine alveolar consonant segments that
are created with a rapid upward action and closure of the tongue tip articulator; these are
among the fastest speech sounds [142], taking place more rapidly than the vowels articulations
elicited in the first stimuli set (Table 3.1). Fifteen sentences having two different stress
patterns were used to elicit the American English alveolar tap/flap consonant [R] (in the
46

Table 3.2: Stimuli design II. Intrinsically fast tongue-tip consonants.
falling stress pattern) and eliciting the alveolar consonant [t] (in the level stress pattern [no
flapping]). These alveolar consonants are articulated by a single constriction action of the
tongue tip, and the tap [R] (e.g., the medial consonant in the word “pita”) can be understood
as a very rapid stop consonant gesture. The canonical stop consonant [t] (e.g. the wordfinal consonant in “peat happily”), while still a rapidly articulated tongue tip consonant,
produces a somewhat longer closure period with greater lingual contact and compression at
the alveolar ridge than for the tap. In addition to word-internal and word-final taps (e.g.
“pita” versus “Pete a”), two other tap variants are included. A rhotic “flap” is elicited
(e.g., in the word “Peter”) in which the tongue tip constriction has a retracted posture and
trajectory, yielding a significant sublingual cavity. Lastly, a laterally released tap is elicited
(e.g., in the word “petal”), in which the tongue profile is narrowed/stretched front-to-back
(assisted by tongue rear retraction) such that the sides of the tongue lower, potentially
pulling away from the palate, so as to allow the lateral airflow required for the following [l]
consonant.
47

3.1.8 Evaluation
Results were compared using intensity vs. time plots during the temporal interval from
the end of the vowel /ei/ in the frame sentence’s “gave” to the /a/ or /æ/ vowel in the
frame sentence’s ‘poppy’ or ‘happily.’ The mean value of Regions-Of-Interest (ROIs) vs.
time was used to visualize the clarity of detected movements [143]. A 4 × 4 pixel ROI was
placed in the mid-sagittal views focusing on tongue gesture excursion during its maximum
consonantal constriction. Intensity normalization among 3D and 2D techniques was done by
taking the reference of middle of the tongue region. Time alignments of different methods
were performed manually based on the time frames with narrowest constriction, e.g. the
tongue tip contacting the alveolar ridge. Note that perfect time alignment across natural
tokens cannot be expected or achieved due to normal variations in speech rate, resulting in
a drift exhibiting slight temporal misalignments later in time from the constrictions’ anchor
point.
48

Figure 3.2: Results of parameter sweep during production of “/loo/ - /lee/ - /na/ - /za/ -
/la/ - /za/.” From top to bottom, spatial TV penalty λs ranges from 0, 0.001 to 0.1. From
left to right, temporal FD penalty λt ranges from 0.0001, 0.01 to 1. A mid-sagittal view (left)
and intensity-time plot (right) with sagittal reference (white dashed line) are shown for each
parameter setting. Spatial and temporal blurring is observed as λs and λt are increased,
respectively. Based on the spatiotemporal visualization of articulator movements, we believe
λs and λt should be chosen in the magnitude of 0.001 and 0.01, respectively (blue dashed
boxes).
3.2 Results
3.2.1 Regularization Parameter Selection
Figure 3.2 illustrates the selection of regularization parameters λs and λt
in Eq. [3.2]
for Case I. The spatial regularization term (TVs) controlled by λs provides denoising as
visualized in the sagittal plane but results in excessive spatial smoothing when large (e.g. λs
= 0.1). The temporal regularization term (FDt) controlled by λt suppresses undersampling
artifact and recovers intensity changes but results in excessive temporal blurring when large
(λt = 1), as best visualized in the intensity-time plots. Based on a consensus of three expert
49

readers, we chose λs = 0.008 and λt = 0.03, and applied this to all sampling schemes for the
simplicity of pairwise comparisons in the remainder of this work.
3.2.2 Experimental Optimization
Figure 3.3 presents pairwise comparisons of the candidate sampling methods during a representative utterance. We extracted four different intensity vs. time plots from different
views that exhibit the lips and tongue movements. I & II: We observed sharper articulator
boundaries in the mid-sagittal plane and clearer lips fluctuations and tongue root movements (pink arrows) for Case II than for Case I. Similar improvement was observed when
comparing III & IV (not shown for brevity). II & IV: We observed improved delineation
of tongue body motion in Case IV in both mid-sagittal and coronal views (blue arrows).
Similar improvement was observed when comparing I & III (not shown for brevity). IV
& V: We observed substantial improvement in boundary sharpness and background noise
suppression for Case V. The region marked with red arrows shows that Case V captures the
fast movements of lips and tongue that are also seen on interleaved 2D RT-MRI (not shown).
We concluded that Case V provided the best image quality among the candidate sampling
strategies.
50

Figure 3.3: Pairwise Comparisons. Four intensity vs. time plots (right) are shown for each
sampling strategy. The line locations (left) are shown over sagittal, axial and coronal views
images. Results from four proposed sampling schemes (I, II, IV and V) are compared. I
& II: The notable differences exist in sagittal views. High frequency lip motion and back
of the tongue movement near the pharyngeal walls are better depicted (pink arrows) in II,
II & IV. The identification of tongue body movements in IV is clarified in both intensity
vs. time displays of sagittal and coronal views (blue arrows). IV & V: The motion-induced
image blurry is successfully removed among sagittal, axial and coronal views in V. The fast
motion of articulators (red arrows), e.g. lips & tongue body, is accurately visualized. Case
V substantially outperforms other methods, providing a high temporal fidelity.
Figure 3.4 illustrates retrospective selection of temporal resolution for Case V, using
51

the results of stimuli set 1 (Table 3.1) with a speeded speaking rate. The use of finer
temporal resolution enables capturing rapid opening and closing of upper and lower lips
(green arrows in c), while a relatively better image quality can be preserved with adequate
temporal resolution shown in a. The proposed method provides a flexible choice of a wide
range of temporal resolutions, while a trade-off exists between image quality and temporal
resolution.
52

Figure 3.4: Retrospective selection of temporal resolution. (Left) sagittal view. (Right)
intensity vs. time plot with reference cut (red dashed lines). We reconstructed stimuli
videos using stimuli set 1 at a speeded speech rate with a 61 ms per frame, b 30 ms per
frame and c 15 ms per frame temporal resolution using λs= 0.008, λt= 0.03. The coarse
sampling rate (a) is inadequate to accurately visualize lip movement from the intensity vs.
time profile (green arrows). However, data acquired with high temporal resolution (b and
c) shows clear lip opening and closing events without temporal blur. On the other hand, the
image quality in the mid-sagittal view at the highly accelerated rate (c) is degraded due to
the coarse spatial resolution.
3.2.3 Application to Speech Production
Figure 3.5 compares the results of the original and the proposed data sampling method,
with the reference interleaved multi-slice method. Six snapshots are selected from stimuli 1
(Table 3.1) at normal and speeded speech rates, along with static images in resting state.
53

Both mid-sagittal and axial views are shown. The high (/i/) and low (/a/, /æ/) vowels
shown in Figure 3.5 are monosyllabic and bounded by labial consonants. 1) Static images
at resting state share the same information among original, proposed and reference methods.
2) For the /p/ sound (b), selected at earliest closure of the lower lip, we observed blurring of
the lower lip and jaw with the original method. The proposed method does not have this issue and is consistent with the interleaved 2D reference. 3) For low vowels (c, d) we observed
blurring at the anterior tongue and tongue root with the original method. The proposed
method shows clearer tongue surface boundaries, which corresponds with the 2D reference.
4) For high vowels (e), we observed blurring at the superior surface of the tongue with the
original method. The quality of same region is improved in the proposed method. The reference method shares the same information as the proposed method. 5) For the speeded speech
rate (f, g), the proposed method demonstrates improved sharpness of articulator boundaries
compared to the original approach. See also Supporting Information Video S2 showing
dynamic tongue movement cycles.
54

Figure 3.5: Demonstration of speech application of the proposed method in stimuli set 1.
Static voiceless images (a) are captured at rest. Six sounds (b - g) are displayed for both
normal and speeded speech rates. The results are compared between the original method
(Case I; 1st column), the proposed method (Case V; 2nd column), and the reference multislice 2D method (3rd column). The qualities of static images (a) are comparable among
the three methods. Blurriness around lower lip and jaw (b) is noticeable in the original
method but is mitigated in the proposed method. Severe motion blurring around the tongue
body is observed in the original method in low tongue postures (c, d, f), whereas clear
tongue surface contours are preserved in the proposed method. The anterior tongue surface
boundaries in the proposed method are resistant to blurring artifacts in high tongue postures
in both normal and speeded rates (e, g). The proposed method and reference method share
similar information in sagittal and axial views.
55

Figure 3.6 demonstrates the improvements of temporal fidelity of the proposed method
compared with the original using the stimuli set 2 (Table 3.2). Six examples are shown:
(a) alveolar consonant [t]; (b-d) word-internal consonant [R; (e-f) rhotic and lateral flap
consonants, respectively. Mid-sagittal views (A) at the time of constriction and the intensity
vs. time plots (B) during a temporal interval of 1.21s are shown. Cut 1 and 2 in (B)
demonstrate different directions of tongue movements. 1) Results for the proposed method
in (B) show clear tongue tip articulations (red dashed boxes) during ‘peat’ (a) and ‘Otto’
(c), demonstrating a single contact of the tongue tip to the alveolar ridge. However, this
critical event cannot be detected in the original method. 2) The sharp boundaries reflecting
the retraction motion (yellow arrows) in (B) during the medial consonants of ‘pita’ (b) and
‘Edda’ (d) were observed in the proposed method. However, this is not clear in the original
method. The mid-sagittal results in (A) are blurry around the tongue tip and the blade
during constriction in the original method. 3) In (e), we observed severe blurring artifacts
(red dashed boxes) for the original method in Cut 1 & 2 (B). However, the corresponding
areas in the proposed method show crisp patterns as well as clean backgrounds. This is
consistent with (A) where the boundaries of sublingual cavity are sharper in the proposed
method. 4) Similar enhancements in the proposed method for stimulus ‘beatle’ (f) were
observed in (B), showing the superior resistance to motion blurring artifacts.
56

Figure 3.6: Demonstration of speech application of the proposed method (Case V) in stimuli
set 2. (A). Mid-sagittal views of six selected alveolar consonants at time of constriction, for
the original (Case I; 1st column) and the proposed (Case V; 2nd column) methods. (a) An
example of consonant [t] showing (in proposed method) the contact of the tongue tip with
the alveolar ridge; (b, c, d) three examples from the word-internal group of [R] (“taps”),
having a rapid and brief linguopalatal contact of the tongue tip. Tongue tip contact in these
four cases is clearly visualized in the proposed method, whereas such events are not captured
in the original method, with blurring artifacts around the tongue surface boundaries. (e,
f) Show two other flap/tap examples with variant 3D tongue shaping of a rhotic flap and
laterally released tap respectively. The sharpness of sublingual cavity boundaries in (e) and
of the tongue tip in (f) are well preserved in proposed method as compared to the results
in original method. (B). Intensity vs. time plots with the reference Cuts 1 & 2 marked as
dashed lines in (A). Red dashed boxes highlight the constriction execution, and yellow arrows
denote when the lingual constriction begins its release. Results in the proposed method
show clear articulator motion patterns in the time series, which have clean backgrounds
and are temporally sharp. However, in the original method these same motion patterns are
‘smoothed’, blurred, and less sharp, indicating limited temporal fidelities in these stimuli.
57

Figure 3.7 demonstrates the tongue shaping geometries for rhotic and lateral flaps in
stimuli set 2. Three displays are selected to assess imaging improvements. (A) Intensity
vs. time plots are shown to demonstrate the temporal fidelities in capturing the sublingual
airway cavity in the rhotic flap and release in “Peter”. (B) The tilted axial views are shown
to demonstrate the tongue’s side channels while articulating the laterally released tap in
“bottle”. (C) The mean value vs. time plots are displayed for the original, the proposed
method, and the reference method. In (A), we observed in both Cut 1 and Cut 2 the creation
of a sublingual cavity as the tongue tip quickly rises and retracts (period within the orange
boxes). However, similar regions are blurred in the original method. From (B), clear dark
channels at the sides of the tongue are detectable along the axis from tongue tip to the
retracting back of the tongue. Channel formation for lateral lingual airflow (along one or
both sides of the tongue) is a critical linguistic feature of this consonant production. Panel
(C) presents a region-of-interest analysis of the formation for the sublingual category also
shown in A. A pixel intensity drop in this vocal tract region is observed for the proposed and
the reference methods during the time of the sublingual cavity formation, but the original
method shows a shallow curve with poorly defined onset and end as compared to the sharper
and greater magnitude reflection of the cavity formation in the proposed method. Overall,
results in the proposed method match closely with the reference. In sum, the results illustrate
the capability of the proposed method for capturing more complex geometries associated with
rapid raising and retraction of the tongue tip and narrowing of the tongue profile created by
anterior-to-posterior stretching for the tongue.
58

Figure 3.7: Illustration of sublingual cavity during a rhotic flap and the tongue side channeling during a laterally release tap. (A) Intensity vs. time displays; (B) axial views; (C)
Mean-of-ROI vs. time plots among original (Case I), proposed (Case V) and 2D reference
methods are shown. (A) Horizontal (cut 1) and vertical (cut 2) lines are selected from midsagittal plane (dashed lines). The “dark” sublingual cavity pattern (see orange boxes) from
cut 1 is blurry in the original method, compared with the proposed and the reference. In cut
2 the cavity (orange arrows) is detectable from the proposed and the reference, in contrast to
in the original method in which it is mitigated. (B) Three tilted axial planes were selected
from the mid-sagittal, indicating the stretched tongue sides. Visible dark channels (white
arrows) indicate airflow pathways in the original and the proposed methods. (C) An ROI
(white dashed box) located at the sublingual cavity were defined in the mid-sagittal plane.
A clear intensity drop and recovered pixel intensity curve is observed in the proposed and
reference method. However, the curve is shallower with less defined onset and offset for the
original method, demonstrating the degraded image quality.
59

3.3 Discussion
We demonstrate improved 3D RT-MRI of speech production by using stack-of-spiral sampling
with advanced data sampling strategies and spatio-temporally constrained reconstruction.
We demonstrate application of the improved 3D RT-MRI to capture alternating high-low
vowels in both normal and speeded speech rates in healthy adults. This indicates that the
proposed method provides adequate temporal resolution to capture natural movement of
articulators, e.g. lips and tongue, in speech production studies, needing no speech slowing,
prolonged or sustained speech or repetitions. Further, we show the capability of the proposed
method to visualize fast alveolar consonant segments during natural speech, along with
capturing complex 3D tongue shape geometries. This approach provides better depiction of
rapidly moving articulators and roughly 2-fold finer temporal resolution as compared to the
prior state-of-the-art [84]. Furthermore, 3D speech RT-MRI now provides roughly 10-fold
greater spatial coverage compared to mature 2D speech RT-MRI technology, with desired
capabilities of capturing rapid temporal articulatory events. We believe this will enable
significant new insights in speech science.
The substantial improvements in temporal resolution and spatial coverage stand to allow more complete characterization of complex vocal tract geometries, such as those found
in lingual grooving, channeling, and concavity found in phonetic sibilant, lateral, retroflex,
and other consonants (e.g., /s, l, ô, ú/), as well as the tongue side bracing found for most
constricted speech sounds including mid and high vowels and the labial protrusion critical to
certain vowels (e.g., /u, y/) and consonants (e.g., /ô, w, S/). It allows for the inclusion of rapid
articulations such as those found in consonants that occur many of the world’s languages
such as the taps, trills, and ejective [144] and for clicks (found for example in many southern
African languages) and in vocal performance such as beatboxing [145]. Finally, a thorough
and accurate spatiotemporal characterization of 3D vocal tract shaping is critical for models
connecting articulation to the acoustic resonance and noise structure that it produces and
60

that characterizes the transmitted speech signal. Beyond typical speech production, complete characterization of 3D vocal tract dynamics may prove important for describing and
remediating speech post-surgery, for example in cancer patients having glossectomies.
Off-resonance is an important consideration in speech MRI due to the susceptibility
effect at air-tissue interfaces [126]. Off-resonance effects are mitigated in this study by
imaging at 1.5 T, using a very short (2.5ms) readout duration and manual center frequency
adjustment. This follows our institutional best practice for speech RT-MRI data collection
that has evolved over 15 years and more than 100 subjects (each with more than one-hour
scan time). Off-resonance should be carefully examined on each new imaging platform,
particularly when imaging at higher field strengths (e.g. 3 T) or with longer readout times.
Recent advances in dynamic off-resonance correction for 2D spiral RT-MRI can be adapted
for 3D stack-of-spiral RT-MRI. These include methods based on model-based reconstruction
[146] and convolution neural networks [147, 148].
3.4 Conclusion
We demonstrate improved 3D RT-MRI of human speech production using innovative kt data sampling. We propose a stack-of-spiral sampling trajectory with variable density
randomized temporal order along kz and golden-angle increment in the kx-ky plane sampling
trajectory, combined with spatio-temporally constrained reconstruction. This scheme shows
superior resistance to motion artifacts and outperforms other proposed sampling schemes.
The improvements in spatiotemporal resolution and in ameliorating blurring allow better
visualization of fast lip and tongue movements of both large and small articulatory magnitude
and are successful at imaging even very rapid consonant segments and complex 3D tongue
shaping geometries such as sublingual cavities and lateral channels. Both dynamic speech
production imaging as well as other imaging applications in which dynamic edge information
and/or mitigating off-resonance are critical stand to benefit from this advance.
61

Chapter 4
Free-breathing 3D pulmonary ventilation mapping at 0.55T using Stack-ofSpiral Out-in bSSFP
4.1 Introduction
Lung MRI provides both structural and functional information for the screening, diagnosis, and longitudinal assessment of lung diseases without ionizing radiation [94, 95, 96, 97,
98]. Proton-based lung MRI suffers from a low signal-to-noise ratio (SNR) at conventional
field strengths (≥1.5T) [59]. However, contemporary 0.55T MR systems offer significant
advantages for pulmonary MRI [124, 89, 125]. At 0.55T, lung parenchyma signals are higher
compared to conventional field strengths due to reduced susceptibility differences between
air (alveoli) and tissue (parenchyma and blood) [126]. The lung parenchyma T
∗
2
is approximately 10ms at 0.55T, compared to ≤2ms at 1.5T [59, 127], resulting in substantially higher
SNR. Additionally, contemporary gradient hardware enables the development of fast imaging sequences, such as spiral [90, 2, 107] and echo-planar imaging (EPI) [128], which can
improve acquisition efficiency.
Structural lung imaging at 0.55T has been recently developed to achieve 3D isotropic
62

high resolution under breath hold and free breathing conditions. The balanced steadystate free precession half-radial dual-echo (bSTAR) [109, 20] lung imaging leverages the
high SNR of balanced steady-state free-precession (bSSFP) at 0.55T with very short TR
(≤2.14ms). Under free breathing, bSTAR has demonstrated superior parenchyma details
(0.9mm isotropic resolution), but with relatively long scan times (13-min). Acquisition time
can be shortened by adopting k-space trajectories with higher SNR efficiency. Recently, a
stack-of-spiral out-in (SOS out-in) bSSFP variant has been explored to provide superior
vessel-to-parenchyma contrasts for structural evaluation under breath-hold conditions. The
spiral-based trajectory offers more efficient sampling and flexible TR selection compared to
radial trajectories [112, 113, 111]. However, respiratory-resolved lung structural imaging
with bSSFP at 0.55T has not yet been explored.
One technical challenge is the design of a respiratory navigator bSSFP. Javed et al.
[107] utilized a superior-inferior (S-I) navigator to estimate respiratory motion every 200ms
at 0.55T using UTE stack-of-spirals imaging. However, inserting a Cartesian based S-I
navigator is likely to disrupt the steady-state of bSSFP, leading to unwanted image artifacts
[149]. FID navigators, which sample the k-space center every TR, capture respiratory-motion
without interrupting the bSSFP steady-state and have been combined with various imaging
trajectories, including 3D cones [118, 150] and WASP [20]. The high sampling rate of the
FID signal is especially advantageous for capturing fast breathing motions, such as those
common in pediatric subjects (e.g., infants, with a resting breathing frequency around 1 Hz
[151]. Therefore, an FID-based respiratory navigator is well-suited for respiratory-resolved
lung structural imaging with bSSFP.
Pulmonary functional analysis is a valuable endpoint for motion-resolved lung MRI. Lung
ventilation can be measured based on either parenchyma deformation [120] or tissue signal intensity changes [132, 135]. PREFUL measures regional lung ventilation from lung
parenchyma signals that are modulated by tissue density changes throughout the respiratory
cycle. This technique generates a full cycle of respiratory phase-resolved images by sorting
63

a continuously acquired, free-breathing image series. At 0.55T, PREFUL has demonstrated
its utility in detecting pulmonary dysfunction in patients with persistent symptoms after
COVID-19, as well as in both adults [113] and pediatric populations [116]. However, 3D
ventilation measurements using PREFUL at 0.55T have not yet been studied. Ventilation
can also be measured based on parenchyma deformation, which primarily relies on detailed
vessel structures and diaphragm movements during breathing. Fei et al. proposed a motioncompensated, low-rank regularized reconstruction using 3D UTE MRI, which provides both
structural and functional lung imaging at 3T [120]. However, the Jacobian-based ventilation, which relies on structural images with low vessel-parenchyma contrast and limited
vessel details, has not been rigorously evaluated.
In this study, we propose a free-breathing high-resolution 3D bSSFP-based lung imaging
method at 0.55T that simultaneously provides respiratory-resolved structural images and
regional ventilation information. The proposed method uses a 3D stack-of-spiral out-in
bSSFP sequence with an embedded FID navigator, constrained reconstruction, and regional
ventilation mapping using Jacobian-based methods [120]. We hypothesize that the proposed
method will be able to reveal detailed vessel features and provide robust, reproducible 3D
ventilation maps.
The proposed technique involves several innovations. 1). Spiral-based trajectories enable
higher sampling efficiency resulting in shorter scan time, compared to 3D radial center-out
trajectories. 2). The structural lung images provide high vessel-to-parenchyma contrast
with short TR, leveraging the high signal intensities of bSSFP at 0.55T, which enables
accurate registration for ventilation estimates [89]. 3). To our knowledge, this is the first
demonstration of 3D ventilation mapping using spiral-based bSSFP at 0.55T, capturing
regional volume changes throughout the entire respiratory cycle.
64

4.2 Method
4.2.1 Experimental methods
Imaging experiments were performed on a whole-body 0.55T system (prototype MAGNETOM Aera; Siemens Healthineers, Erlangen, Germany) equipped with high-performance
shielded gradients (45 mT/m amplitude, 200 T/m/s slew rate). A 6-channel body array
(anterior) and 6 elements from an 18-channel spine array coil (posterior) were used for phantom and in-vivo experiments. The study was approved by the Institutional Review Board of
the University of Southern California. All subjects provided written informed consent. Six
healthy volunteers (3F/3M, 28-35 years old) were scanned in both supine and prone positions,
following a protocol described in Table 4.1 (see Section 4.2.6 for detailed description).
All subjects were scanned using a free-breathing 2-mm isotropic-resolution SOS out-in
imaging protocol with the following imaging parameters: flip angle = 25°, FOV = 432 mm
×432mm×240 mm, matrix size = 216×216×120, T E1/T E2/TR = 0.64/2.60/3.25ms, readout duration = 1.86 ms, peak gradient amplitude = 23 mT/m, maximum gradient slew
rate = 175 T/m/s, total acquisition time = 4min 49sec, and coronal acquisition. Center
frequency was monitored and readjusted before and after each SOS out-in sequence.
To determine the optimal flip angle to maximize lung parenchyma and vessel signals, one
healthy volunteer (F, 28) was asked to perform eight 47sec end-of-inhalation breath-hold
scans with different flip angles in the supine position. The sequence used for the breathhold scans was the clipped free-breathing SOS out-in sequence with a total number of 14,
400 interleaves and 120 unique spiral arms in the kx-ky plane. A 5min pause was inserted
between breath-hold scans.
2D PREFUL images were acquired in six healthy volunteers as a comparison to the 3D
ventilation maps. Two coronal slices were obtained. The anterior slice was positioned at the
tracheal bifurcation and the posterior slice was placed to capture both the main descending
65

aorta and the spine. The imaging parameters were as follows: coronal plane, flip angle = 30◦
,
resolution = 1.7×1.7×15 mm3
, matrix size = 128×128 (interpolated to 256×256), bandwidth
= 1149 Hz per pixel, TR/TE= 1.17/2.33 ms, parallel imaging acceleration factor=2, 250
frames per scan (1 min 2 sec) and temporal resolution = 248 ms/frame.
4.2.2 Pulse sequences
Figure 4.1 (A) illustrates the proposed free-breathing SOS out-in pipeline. Non-selective
excitation was generated with a 200 µs hard pulse, and the slice-selection direction (kz) is in
the anterior-posterior direction. (Top) For each kz partition, a rapid spiral out-in readout
created two echoes in the kx-ky plane. It took approximately 14, 400 interleaves (46.8 sec)
to fully sample the 3D image. (Bottom) Cartesian sampling was used along the kz direction,
following a “ping-pong” order to avoid rapid switching of kz phase-encoding in balanced
SSFP, which induces eddy currents and resulting in image artifacts [152, 153]. The kx-ky
spiral takeoff angle started at 0◦
, and was increased by a tiny golden angle (GA = 23.6281◦
)
for every 120 TRs (120 kz encoding steps) to minimize the eddy currents [154]. A KaiserBessel windowed ramp preparation [155] using 7 TRs was added to minimize oscillatory
transients in bSSFP imaging. The pulse sequence was implemented in the Pulseq framework
[156].
4.2.3 Respiratory navigator and data binning
An FID navigator was designed and embedded into a 3D SOS out-in sequence for respiratory signal extraction. Specifically, a 20-sample ADC event of duration 0.05ms was inserted
every TR between a non-selective RF pulse and a kz partition-encoding gradient. The FID
navigator module has a duration of 0.07ms including two ADC dead times, 0.01ms each,
before and after the ADC event. The respiratory signals were then passed through angular
filtering [157] in the [0-0.04] rad−1 angular frequency domain, followed by bandpass filtering
with a passband [0.05–0.6] Hz. Coil clustering was performed to weight the signals based
66

Figure 4.1: Free-breathing SOS out-in pipeline. (A) 3D sampling scheme. Spiral outin interleaves are acquired in the kx − ky plane, with kz sampling following a ‘ping-pong’
ordering. The spiral rotation angles increase by a tiny GA after each complete set of kz
steps. (B) Reconstruction and ventilation analysis. Respiratory data are sorted according
to the respiratory-direction dependent manner, with raw k-space data binned accordingly.
A constrained reconstruction is applied to generate respiratory-resolved images. Regional
3D ventilation maps are generated following motion field extraction and image registration.
67

on highly correlated coils, resulting in a 1-D respiratory signal with a length equal to the
total number of TRs. Figure 4.1 (B) demonstrates the data-binning strategy and remaining pipeline. The extracted navigator was used to bin the corresponding k-space data in
a respiratory-direction-dependent manner. Respiratory directions were defined by inhalation (from end-of-exhalation to end-of-inhalation) and exhalation (from end-of-inhalation to
end-of-exhalation) periods, based on the sign of the derivatives of the respiratory signals.
Different motion stages were grouped to have an equal displacement range of the respiratory
signal. The sorted signals from the two extreme displacements were combined, representing end-of-exhalation and end-of-inhalation states. k-space data were binned into different
number of motion stages based on the motion estimates.
4.2.4 Image reconstruction
All reconstructions were performed offline using MATLAB R2021a (MathWorks, Natick,
MA). Raw datasets were transferred and processed using ISMRMRD format [158]. The
spiral out-in trajectory was measured using the method proposed by Zhao et al. [159]
using a spherical ball phantom of 14 cm in diameter. The k-space data were pre-whitened
using the correlation matrix calculated from noise measurements [53]. Coil sensitivity maps
were estimated with ESPIRiT [141] and were assumed to be time-invariant. An estimationsubtraction method [160] was used to mitigate the spiral aliasing artifacts originating from
arm regions outside the FOV.
Motion-resolved images X were reconstructed with an iterative constrained reconstruction using spatial L1-based wavelets (Φ) and temporal finite difference (FD) along the motion
dimension. The cost function is given by:
arg min
X
X
M
m=1
||AXm − bm||2
2 + λs||ΦX||1 + λm||FDm(X)||1 (4.1)
where Xm is vectorized motion-resolved 3D images parametrized by the motion state
68

index m (m = 1,2,. . . ,M). A is the encoding function (including coil sensitivity and nonCartesian Fourier transform), bm represents multi-coil k-space data for each sorted motiondimension. Φ represents a spatial 3D L1-wavelet transform, and F Dm represents a firstorder finite difference along the motion dimension. λs and λm are the spatial and motion
regularization parameters, respectively.
This optimization problem was solved using the alternating direction method of multipliers (ADMM) algorithm [161], as implemented in the Berkeley Advanced Reconstruction
Toolbox (BART) [140]. The spiral out and in acquisitions were reconstructed separately and
were then combined using the square root of the sum-of-squares. To determine the suitable
number of motion states which provides sufficient aSNR and minimize motion blur, one 4min
49sec dataset was reconstructed using 4, 8, 12, 16, 20, and 30 states. Parenchyma aSNR
and sharpness (using MD) were evaluated for each reconstruction. Gridding reconstruction,
followed by coil combination using ESPIRIT coil sensitivity maps [141], was performed to
reconstruct breath-hold artifact-free images.
We analyzed the sensitivity of regularization parameters λs and λm in the constrained
reconstruction using k-space data acquired from one healthy volunteer (31,F). The regularization parameters λm and λs were selected visually based on image quality, ventilation
results, parenchyma aSNR, MD, and the mean of the ventilation values. A range of λs values
from 0, 0.001, 0.001 to 0.01 and λm values from 0, 0.001, 0.001 to 0.01 were tested.
4.2.5 Ventilation analysis
All quantitative measurements were performed in MATLAB (MathWorks, Natick, MA).
Lung parenchyma was segmented using the Total Segmentator [162]. The end-of-exhalation
and end-of-inhalation states of each subject were segmented. A 3D tidal volume change
(in percentage) was then calculated by diving the difference in volumes between the two
extreme stages by the end-of-exhalation stage. PREFUL images were processed using the
LungMRI 2.2.0 research package (Siemens Healthineers, Erlangen, Germany) [132], resulting
69

in a fractional ventilation map (in percentage).
The end-of-exhalation frame, typically considered as the stable state, was selected as the
reference frame for registration. All other motion state images were registered to this reference via Demons nonrigid registration [163] (4 pyramid levels, and 100 iterations). Jacobianbased regional ventilation maps (Rvent) were then calculated from the estimated motion fields
using:
Rvent(%) = Vi − Ve
Ve
× 100% = Vi
Ve
− 1 = | det(∇M + Id)| − 1 (4.2)
where Vi and Ve are the lung volumes at i
th respiratory state and the end-of-exhalation
state, respectively. ∇M, a 3 × 3 matrix for each voxel, represents the gradient (∇) of the
motion fields (M). Id is a 3×3 identity matrix. det demotes the determinant of the matrix.
Rvent 0 indicates volume expansion, while Rvent < 0 indicates volume contraction. Since
the reference frame is the end-of-exhalation, where the volume is often at a minimum, Rvent
is typically positive, representing voxel-wise expansion throughout the respiratory cycle.
To evaluate the gravitational effects on ventilation between supine and prone positions,
gravitational profiles were generated by averaging the ventilation maps along the superiorinferior and left-right axes for all subjects. Anterior-Posterior (A-P) profiles were then
interpolated and normalized to a common width (0-1). The 20%-80% range of the profiles
was used to avoid the outlines on the edges. The slope of a linear fit was calculated to
evaluate the rate of ventilation changes along the A-P direction in both supine and prone
positions.
4.2.6 In-vivo experimental study design
Table 4.1 illustrates the experimental design. The total scan included three sections (A,
B and C), lasting 1hr 20min. Sections A vs. B were used to test intra-scan test-retest
repeatability for 3D SOS out-in and PREFUL. Sections A or B versus C were used to
70

evaluate posture-related gravity dependence of ventilation. Sections A and B were scanned
in the supine position, with the volunteer remaining on the patient table, with a 10min
break between the two sections. Section C was scanned in the prone position with the arms
aside the torso. Each scan section lasted approximately 9 minutes and included both a 3D
SOS out-in and 2D PREFUL with each coronal slice measured twice. The order of the 3D
SOS out-in and 2D PREFUL acquisitions was randomly alternated within each section. A
15-min break was taken before each repositioning, which includes the time for changing of
the position (approximately 3min) and allowing the volunteer to stabilize.
Table 4.1: In-vivo experimental design. The scan consisted of 3 sections (A, B and C).
Sections A vs. B were scanned in the supine position, with a 10-minute break between the
two sections to prevent gradient heating. Each section included both a 3D SOS out-in scan
and a 2D PREFUL scan, with two coronal slices, each repeat twice. Section C was scanned
in the prone position. After positioning, there was a 15-minute break to allow volunteers to
stabilize. The center frequency was measured and adjusted before and after the SOS out-in
acquisition to correct for any frequency shifts caused by gradient heating.
71

4.2.7 Image quality analysis
Apparent SNR (aSNR) was calculated for the structural echo-combined images. Regions
of interests (ROIs) of the lung parenchyma and background were manually selected. aSNR
was computed as the mean signal intensity in each ROI divided by the standard deviation
of the signal intensities in the background ROI.
The maximum derivative (MD) of the diaphragm was also used to evaluate the boundary sharpness of the reconstructed images. Two rectangular ROIs were first selected: one
covering the diaphragm boundary in the sagittal view and the other covering the liver area.
MD was then calculated by dividing the maximum diaphragm gradient by the mean signal
intensity in the liver ROI [120]. Similar ROIs were manually selected for all subjects.
4.2.8 Statistical analysis
Bland–Altman analysis of the mean of the ventilation maps was performed between the
anterior and posterior slices of PREFUL within 2 minutes and within 15 minutes, as well
as between repeated SOS out-in scans within 15 minutes. Correlation plots of regional
ventilation were generated for comparisons between 2D PREFUL and 3D SOS out-in, and
between 3D SOS out-in and 3D volume changes (%). Relative differences were calculated
by dividing the difference between two measurements by their mean, and these were used
to assess repeatability, with a threshold of less than 10% defining good repeatability. A
significant level of p < 0.001 was considered for correlation analyses.
4.3 Results
4.3.1 Sequence and reconstruction parameter optimizations
Figure 4.2 illustrates the flip angle optimization for blood and parenchyma signal intensities.
The mean signal intensity curves for (A) blood and (B) parenchyma regions were plotted
72

Figure 4.2: Flip angle optimization. (A) Averaged signal curves for blood and parenchyma
(mean signal intensity ± standard error; solid lines) versus flip angle. Simulated curves are
plotted as dashed lines. (B) Regions of interest (ROIs) from blood and parenchyma areas
(green boxes) . Eight flip angles were tested, ranging from 5◦
to 75◦
in 10◦
increments. A flip
angle of approximately 45◦ provided the strongest blood signal, while a FA in 25◦ yielded
the highest mean parenchyma signal. The parenchyma signals deviating from the simulation
curve (at 35◦ and 45◦
) may be due to inconsistent breath holds.
73

across flip angle ranges from 5◦
to 75◦ based on 8 breath-hold scans. For blood signals, the
measured signal closely matches the simulated curve, with peak signal intensity observed with
FA=45◦
. For parenchyma signals, the trends generally align with the simulation, deviations
were noted at FA=35◦ and 45◦
. The highest parenchyma signal was achieved at FA = 25◦
,
which was selected for subsequent ventilation studies.
To assess the impact of regularization parameters, 8 respiratory phases (with undersampling factor of 1.6) were reconstructed for all combinations of λm and λs, as shown in Figure
4.3. When both λm and λs were set to trivial values (λm = 0 and λs = 0;upper left), noise
dominated the reconstructed images, leading to an underestimation of ventilation values.
In this case, the MD was the highest, while the aSNR of the parenchyma was the lowest.
As λm increases, aSNR improved and the ventilation patterns became more homogeneous.
However, when λm reached 0.01, diaphragm boundary became noticeably blurred (green arrow), resulting in overestimated regional ventilation values (red arrow), although the mean
of ventilation values remained relatively small. Similarly, extreme values of λs caused spatial blurring in the vessels (green box), leading to blurry ventilation maps (e.g., λs = 0.01).
When both λm and λs were within the range of 0.0001 to 0.001, a good balance of denoising
and image sharpness was achieved. Considering the minimal motion blurring (MD>0.350),
realistic ventilation values (3D Rvent = 11.29%), and a decent parenchyma aSNR (>4), a
combination of λs = 0.0001 and λm = 0.001 (red dashed box) were chosen for the remainder
of the analysis.
Figure 4.4 illustrates the process used to determine the suitable number of reconstructed
motion states. As the number of motion states increased from 4 to 30, diaphragm MD
increased, plateauing at around 8 motion states, while the aSNR decreased and plateaued
at around 16 motion states. Ventilation values are closest to the tidal volume reference (red
dashed line in (C)) when 8 motion states were used. Based on these findings, eight motion
states were selected for all subsequent processes.
74

Figure 4.3: Structure and ventilation results for different choices of regularization parameters.
From top to bottom, the spatial L1-wavelet regularization parameter λs varies from 0, 0.0001,
0.001 to 0.01. From left to right, the finite difference regularization term λm along the motion
dimension ranges from 0, 0.0001, 0.001 to 0.01. For each combination of λs and λm, a coronal
view (left) and the corresponding ventilation map (right) are shown. Diaphragm MD and
parenchyma aSNR are reported at bottom, along with the mean values of the 2D ventilation
maps. With trivial values of λm and λs (λm = 0, and λs = 0), the structural images are
dominated by noise, leading to underestimated ventilation. Increasing λm or λs improves the
parenchyma images by reducing noise. However, extreme values of λm (e.g, λm = 0.01) cause
diaphragm blurring (green arrow) and result in overestimated regional ventilation values (red
arrow), while the mean ventilation value becomes underestimated. Similarly, extreme values
of λs cause spatial blurring in the vessels (green box) and result in blurry ventilation maps
(e.g., λs = 0.01). Overall, structural images with λm = 0.001-0.01 and λs = 0.001-0.01 yield
the best image qualities with parenchyma aSNR > 3.00 and diaphragm MD > 0.28. The
suitable parameter set (λm = 0.01 and λs = 0.001; pink dashed box), which provides the
most accurate ventilation estimates compared to tidal volume, was chosen for all subsequent
processes.
75

Figure 4.4: Impacts of motion-states on structure and ventilation. (A) Reconstructed images (end-of-exhalation frame) and ventilation maps in coronal views. The ventilation values
gradually increase as the number of motion states increases from 4 to 16. (B) Image quality matrices plot showing parenchyma aSNR and diaphragm MD. Both the aSNR of the
end-of-exhalation and end-of-inhalation frames decrease as the number of respiratory states
increases, and plateau when number of states exceeds 16. Diaphragm MD increases as the
number of motion states ranges from 4 to 8; and fluctuates when number of motion states
exceeds 10. (C) Violin plots (zoomed-in) of 3D ventilation values for different numbers of
motion states. Red dashed lines indicate the reference of tidal volume change. The overall
ventilation distribution increases as the number of motion states increases from 4 to 16, with
the median of the 8-motion-state configuration being closest to the reference value. Based
on these findings, 8 motion states were selected for use in subsequent analyses.
76

Figure 4.5: Example structural images and ventilation maps from one healthy volunteer
(31/M). Coronal (1st row), axial (2nd row) and sagittal (3rd row) views of images reconstructed from the first echo, second echo and combined echoes in (A) exhalation state and
(B) inhalation states. (C) representative ventilation maps in the corresponding coronal,
axial and sagittal views. Diaphragm positions are illustrated using red dashed lines. Diaphragm boundaries are slightly blurry in inhalation stage. In both end-of-exhalation and
end-of-inhalation stages, echo combined images show superior image qualities with less noise
in the parenchyma, and clear vessel delineation.
4.3.2 Structure and ventilation maps
Overall, the structural images exhibit clear diaphragm boundaries, and detailed vessel structures (Figure 4.5). Echo-combined images demonstrate superior image qualities, with
higher parenchyma aSNR and sharp boundaries compared to single echo images. The second
echo images are nearly identical to the first echo images in terms of the contrast, without
noticeable artifacts. End-of-inhalation state images typically show lower SNR, as noticed
in diaphragm regions, resulting in less distinct diaphragm boundaries compared to the endof-exhalation images. The ventilation maps reveal overall homogeneous ventilation patterns
across all three views, with higher ventilation values observed in the posterior regions of the
sagittal views.
77

4.3.3 Repeatability and validation
Figure 4.6 represents the results of the repeatability studies for regional ventilation in 3D
SOS out-in and 2D PREFUL. (A) Bland-Altman analysis of the mean regional ventilation
values is shown for PREFUL scans within (a) 2-minute, (b) 15-minute, and (c) two 15-minute
separated SOS out-in scans. The mean bias between PREFUL scans within 2 minutes was
0.3% (-1.69, 2.30) in the anterior slice, and 0.74% (-1.54, 3.02) in the posterior slice. The bias
between PREFUL scans within 15 minutes was 0.08% (-1.48, 1.63) for the anterior slice and
0.21% (-2.92, 2.5) for the posterior slice. For SOS out-in, the bias within 15-minute repeated
measurements was 0.19% (-0.88, 1.26) for the anterior slice and 0.78% (1.63, 3.18) for the
posterior slice. Overall, posterior slices showed slightly higher bias compared to the anterior
slices for both PREFUL and SOS out-in. The relative differences in all repeated experiments
(a, b and c) were less than 10%, indicating good repeatability. (B) Representative SOS outin ventilation maps for single-slice measurements are shown in coronal, sagittal and axial
views for two healthy volunteers (30/M, 31/M). The ventilation patterns exhibit strong
agreements between the two scans. Orange boxes highlight similar ventilation patterns in
the coronal and axial views in Volunteer 1 and 2, respectively.
Figure 4.7 compares the ventilation results between SOS out-in and PREFUL. Ventilation quantification using SOS out-in was significantly correlated with global tidal volume
change (R2 = 0.8738, p < 0.0001) and PREFUL measurements (R2 = 0.8967, p < 0.0001) at
matched slice locations. Regional ventilation patterns generally match well between PREFUL and SOS out-in. In particular, posterior slices show consistently high ventilation values.
Orange boxes highlighted regions of higher regional ventilation in same coronal areas in Volunteer 5, as seen in both PREFUL and SOS out-in. However, unmatched regions, such as
the corners of the diaphragm in Volunteer 1 and 2, may result from registration errors.
78

Figure 4.6: Repeatability studies of ventilation. (A) Bland–Altman analysis of the mean
regional ventilation for PREFUL in (a) 2 minutes and (b) 15minutes, and (c) SOS outin in 15 minutes. Anterior and Posterior slices are represented by different colors. Solid
lines represent the mean bias and dashed lines indicate the 95% confidence interval. (B)
Representative SOS out-in ventilation maps from 15-minute repeated scans in three views
for two healthy volunteers (30/M and 31/M). The mean ventilation values for the 2D slices
are reported.
79

Figure 4.7: Ventilation comparison with PREFUL. (A) Correlation plots showing (left) the
correlation between the mean ventilation values for 3D SOS out-in and 3D volume changes,
and (right) the correlation between the mean ventilation values for 3D SOS out-in and 2D
PREFUL. (B) Representative slice-matched regional ventilation maps (anterior and posterior
slices) for PREFUL (1st row) and SOS out-in (2nd row) from all six volunteers. The means
of the ventilation values are reported. Regional ventilation patterns are well-matched in the
indicated areas. A notable increase in ventilation values is observed at the corner of the
diaphragm in the SOS out-in results (red arrow), probably due to the mis-registration.
80

Figure 4.8: Posture-related ventilation differences. Averaged gravitational profiles of mean
ventilation values along the A-P directions for all healthy volunteers in (A) supine and (B)
prone positions, respectively. The values in the plots are represented as means ± standard
deviations. Dashed lines represent the linear fit of the gravitational profiles. In the supine
position (A), there is a clear increase in ventilation values (slope = 5.76), showing an average
ventilation increase of 5.76% across the whole lung. In contrast, the prone position (B)
shows a decreasing trend in ventilation (slope = -11.56), with an average of 11.56% decrease
in ventilation.
4.3.4 Gravity effects
Figure 4.8 illustrates the gravity dependence of ventilation along the A-P direction for all
healthy volunteers in both supine and prone positions. In the supine position (A), a clear
increasing trend in ventilation is observed from anterior to posterior, with an estimated
ventilation difference of 5.76%. In contrast, in the prone position (B), a decreasing trend is
noted, with an approximate ventilation difference of 11.56% across the lung. These findings
highlight how regional ventilation is influenced by gravity, where the higher ventilation values
seen in the dependent regions of the lungs.
81

4.4 Discussion
In this study, we demonstrate a 5-minute free-breathing lung imaging method that provides
2-mm isotropic structural images and 3D regional ventilation maps. The structural lung
images exhibit good vessel-parenchyma contrasts, enabled by the short TR (3.25ms) and
bSSFP readout [149] at 0.55T, which supports ventilation estimation. The ventilation results
show good intra-scan repeatability (relative differences<10%) within a 15-minute time gap
and are strongly correlated with 2D PREFUL and global tidal volume (R2 >0.7, p<0.001).
The proposed method also allows capturing the breathing dynamics of the entire respiratory
cycle, allowing generating more intermediate motion states (Figure 4.4). To our knowledge,
this is the first demonstration of 3D regional ventilation mapping using spiral-based bSSFP,
showing both feasibility and repeatability. This method warrants further investigation in
a larger cohort, including in patients with pulmonary dysfunction and a greater range of
regional ventilation defects, such as chronic obstructive pulmonary disease (COPD) and
asthma.
Several practical considerations were addressed in the in-vivo experimental design. Volunteers were instructed to rest for 15 minutes after positioning and before the starting of
the scan to stabilize lung volume [164] and minimize the effects of lung water redistribution
due to positioning [165]. A 10min break between section A and B was included primarily to
allow the system to cool down, preventing center frequency shifts caused by gradient heating
[166], which is common in spiral acquisitions. We observed a frequency drift of approximately
60Hz during the 5-minute scan, which could lead to banding artifacts in the shoulder/neck,
and inferior liver regions, however not affecting the parenchyma regions. In the repeatability studies, we also found that the changes in volunteers’ breathing patterns, lung volumes
affect the ventilation measurements in both PREFUL and SOS out-in [132]. Therefore, volunteers were instructed to maintain a resting state, without being disturbed or interrupted
by the scan performer, throughout the scan. This ensures that their physical condition was
82

consistent during the whole experiment. This will be crucial for future patient scans. Our
results reveal a gravitational gradient in regional ventilation, with increased ventilation in
the dependent lung regions. This gravity dependence in ventilation has been observed using
PREFUL as well as other imaging modalities such as single-photon emission computed tomography (SPECT), and positron emission tomography (PET) [133, 167, 168]. All prior
work has reported physiologic vertical gradients with increased ventilation with decreasing
lung height. Our results show a larger gradient in the prone position (11.56% difference in
ventilation across the lung) compared to the supine position (5.76% difference), which is
consistent with the studies using PREFUL [133]. Similar results have been reported in PET
studies, which found an averaged of 8% ventilation heterogeneity [168].
This work has limitations. First, the ventilation results rely on a Jacobian-based method,
which uses the calculated motion fields during non-rigid registration. However, intensitybased ventilation measurements were not considered in this study, as the signal differences
between end-of-inhalation and end-of-exhalation were found to be less sensitive. Further
investigation into the reconstruction is needed to account for signal differences between motion states. Additionally, this study provides global repeatability and accuracy comparisons
between SOS out-in and PREFUL in a slice-matched manner, there are limitations in the
registration process. Specifically, PREFUL is limited in 2D registration and doesn’t account
for 3D through-plane motion, which may lead to subtle regional differences in ventilation
between 2D PREFUL and 3D SOS out-in, especially in regions near distal vessels and the
heart. Future studies should include a complete comparison in 3D ventilation, potentially
using hyperpolarized gas [169] MRI, to further validate the 3D ventilation maps.
4.5 Conclusion
We present a free-breathing 3D pulmonary structural and ventilation method at 0.55T that
combines a 3D stack-of-spiral out-in bSSFP acquisition with an embedded respiratory nav83

igator, a constrained reconstruction, and a non-rigid registration. This method provides
motion-resolved lung images with 2-mm isotropic resolution and generates 3D regional ventilation maps in just 5 minutes. The proposed ventilation analysis demonstrates good repeatability, and strongly correlations in both 2D PREFUL in matched slices and 3D tidal
volumes during free-breathing. Additionally, we observe a gravitational gradient in ventilation in the dependent lung regions, consistent with the literature. Overall, 3D SOS out-in
ventilation offers a feasible, repeatable and sensitive method of measuring regional ventilation
in healthy volunteers.
84

Chapter 5
Multidimensional RF Pulse Design with
Consideration of Concomitant Field Effects
5.1 Introduction
Multidimensional radiofrequency (RF) pulses have many common uses including localization
[170], motion tracking[171], reduced field of view (FOV) [172, 173], inner volume imaging
[174], parallel excitation using multiple transmit coils [175, 176], and improving flip angle
homogeneity [177]. The small-tip-angle (STA) approximation to the Bloch equations by
Pauly et al. [170] has facilitated the design of multidimensional RF pulses using an excitation
k-space formalism. This provides an elegant linear approximation that links time-varying
gradients and an RF waveform to the transverse excitation pattern. When large flip angles
9 are desired, certain symmetric conditions should be obeyed for this formalism.
There are several important considerations in multidimensional RF pulse design, each
with practical solutions. Multidimensional RF pulses [174] often have long pulse durations
that can lead to artifacts stemming from T2 decay, and require high peak RF power that
leads to high specific absorption rate (SAR) [178]. Parallel RF transmission is a powerful
85

tool for addressing both of these [175, 176, 179, 180]. Variable rate selective excitation
(VERSE) design [181] and its variants [182, 183, 184, 185, 186, 187, 188] also help to limit
peak RF power and reduce SAR. Another consideration is the accuracy of an excitation
k-space trajectory. Gradient nonidealities [189, 190, 191, 192] cause actual excitation kspace trajectories to deviate from what is intended. Linear and time-invariant gradient
nonidealities can be addressed by considering a gradient impulse response function (GIRF)
[193, 194, 195, 196] during RF pulse design [188, 197, 198].
In this chapter, we describe an extension of the excitation k-space formalism that fully
incorporates concomitant field effects. The concomitant fields, nonrotating transverse magnetic fields in the laboratory frame, were approximated as a Bloch-Siegert (BS) shift [199,
200] in the rotating frame to reduce the computational complexity of Bloch simulation. To
the best of our knowledge, the BS shift approximation is the only theoretical way to derive
additional phase terms due to concomitant fields. The BS shift is an NMR phenomenon that
describes a frequency shift of the spin when an off-resonance RF field is applied, and it is
often used for B
+
1 mapping [201]. After this BS shift approximation, concomitant fields were
modeled as spatiotemporal phase terms consisting of gradient fields and spatial variables.
This expression is identical to that derived for MR acquisition in the presence of concomitant
fields [61], demonstrating yet another symmetry between excitation k-space under the STA
approximation and acquisition k-space [202].
5.2 Theory
5.2.1 Concomitant fields
According to Maxwell’s equations, a vector magnetic field must satisfy two conditions: negligible curl, ∇ × B ≈ 0, assuming that the true and displacement current densities are
negligible [203], and zero divergence, ∇ · B = 0. As a result, additional transverse magnetic
fields, denoted as concomitant fields Bc(r, t), are created when gradient fields are applied.
86

Based on Bernstein et al. [61], a vector Taylor expansion of a vector magnetic field Blab(r, t)
in the laboratory frame neglecting gradient nonlinearity is given by:
Blab(r, t) =
0 0 B0

+


−αGz(t) g Gx(t)
g (α − 1)Gz(t) Gy(t)
Gx(t) Gy(t) Gz(t)




x
y
z


=


−αGz(t)x + gy + Gx(t)z
gx + (α − 1)Gz(t)y + Gy(t)z
0


+


0
0
B0 + G(t) · r


(5.1)
where G(t) = [Gx(t), Gy(t), Gz(t)]T
is a vector of applied gradient fields along the physical
x-, y-, and z-axes, respectively, α is the symmetry parameter which depends on coil geometry.
g is a constant that represents two identical unspecified field derivatives (i.e., g =
∂By
∂x =
∂Bx
∂y )
that are not directly related to the desired gradients, therefore g is often assumed to be zero.
r =


x
y
z


∈ R
3 denotes spatial coordinates in the physical coordinate system, and B0 is the
main magnetic field strength.
87

Figure 5.1: Derivation of concomitant gradient terms from concomitant fields using the
Bloch-Siegert shift. (A) The vector magnetic field Blab(r, t) in the laboratory frame consists of a rotating on-resonance RF field B1(r, t) nonrotating concomitant fields Bc(r, t), a
static B0 field, and gradient fields G(t) · r. (B, C) Bloch-Siegert shift approximation. The
concomitant fields are considered as an off-resonance RF field with a carrier frequency of
ωrf = 0ko. To calculate the BS shift caused by concomitant fields, the vector magnetic field
Blab(r, t) excluding other RF fields is transformed into a rotating frame with an angular
frequency of ωrot = ωrf = 0ko. The difference between the magnitude of the effective magnetic field Beff(r, t)|ωrot=0 and the magnitude of its z component is known as the BS shift
caused by the concomitant fields (indicated by the red curly bracket). (D) Once the BS shift
approximation is used, concomitant fields are treated as an additional off-resonance term in
a rotating frame with a frequency of ωrot = −ω0ko.
Consider a general scenario where an RF field with an angular frequency of ωrf = −ωrfko
and gradient fields are applied to a spin system in the laboratory frame, as shown in Figure
5.1A. Assuming uniform transmit sensitivities, the vector magnetic field Blab(r, t) in the
laboratory frame can be described as a superposition of an RF field, the static B0 field,
gradient fields, and concomitant fields:
88

Blab(r, t) = Rz(−ωrft)B(1,xy)(t) + Bc(r, t) +


0
0
B0 + G(t) · r


with B(1,xy)(t) =


B(1,x)(t)
B(1,y)(t)
0


and Rz(α) =


cos α − sin α 0
sin α cos α 0
0 0 1


(5.2)
where B1,xy(t) is the RF field and Rz(α) represents a right-handed rotation matrix that
rotates a vector about the z axis by an angle α. It is important to note that concomitant
fields Bc(r, t) are stationary (i.e., nonrotating) in the laboratory frame and are considered
as an off-resonance RF field in the context of RF pulse design.
5.2.2 Approximation of concomitant fields as a Bloch-Siegert shift
The BS shift [199, 200] describes a frequency shift of a spin when an off-resonance RF field
is applied. When an on-resonance RF field and concomitant fields are given (Figure 5.1A),
concomitant fields are considered as an off-resonance RF field. To calculate the BS shift
caused by an off-resonance RF field, a vector magnetic field excluding any other RF fields
in the laboratory frame is first transformed into a rotating frame with an angular frequency
equal to the carrier frequency of the off-resonance RF field. In general, the effective field
Beff(r, t)|ωrot describes the laboratory magnetic field Blab(r, t) in a rotating frame with an
angular frequency of ωrot = −ωrotko 47,48:
Beff(r, t)|ωrot = Rz(ωrott)Blab(r, t) + ωrot
γ
. (5.3)
Since the carrier frequency of concomitant fields is zero (i.e., ωc = 0), the effective magnetic
field Beff(r, t)|ωrot=0 in a rotating frame with an angular frequency of ωrot = ωc = 0 ko
89

including static off-resonance ∆f(r) is given by
Beff(r, t)|ωrot=0 =


Bc,x(r, t)
Bc,y(r, t)
B0 + G(t) · r +
2π∆f(r)
γ


(5.4)
This is graphically described in Figure 5.1B. The magnitude of Beff(r, t)|ωrot=0 is slightly
increased relative to the z component of Beff(r, t)|ωrot=0 and this additional contribution is
approximated as the BS shift caused by the off-resonance RF field (Figure 5.1C). Following
the approach by Bernstein et al. [61], the magnitude of Beff(r, t)|ωrot=0 is given by
||Beff(r, t)||l2 =
s
(Bc,x(r, t))2 + (Bc,y(r, t))2 +

B0 + G(t) · r +
2π∆f(r)
γ
2
= B0
s
1 +
(Bc,x(r, t))2 + (Bc,y(r, t))2
B
2
0
+ 2
(G(t) · r + · · ·)
B0
+

(G(t) · r + · · ·)
B0
2
(5.5)
Using the Taylor series expansion √
1 + u = 1 + 1
2
u −
1
8
u
2 +
1
16u
3 − · · · and substituting
Bc(r, t) from Equation 1 into Equation 5, the magnitude consists of the static B0 field,
gradient fields, static off-resonance, and a spatiotemporal BS shift (known as concomitant
gradient terms [61]):
||Beff(r, t)|ωrot=0||l2 = B0 + G(t) · r +
2π∆f(r)
γ
+
X
Nl
l=4
hl(t)pl(r), (5.6)
where l counts the concomitant gradient terms, hl
is the l
th dynamic coefficient, and pl
is
the l
th concomitant field basis function. See Table 1 of reference [64] for full expressions of
concomitant gradient terms. Assuming a symmetric gradient system (α=0.5) and neglecting
g, the concomitant gradient terms (up to order 1/B0) are expressed as
N
Xl(=9)
l=4
hl(t)pl(r) = 1
2B0



Gx(t)
2
z
2 + Gy(t)
2
z
2 +
1
4Gz(t)
2
(x
2 + y
2
)
−Gx(t)Gz(t)xz − Gy(t)Gz(t)yz



(5.7)
90

Concomitant gradient terms of order (1/B2
0
) are ignored at 0.55T because of their negligible
effects [64]. Note that in RF pulse design both concomitant fields and static off-resonance
should be compensated simultaneously to mitigate off-resonance effects on the excitation
profile.
5.2.3 Solving Bloch equations with concomitant fields in a rotating frame
Approximating concomitant fields as the BS shift and with the expression from Equation
2 and Equation 3, the effective field ||Beff(r, t)||ωrot=−ω0 of Blab(r, t) in an on-resonance
rotating frame with an angular frequency of ωrot = −ω0ko = −ωrfko = −γB0ko including
static off-resonance can be written as
Beff(r, t)|ωrot=−ω0 ≈ Rz(ω0t)Rz(−ωrft)B1,xy(t)+
Rz(ω0t)


0
0
B0 + G(t) · r +
2π∆f(r)
γ +
PNl
l=4 hl(t)pl(r)


+


0
0
−ω0
γ


≈


B1,x(t)
B1,y(t)
G(t) · r +
2π∆f(r)
γ +
PNl
l=4 hl(t)pl(r)


(5.8)
This is graphically described in Figure 5.1D. The Bloch equations in this rotating frame
neglecting relaxation can be written as
∂Mrot(r, t)
∂t |ωrot=−ω0 = Mrot(r, t)|ωrot=−ω0 × γBeff(r, t)|ωrot=−ω0
. (5.9)
91

For notational simplicity, we omit subscripts “rot”, “eff”, and “ωrot = −ω0”. Equation 9
can be recast into a matrix-vector form as:


dMx(r,t)
dt
dMy(r,t)
dt
dMz(r,t)
dt


= γ


0 Bz(r, t) −By(r, t)
−Bz(r, t) 0 Bx(r, t)
By(r, t) −Bx(r, t) 0




Mx(r, t)
My(r, t)
Mz(r, t)


. (5.10)
The STA approximation (Mz(r, t) ≈ M0,
dMz(r,t)
dt ≈ 0) yields
dMxy(r, t)
dt + jγBz(r, t)Mxy(r, t) = jγM0Bxy(r, t) (5.11)
where Bxy(r, t) = Bx(r, t) + jBy(r, t), and Mxy(r, t) = Mx(r, t) + jMy(r, t). Note that the
Bz(r, t) term does not include additional concomitant gradient terms in the standard STA
approximation.
The above differential equation can be solved by multiplying both sides of Equation 11
by the integrating factor e
j
R t
−∞ γBz(r,s) ds and with the product rule of partial derivatives. The
detailed derivation can be found in reference [204]. The final expression can be expressed as
Mxy(r, t) = jγM0(r, t)
Z t
−∞
Bxy(r, τ ) exp(k(t) · r)
× exp (j2π∆f(r)(τ − t)) exp
−jγX
Nl
l=4
Z t
τ
hl(s) ds
pl(r)
!
dτ.
(5.12)
where the excitation k-space k(t) is defined as k(t) = −γ
R t
τ G(s) ds. Here, the concomitant
fields are considered as an additional off-resonance term during the STA approximation.
Note that the expression in Equation 7 is identical to the expression in Nielsen’s work [62]
when α = 0.5. While our derivation considers single-channel RF pulse design for simplicity,
Nielsen’s work [62] includes multi-channel RF pulse design.
92

5.3 Methods
5.3.1 Validation of Bloch-Siegert Approximation
We evaluated the estimation accuracy across simulated field strengths (from 0.005T to 7T)
using an 8-ms 2D excitation RF pulse [170] with a single-shot spiral-in k-space trajectory.
The dwell time was set to 0.7ns (1.43 GHz sampling frequency). Concomitant fields are
stationary in the laboratory frame but are rotating in a counterclockwise direction in the
Larmor frequency rotating frame. The dwell time must be chosen sufficiently small to resolve
the highest simulated Larmor frequency in the rotating frame, which for this study was 298.06
MHz at 7 Tesla. The target pattern was placed at 20cm from isocenter along the physical
z axis (B0 direction). The ‘laboratory frame’ Bloch simulation refers to Bloch simulations
where concomitant fields were realized as an off-resonance RF field, i.e., Rz(ω0t)Bc(r, t).
The ‘rotating frame’ Bloch simulation refers to Bloch simulations where concomitant fields
were approximated as additional off-resonance terms using Equation 8. Both magnitude
and phase of excitation patterns were reported. Normalized root-mean-square excitation
errors (NRMSEs) = ||mlab−mrot||l2
||mlab||l2
were calculated where mlab and mrot denote the excitation
profiles by the ‘laboratory frame’ and ‘rotating frame’ Bloch simulations, respectively.
5.3.2 Evaluation of Bloch-Siegert Approximation - Multi-channel
RF design
To evaluate the approximation of the concomitant fields using the BS shift, we designed
spatially selective multi-channel RF pulses with iterative VERSE-guided RF pulse design
(i.e., reVERSE) proposed by C¸ avu¸so˘glu et al. [197] This method integrates GIRF prediction of k-space trajectories at each iteration during RF pulse design. We compared the
performance of RF pulses designed with reVERSE methods using the original method and
the proposed method. The original method refers to reVERSE using the system matrix
93

without concomitant gradient terms (i.e., reVERSE w/o. concomitant fields). The proposed
method refers to the reVERSE method using the system matrix with concomitant gradient
terms (i.e., reVERSE w/ concomitant fields) [197]. Excitation profiles by both methods were
simulated using the ‘rotating frame’ Bloch simulation that approximates concomitant fields
as off-resonance terms using the BS shift. A target pattern was a 52 × 52mm2
square beam
and was placed at a 6-mm distance from isocenter along the physical y-axis. A single-shot
spiral-in trajectory was used with 26ms duration and 0.25cm resolution, supporting an FOV
of 128 × 128mm2
. The RF and gradient raster times were set to 6.4 µs. The maximum
peak RF amplitude was set to 12 µT. The maximum gradient amplitude and maximum slew
rate were set to 30 mT/m and 180 T/m/s, respectively. All computations were performed
in MATLAB (MathWorks, Natick, MA, USA).
We utilized published GIRF measurements (Setup A by C¸ avu¸so˘glu et al. [197]), B
+
1
and ∆B0 maps (obtained at 3T with 8-element parallel transmit by Malik et al. [188]).
The target flip angle was set to 90°. The impact and correction of concomitant fields were
demonstrated at various B0 field strengths (0.2T, 0.55T, 1.5T, 3T, and 7T). The effect of
a distance from isocenter (0cm – 30cm) along the physical z-axis on excitation profiles was
demonstrated in Supporting Figure 5.6. The effects of T2 relaxation as a function of pulse
duration were demonstrated in Supporting Figure 5.7. The real (Mx
M0
) and imaginary (My
M0
)
parts of excitation profiles were reported along with NRMSE relative to the target pattern.
5.3.3 Evaluation using single-channel RF design and phantom experiment at 0.55T
Phantom experiments were performed on a whole-body 0.55T scanner (prototype MAGNETOM Aera; Siemens Healthcare, Erlangen, Germany) with gradients capable of 45mT/m
maximum amplitude and 200T/m/s maximum slew rate 50 and a single-channel body transmit coil. We re-performed the multi-channel simulation for this system. The target pattern
and k-space trajectory were the same as those used in the multi-channel simulation. The
94

maximum RF amplitude was set to 13.3 µT. The maximum gradient amplitude and maximum slew rate were set to 24 mT/m and 180.18 T/m/s, respectively. GIRF measurements
were obtained with Duyn’s method [193] as described by reference [205]. The impact and
correction of concomitant fields were demonstrated at various off-isocenter distances (0cm,
5cm, 10cm, and 15cm) along the physical z-axis.
A small uniform ball phantoms with 10cm diameter were chosen to match the supported
FOV in the simulation. A 16-channel head/neck receive coil was used for phantom experiments. The center slice of a ball phantom was imaged at various off-isocenter distances
from isocenter (0cm, 5cm, 10cm, 15cm) along the physical z-axis. Standard B0 shimming
was used. A magnitude B
+
1 map at each distance was acquired using saturation recovery
TurboFLASH [206] with imaging parameters: single-slice, axial, FOV = 220×220mm2
, slice
thickness = 6 mm, flip angle = 80◦
, matrix size = 64 × 64, TE = 1.92 ms, and total TR
= 6040 ms. A single-echo 2D Cartesian FLASH scan was repeated to acquire datasets at
different TEs (3.7, 4.7, 5.7, 6.7, and 7.7ms) with imaging parameters: single-slice, axial,
FOV = 220 × 220mm2
, slice thickness = 6 mm, flip angle = 20°, matrix size = 64 × 64, and
TR = 11ms. A linear least square fitting was performed to estimate ∆B0 in a pixel-by-pixel
manner. The obtained B
+
1
and ∆B0 maps were incorporated in the RF pulse design for both
numerical simulations and actual experiments.
The designed 2D RF pulses were used as excitation in a 2D spin-echo sequence with a
pair of selective 180° adiabatic pulses as a refocusing pulse. Other imaging parameters were:
single-slice, FOV = 256 × 256mm2
, slice thickness = 5 mm, flip angle = 90◦
, matrix size =
256 × 256, and bandwidth = 122 Hz/pixel. We used the open-source Pulseq framework
[156] for pulse sequence implementation.
95

5.4 Results
Figure 5.2 shows the accuracy of the approximation of concomitant fields as concomitant
gradient terms (i.e., off-resonance terms) using the BS shift. Excitation profiles produced by
the ‘laboratory frame’ Bloch simulation (top/bottom: magnitude/phase) are shown in Figure 5.2A. The magnitude (row i, ii) and phase (row iii, iv) of excitation profiles produced
by the ‘laboratory frame’ Bloch simulation (i, iii) and ’rotating frame’ Bloch simulation
(ii, iv) at various field strengths (0.2T, 0.55T, 1.5T, 3T, and 7T) are shown. NRMSEs
are 4.30/1.60/0.64/0.33/0.16% for 0.2T to 7T, respectively. NRMSEs as a function of field
strength are shown in Figure 5.2B. The estimation error is smaller for larger B0. The
approximation yields NRMSEs less than 5% for B0 greater than 0.2T.
96

Figure 5.2: Evaluation of the accuracy of the approximation of concomitant fields as concomitant gradient terms using the BS shift at various field strengths. Note that the ‘laboratory
frame’ Bloch simulation realizes the concomitant fields as a rotating off-resonance RF field
(counter-clockwise) in the Larmor frequency rotating frame without the BS shift approximation, while the ‘rotating frame’ Bloch simulation approximates the concomitant fields as
off-resonance terms using the BS shift in the rotating frame. (A) (Left) Excitation profile
produced by a 2D RF pulse without concomitant fields (top/bottom: magnitude/phase).
(Right) Magnitude (row i, ii) and phase (row iii, iv) of excitation profiles simulated (i,
iii) without and (ii, iv) with the BS shift approximation at 0.2T, 0.55T, 1.5T, 3T and 7T.
NRMSEs are 4.30/1.60/0.64/0.33/0.16% for 0.2T to 7T, respectively. (B) NRMSEs as a
function of field strength. Note that NRMSEs are < 5% for B0 ≥ 0.2T.
Figure 5.3 demonstrates the off-isocenter dependence of concomitant field effects at
0.55T. The (A: i, iii) imaginary and (A: ii, iv) real parts of excitation profiles produced
by (i, ii) the original method and (iii, iv) the proposed method are shown. The excitation
profiles were produced with the ‘rotating frame’ Bloch simulation. The remaining signals
97

on the backgrounds of the excitation profiles are reduced in the proposed method when offisocenter distance are less than 15cm. There are some remaining artifacts on the edge of
the excitation profile at z = 15cm. NRMSEs are plotted, which increases with off-isocenter
distance. The proposed method completely mitigates NRMSE by 62.69% for shifts up to
10cm.
Figure 5.3: Off-isocenter dependence of concomitant field effects at 0.55T without (A, 1-2
rows) and with (A, 3-4 rows) consideration and correction of the concomitant fields. (i, iii)
Imaginary and (ii, iv) real parts of excitation profiles produced by the (i, ii) original and (iii,
iv) proposed methods. The original and proposed methods refer to reVERSE method with
and without consideration of the concomitant fields during RF pulse design, respectively.
Shown are excitation profiles at various distances from isocenter ranging from 0 to 30 cm
with an increment of 5 cm (1st – 7th columns). (B) NRMSEs are plotted on the right. Notice
that NRMSE increases with off-isocenter distance. The proposed method is able to improve
the NRMSE of the excitation patterns by more than 50% when off-isocenter up to z = 10cm.
Figure 5.4 demonstrates the B0 dependence of concomitant fields. The (A: i, iii) imaginary and (A: ii, iv) real parts of excitation profiles were generated using the (i, ii) original
and (iii, iv) proposed methods. The original and proposed methods refer to reVERSE
method with and without consideration of the concomitant fields during RF pulse design,
98

respectively. The excitation profiles were produced with the ‘rotating frame’ Bloch simulation. The target pattern (1st column) and excitation profiles at 0.2T, 0.55T, 1.5T, 3T,
and 7T are shown (2nd – 6th columns). The effects of concomitant fields are large at low
fields (0.2T and 0.55T) and diminish with higher field strengths (1.5T, 3T and 7T) with the
original method (A: row i, ii). The NRMSEs of the proposed method are further reduced by
30.95%, 55.57%, 74.86%, 70.77% and 38.53% at 0.2T, 0.55T, 1.5T, 3T and 7T, respectively.
The reduction of NRMSE is less at 7T is mainly due to the large B0 variations used in
the simulation. Note that the proposed method failed to completely mitigate the excitation
errors at 0.2T, showing noticeable spatial blurring on the edge of the target pattern (orange
arrow). The proposed method corrected the concomitant fields effects when field strength
> 0.2T, reducing NRMSEs by more than 50% in other cases (B0 > 0.2T).
99

Figure 5.4: Simulated B0 dependence of concomitant fields obtained without (A, i-ii) and
with (A, iii-iv) modeling of concomitant fields effects. (i, iii) Imaginary and (ii, iv) real
parts of excitation profiles produced by the (i, ii) original and (iii, iv) proposed methods.
The proposed method refers to the reVERSE method using the system matrix including
concomitant gradient terms. The excitation profiles were simulated with the ‘rotating frame’
Bloch simulation. Shown are the target pattern (1st column) and excitation profiles at 0.2T,
0.55T, 1.5T, 3T and 7T (2nd–6th columns). (B) NRMSEs are plotted as a function of field
strength. The NRMSE decreases as B0 increases, and excitation errors are reduced with the
proposed method. The NRMSEs are further reduced by 30.95%, 55.57%, 74.86%, 70.77%
and 38.53% at 0.2T, 0.55T, 1.5T, 3T and 7T, respectively. Note that the proposed method
fails in mitigating concomitant fields using the BS shift approximation at 0.2T, showing
noticeable spatial blurring on the edge of the target pattern (orange arrow).
Figure 5.5 shows the results of single-channel (A) simulation and (B) phantom validation at 0.55T. Excitation profiles at off-isocenter distances of z = 0cm, 5cm, 10cm, and
15cm (1st – 4th columns) are shown using the (i, iii) original and (ii, iv) proposed methods. The excitation profiles produced by both methods were simulated using the ‘rotating
frame’ Bloch simulation. The simulation study shows an increase in background signals of
the excitation profiles created by the original with a larger distance from isocenter. This is
100

especially noticeable at z = 15cm. In addition, the edges of the excitation profile are significantly blurred. In contrast, the proposed method retains sharp boundaries of the excitation
profiles (pointed by orange arrows) at all off-isocenter distances. The results from phantom
experiments match well with those from numerical simulations, showing the similar level of
blurring near the boundaries of the excitation profiles for off-isocenter cases. Red-dashed
boxes indicate the target pattern. In the original method, the excitation profile matches well
with the target pattern at the isocenter, and becomes dilated and blurred for off-isocenter
cases. However, the proposed method shows excitation profiles matching reasonably well
with the target profile in all cases. A green arrow indicates a region with ringing artifacts on the edge of a ball phantom, which may be caused by inaccuracy in inhomogeneity
measurements, and/or the assumption of the BS shift failed.
101

Figure 5.5: Single-channel (A) simulation and (B) phantom validation results at 0.55T. (i,
ii) Simulated and (iii, iv) experimental excitation profiles at various off-isocenter distances
(1st–4th columns: 0/5/10/15 cm) are produced by the (i, iii) original and (ii, iv) proposed
methods. The excitation profiles were simulated with the ‘rotating frame’ Bloch simulation.
The excitation profiles produced by the original method show noticeable unwanted, background signals at z =10cm and z=15cm, which are mitigated with the proposed method.
These artifacts are reproduced in the real phantom results with the original method. The
proposed method provides sharper profiles (emphasized using orange arrows) compared to
the original method. A green arrow indicates remaining ring artifacts near the boundary of
a ball phantom, which may be caused by inaccuracies in field inhomogeneity measurements.
102

5.5 Discussion
This work has many potential applications ranging from cylindrical excitations used for respiratory navigators to 3D selective excitations for tumor spectroscopy. We chose cylindrical
excitations as an example because they are widely used in cardiac and thoracic imaging,
which appears to be promising at low field strengths [207, 208, 89] where concomitant fields
are substantial.
The proposed approach relies on modeling nonrotating concomitant fields in the laboratory frame as the BS shift in the rotating frame. This BS shift approximation should be
valid whenever the static B0 field is substantially larger in magnitude than the concomitant
fields (in the transverse plane) induced by gradient fields. For the scenario investigated in
this manuscript, the BS shift approximation appears to be valid when the field strength is
larger than 0.2T, producing NMRSEs less than 5% (see Figure 5.2B).
This work has several limitations. First, the findings are primarily based on numerical
simulations and a validation experiment was only performed at 0.55T with a single-channel
body transmit coil. Additional evaluation is needed on different MRI platforms. Second,
∆B0 variations used in the simulation were scaled by the B0 field strength, but the B
+
1
patterns were taken to be consistent across all field strengths for simplicity, in Figure 5.4.
B
+
1
is primarily dependent on RF coil design, and it is limited by SAR, especially at high
field strength. As field strength increases, the field patterns become spatially more complex
[209, 210]. Third, we only considered 2D selective excitations using spiral trajectories because
these typically have a small excitation field-of-view and benefit from sharp excitation profiles.
We provide source code in the spirit of reproducible research so that motivated readers can
explore the benefits of the proposed formalism for other applications including 3D selective
excitations.
103

5.6 Conclusion
We demonstrate a new small-tip multidimensional RF pulse design procedure that incorporates concomitant fields effects. This new procedure provides increased accuracy in excitation
profiles for all scenarios tested and was validated at 0.55T with a single-channel body transmit coil. In many cases, the errors in excitation profiles were reduced by more than 50%.
The impact is greatest in scenarios where concomitant fields effects are substantial, such as
low field strengths, strong gradient fields, and off-isocenter. The proposed method provides
substantially sharper excitation profiles for off-isocenter RF pulse design at 0.55T.
104

Figure 5.6: Representative optimized waveforms obtained with the proposed method (i.e.,
reVERSE method with the system matrix including concomitant gradient terms) at 0.55T
with parallel transmission. (Top) Optimized 8-channel RF waveforms. The resulting peak
RF amplitude is below the target magnitude (12 µT). The RF duration was optimized to be
17.66 ms. (Bottom left) GIRF-predicted gradient waveforms were used and updated at each
iteration of reVERSE. Nominal gradient waveforms used in the last iteration of reVERSE are
provided as an output. (Bottom right) Excitation k-space trajectory obtained with nominal
gradient waveforms.
105

Figure 5.7: Pulse duration dependence of T2 effects in conjunction with concomitant field
effects at 0.55T. (A: 1st - 2nd columns) without and (A: 3rd−4
th columns) with T2 relaxation
using T2 = 40ms during pulse simulation. (i, iii) Imaginary and (ii, iv) real parts of
excitation profiles produced by the (i, ii) original and (iii, iv) proposed methods. Original
and proposed methods refer to the reVERSE method without and with consideration of
concomitant fields, respectively. The excitation profiles were simulated with the ‘rotating
frame’ Bloch simulation. Shown are (1st and 3rd columns) excitation profiles obtained with
a 17.66ms pulse duration and (2nd and 4th columns) excitation profiles obtained with a
52.07ms pulse duration corresponding to a larger unaliased FOV (38.4 × 38.4 mm2
). The
excitation profiles show blurring near the edges, as expected, for longer excitation pulses
and shorter T2. These artifacts are reduced with the proposed method compared to the
original one. Comparing the results obtained without and with T2 relaxation, we observe
that both components of magnetization (Mx
M0
and My
M0
) are reduced in the excitation profiles.
(B) NRMSEs are plotted on the right. Notice that NRMSE increases with pulse duration,
and these errors are substantially reduced using the proposed method. The errors induced
by concomitant fields increase with pulse duration. The proposed method reduced the errors
by more than 50% when RF duration is less than 20ms.
106

Chapter 6
Concluding remarks
6.1 Summary
This dissertation presents several MRI technologies focused on the applications of 3D dynamic MRI for speech and lung imaging, using field strengths of 1.5T and 0.55T, respectively.
It also introduces an advanced excitation method designed to improve imaging of small volumes over large body areas. The following list summarizes my key contributions:
1) The development of 3D RT-MRI of speech production with improved spatiotemporal
resolution for better visualization of fast lip and tongue movements during normal and 1.5-
times speeded speech. It combines several novel techniques: 3D stack-of-spiral trajectory,
efficient k-t space sampling, parallel imaging and constrained reconstruction. It enables
better capturing alternating high, low tongue postures during vowels during natural speech
tasks, and capture complex tongue shapes during fast alveolar consonant segments.
2) The development of high-resolution 3D free-breathing lung ventilation mapping 0.55T
that achieves motion resolved lung images with 2mm isotropic resolution, along with 3D
regional ventilation maps within a 5min scan time. This proposed ventilation analysis shows
good repeatability, and values that are strongly correlated with 2D reference method in
matched 2D slices, and with 3D tidal volumes in healthy volunteers, offering a feasible,
repeatable and sensitive measurement.
107

3) The development of a multidimensional RF pulse design procedure that incorporates
linear time-invariant gradient imperfections and concomitant field effects. This new procedure provides increased accuracy in excitation profiles for all scenarios tested and was
validated at 0.55T with a single-channel body transmit coil. In many cases, the errors in
excitation profiles were reduced by more than 50%. The impact is greatest in scenarios
where concomitant fields effects are substantial, such as low field strengths, strong gradient
fields,and off-isocenter. It provides substantially sharper excitation profiles for off-isocenter
RF pulse design at 0.55T.
6.2 Future directions
Looking ahead, there are promising directions that warrant further investigation for each
work.
Beyond the improved 3D Real-time MRI of speech imaging, data sampling and constrained reconstruction are intertwined. It is likely that there are tailored constrained reconstruction approaches that would work better for each sampling strategy. Possibilities include
regularized nonlinear inversion [211] and partial separability model-based reconstruction.
[82] Reconstruction time may be reduced by using parallelized GPU computation, as individual time segments can be processed independently. The GPU-based reconstruction has
been investigated by multiple groups and has provided 3-200-fold reduction in computation
time for MRI-constrained re-construction. [212, 213] This would be a valuable direction for
future work, especially for interactive RT-MRI applications that require low-latency reconstruction.
The proposed data sampling and reconstruction method may benefit other applications,
particularly those in which dynamic edge information is critical—for example: 1) dynamic
MRI of cardiac function [214, 215, 216] where the endo- and epicardial contours are of
greatest importance, 2) dynamic MRI of the upper and lower airway [217, 218] where the
108

pharyngeal airway and trachea are most critical, and 3) dynamic MRI of joint motion [219,
220] where the movement of bones and cartilage are most critical. The proposed sampling
strategy may also benefit a 3D stack-of-radial (a.k.a. stack-of-stars) approach due to its
similar distribution in 3D k-space [221, 222]. This data sampling scheme has benefits of
mitigating motion-induced temporal blurring, especially with respect to radial trajectory
that has been proved to be robust to motion artifacts [223].
In the development of the 3D lung ventilation at 0.55T, there are great potentials to
include larger cohort for repeatability, reproducibility and accuracy studies, including patients with chronic lung diseases, such as COPD, asthma and cyclist fibrosis. The high scan
efficiency of the stack-of-spiral out-in enables the proposed method to be augmented with
magnetization preparation acquisition to create desired tissue contrasts, such as T2 preparation and/or diffusion preparation, which will be beneficial for lung cancer screening and
nodule characterization. The inherent dual echo further allows additional T
∗
2 mapping, that
may be beneficial for the evaluation of articular cartilage osteoarthritis, and iron overload
early detection, monitoring, and treatment in the heart.
The 3D spiral out-in bSSFP sequence can be combined with variant magnetization preparation pulses, expanding the possibilities for lung imaging at 0.55T. The high vessel-toparenchyma contrast provided by this sequence allows for detailed exploration of pulmonary
perfusion, such as through arterial spin labeling (ASL) with spiral sampling [224]. By
combining SOS out-in ventilation mapping with perfusion assessments, we can achieve a
comprehensive 3D V/Q analysis of the lung at 0.55T, offering enhanced diagnosis capabilities and deeper insights into pulmonary physiology. In addition, the high scan efficiency of
the proposed method the integration of T2 preparation and/or diffusion preparation, create
desired tissue contrasts that are particularly useful for lung cancer screening and nodule
characterization.
In the work of multi-dimensional RF excitation with the consideration of concomitant
field effects, the current work includes preliminary simulations across field strengths, and
109

initial phantom results. Further investigations of the potential applications should be the
next step, ranging from cylindrical excitations used for respiratory navigators to 3D selective
excitations for tumor spectroscopy. The cylindrical excitation used in the simulation are
widely used in cardiac and thoracic imaging, which appears to be promising at low field
strengths [207, 208, 89] where concomitant fields are substantial. It is suitable for design
navigators to keep track of the motion within a limited region-of-interests, and away from
isocenter. It is also well-suited for 3D static imaging with large FOV, so that to increase the
scan efficiency by zooming-in a smaller volume-of-interests, and help diagnosis, the treatment
during radiation therapy in the clinic.
With this thesis, and other ongoing research, MRI has been shown to be a great imaging tool to provide detailed volumetric images, additional functional information, and pronounced spatial-temporal capabilities, for better understanding of human natural speech and
diagnosing of pulmonary dysfunction. Although faced with challenges, I am convinced that
using 3D dynamic MRI has its unique values and has a promising future.
110

Bibliography
[1] Dynamic Imaging Science Center. url: https://sites.usc.edu/disc/ (visited on
12/20/2024).
[2] Ye Tian and Krishna S Nayak. “Real-time water/fat imaging at 0. 55T with spiral
out-in-out-in sampling”. In: Magnetic Resonance in Medicine 91.2 (2024), pp. 649–
659.
[3] Lukasz Mik et al. “Fusing the electromagnetic articulograph, high-speed video cameras and a 16-channel microphone array for speech analysis”. In: Bulletin of the Polish
Academy of Sciences. Technical Sciences 66.3 (2018).
[4] The Saxophonist’s Anatomy. url: https : / / www . stevenjordheim . com / Larynx _
Pharynx_Anatomy.html (visited on 12/20/2024).
[5] Korin Richmond and Steve Renals. “Ultrax: An animated midsagittal vocal tract
display for speech therapy”. In: Proceedings of the 13th Annual Conference on the
International Speech Communication Association (INTERSPEECH). International
Speech Communication Association. 2012, pp. 74–77.
[6] Ian Stavness et al. “Coupled hard–soft tissue simulation with contact and constraints
applied to jaw–tongue–hyoid dynamics”. In: International Journal for Numerical
Methods in Biomedical Engineering 27.3 (2011), pp. 367–390.
[7] Yoon-Chul Kim. “Fast upper airway magnetic resonance imaging for assessment of
speech production and sleep apnea”. In: Precision and Future Medicine 2.4 (2018),
pp. 131–148.
[8] Joseph S Perkell et al. “Electromagnetic midsagittal articulometer systems for transducing speech articulatory movements”. In: The Journal of the Acoustical Society of
America 92.6 (1992), pp. 3078–3096.
[9] Mark I Singer and Eric D Blom. “An endoscopic technique for restoration of voice
after laryngectomy”. In: Annals of Otology, Rhinology & Laryngology 89.6 (1980),
pp. 529–533.
[10] Pascal Perrier, Louis-Jean Bo¨e, and Rudolph Sock. “Vocal tract area function estimation from midsagittal dimensions with CT scans and a vocal tract cast: modeling the
transition with two sets of coefficients”. In: Journal of Speech, Language, and Hearing
Research 35.1 (1992), pp. 53–67.
111

[11] Jeffrey J Fredberg et al. “Airway area by acoustic reflections measured at the mouth”.
In: Journal of Applied Physiology 48.5 (1980), pp. 749–758.
[12] Maureen Stone et al. “Cross-sectional tongue shape during the production of vowels”.
In: The Journal of the Acoustical Society of America 83.4 (1988), pp. 1586–1596.
[13] Erik Bresch et al. “Seeing speech: Capturing vocal tract shaping using real-time magnetic resonance imaging [Exploratory DSP]”. In: IEEE Signal Processing Magazine
25.3 (2008), pp. 123–132.
[14] Didier Demolin et al. “Real-time MRI and articulatory coordination in speech”. In:
Comptes Rendus. Biologies 325.4 (2002), pp. 547–556.
[15] Sandra Rua Ventura, Diamantino Rui Freitas, and Jo˜ao Manuel RS Tavares. “Application of MRI and biomedical engineering in speech production study”. In: Computer
Methods in Biomechanics and Biomedical Engineering 12.6 (2009), pp. 671–681.
[16] Respiratory System. url: https://www.webmd.com/lung/picture-of-the-lungs
(visited on 12/20/2024).
[17] Yuan Lei. “Lung Ventilation: Natural and Mechanical”. In: Medical Ventilator System
Basics: A clinical guide. Oxford University Press, June 2017. isbn: 9780198784975.
[18] Sara M Tony and Mohamed EA Abdelrahim. “Respiratory Tract, Ventilation, and
Pulmonary Functions Tests”. In: Journal of Clinical and Nursing Research 6.3 (2022),
pp. 209–215.
[19] V/Q Scan. url: https : / / radiopaedia . org / articles / vq - scan - 2 (visited on
12/20/2024).
[20] Grzegorz Bauman et al. “Submillimeter lung MRI at 0.55 T using balanced steadystate free precession with half-radial dual-echo readout (bSTAR)”. In: Magnetic Resonance in Medicine 90.5 (2023), pp. 1949–1957.
[21] Adrienne E Campbell-Washburn. “Lung imaging using high-performance low-field
magnetic resonance imaging”. In: American Journal of Respiratory and Critical Care
Medicine 201.11 (2020), pp. 1333–1336.
[22] Robert O Crapo. “Pulmonary-function testing”. In: New England Journal of Medicine
331.1 (1994), pp. 25–30.
[23] Mario Silva et al. “Low-dose CT for lung cancer screening: Position paper from the
Italian college of thoracic radiology”. In: La Radiologia Medica 127.5 (2022), pp. 543–
559.
112

[24] Jing Liu et al. “Neonatal lung diseases: lung ultrasound or chest x-ray”. In: The
Journal of Maternal-Fetal & Neonatal Medicine 34.7 (2021), pp. 1177–1182.
[25] George R Washko, Grace Parraga, and Harvey O Coxson. “Quantitative pulmonary
imaging using computed tomography and magnetic resonance imaging”. In: Respirology 17.3 (2012), pp. 432–444.
[26] Diana L Miglioretti, E Johnson, AE Williams, et al. “Pediatric computed tomography
and associated radiation exposure and estimated cancer risk”. In: JAMA Pediatrics
167 (2017), pp. 6–14.
[27] Per Thunqvist et al. “Lung function at 6 and 18 months after preterm birth in relation
to severity of bronchopulmonary dysplasia”. In: Pediatric Pulmonology 50.10 (2015),
pp. 978–986.
[28] John Palmer et al. “Comprehensive ventilation/perfusion SPECT”. In: Journal of
Nuclear Medicine 42.8 (2001), pp. 1288–1294.
[29] Zhi-Pei Liang and Paul C Lauterbur. Principles of magnetic resonance imaging. SPIE
Optical Engineering Press Bellingham, 2000.
[30] Robert W Brown et al. Magnetic resonance imaging: physical principles and sequence
design. John Wiley & Sons, 2014.
[31] MA Bernstein. Handbook of MRI pulse sequences. Elsevier, 2004.
[32] Dwight G Nishimura. Principles of magnetic resonance imaging. Standford University, 1996.
[33] Isidor Isaac Rabi. “Space quantization in a gyrating magnetic field”. In: Physical
Review 51.8 (1937), p. 652.
[34] Felix Bloch. “Nuclear induction”. In: Physical Review 70.7-8 (1946), p. 460.
[35] Edward M Purcell, Henry Cutler Torrey, and Robert V Pound. “Resonance absorption
by nuclear magnetic moments in a solid”. In: Physical Review 69.1-2 (1946), p. 37.
[36] Richard R Ernst and Weston A Anderson. “Application of Fourier transform spectroscopy to magnetic resonance”. In: Review of Scientific Instruments 37.1 (1966),
pp. 93–102.
[37] Kurt W¨uthrich et al. “Sequential resonance assignments as a basis for determination
of spatial protein structures by high resolution proton nuclear magnetic resonance”.
In: Journal of Molecular Biology 155.3 (1982), pp. 311–319.
113

[38] Paul C Lauterbur. “Image formation by induced local interactions: examples employing nuclear magnetic resonance”. In: Nature 242.5394 (1973), pp. 190–191.
[39] Peter Mansfield. “Multi-planar image formation using NMR spin echoes”. In: Journal
of Physics C: Solid State Physics 10.3 (1977), p. L55.
[40] Paul A Bottomley et al. “A review of normal tissue hydrogen NMR relaxation times
and relaxation mechanisms from 1–100 MHz: dependence on tissue type, NMR frequency, temperature, species, excision, and age”. In: Medical physics 11.4 (1984),
pp. 425–448.
[41] MR Scanner Key Components. url: https://www.ee.duke.edu/jshorey (visited
on 12/20/2021).
[42] Govind B Chavhan et al. “Principles, techniques, and applications of T2*-based MR
imaging and its special applications”. In: Radiographics 29.5 (2009), pp. 1433–1449.
[43] Jeffrey A Fessler. “On NUFFT-based gridding for non-Cartesian MRI”. In: Journal
of Magnetic Resonance 188.2 (2007), pp. 191–195.
[44] James G Pipe, Ergun Ahunbay, and Padmanabhan Menon. “Effects of interleaf order for spiral MRI of dynamic processes”. In: Magnetic Resonance in Medicine 41.2
(1999), pp. 417–422.
[45] Weiran Deng, Benjamin Zahneisen, and V Andrew Stenger. “Rotated stack-of-spirals
partial acquisition for rapid volumetric parallel MRI”. In: Magnetic Resonance in
Medicine 76.1 (2016), pp. 127–135.
[46] Klaus Scheffler and J¨urgen Hennig. “Frequency resolved single-shot MR imaging using
stochastic k-space trajectories”. In: Magnetic Resonance in Medicine 35.4 (1996),
pp. 569–576.
[47] Krishna S Nayak and Dwight G Nishimura. “Randomized trajectories for reduced
aliasing artifact”. In: Proceedings of the 6th Annual Meeting of International Society
for Magnetic Resonance in Medicine, Sydney, Australia. 1998, p. 670.
[48] Chi-Ming Tsai and Dwight G Nishimura. “Reduced aliasing artifacts using variabledensity k-space sampling trajectories”. In: Magnetic Resonance in Medicine 43.3
(2000), pp. 452–458.
[49] Jan-Ray Liao et al. “Reduction of motion artifacts in cine MRI using variable-density
spiral trajectories”. In: Magnetic Resonance in Medicine 37.4 (1997), pp. 569–575.
114

[50] Jin Hyung Lee et al. “Fast 3D imaging using variable-density spiral trajectories with
applications to limb perfusion”. In: Magnetic Resonance in Medicine 50.6 (2003),
pp. 1276–1285.
[51] RAD 229 Lecture Notes. url: https://web.stanford.edu/class/rad229/Notes/
Lecture-11/ (visited on 12/20/2024).
[52] Anagha Deshmane et al. “Parallel MR imaging”. In: Journal of Magnetic Resonance
Imaging 36.1 (2012), pp. 55–72.
[53] Klaas P Pruessmann et al. “SENSE: sensitivity encoding for fast MRI”. In: Magnetic
Resonance in Medicine 42.5 (1999), pp. 952–962.
[54] Mark A Griswold et al. “Generalized autocalibrating partially parallel acquisitions
(GRAPPA)”. In: Magnetic Resonance in Medicine 47.6 (2002), pp. 1202–1210.
[55] Ze Wang, Jiongjiong Wang, and John A Detre. “Improved data reconstruction method
for GRAPPA”. In: Magnetic Resonance in Medicine 54.3 (2005), pp. 738–742.
[56] Felix A Breuer et al. “Dynamic autocalibrated parallel imaging using temporal GRAPPA
(TGRAPPA)”. In: Magnetic Resonance in Medicine 53.4 (2005), pp. 981–985.
[57] Michael Lustig, David Donoho, and John M Pauly. “Sparse MRI: The application of
compressed sensing for rapid MR imaging”. In: Magnetic Resonance in Medicine 58.6
(2007), pp. 1182–1195.
[58] Wa Edelstein et al. “The intrinsic signal-to-noise ratio in NMR imaging”. In: Magnetic
Resonance in Medicine 3.4 (1986), pp. 604–618.
[59] Jiangsheng Yu, Yiqun Xue, and Hee Kwon Song. “Comparison of lung T2* during
free-breathing at 1.5 T and 3.0 T with ultrashort echo time imaging”. In: Magnetic
Resonance in Medicine 66.1 (2011), pp. 248–254.
[60] Xiaohong Joe Zhou et al. “Concomitant magnetic-field-induced artifacts in axial echo
planar imaging”. In: Magnetic Resonance in Medicine 39.4 (1998), pp. 596–605.
[61] Mat A Bernstein et al. “Concomitant gradient terms in phase contrast MR: analysis
and correction”. In: Magnetic Resonance in Medicine 39.2 (1998), pp. 300–308.
[62] Tim Nielsen et al. “Compensation of concomitant Maxwell gradient effects in 3D
multi-element spatially selective RF excitation”. In: Proceedings 17th Scientific Meeting, International Society for Magnetic Resonance in Medicine, Honolulu. 2009, p. 2588.
115

[63] Ziwei Zhao, Nam G Lee, and Krishna S Nayak. “Multidimensional RF Pulse Design
with Known Spatial Encoding Imperfections”. In: Proceedings 29th Scientific Meeting,
International Society for Magnetic Resonance in Medicine, Online (2021), p3954.
[64] Nam G Lee et al. “MaxGIRF: Image reconstruction incorporating concomitant field
and gradient impulse response function effects”. In: Magnetic Resonance in Medicine
88.2 (2022), pp. 691–710.
[65] Zhixing Wang et al. “Concomitant magnetic-field compensation for 2D spiral-ring
turbo spin-echo imaging at 0.55 T and 1.5 T”. In: Magnetic Resonance in Medicine
90.2 (2023), pp. 552–568.
[66] Kevin F King et al. “Concomitant gradient field effects in spiral scans”. In: Magnetic
Resonance in Medicine 41.1 (1999), pp. 103–112.
[67] Paul T Weavers et al. “B0 concomitant field compensation for MRI systems employing asymmetric transverse gradient coils”. In: Magnetic Resonance in Medicine 79.3
(2018), pp. 1538–1544.
[68] Shengzhen Tao et al. “Gradient pre-emphasis to counteract first-order concomitant
fields on asymmetric MRI gradient systems”. In: Magnetic Resonance in Medicine
77.6 (2017), pp. 2250–2262.
[69] Yiping P Du, Xiaohong Joe Zhou, and Matt A Bernstein. “Correction of concomitant
magnetic field-induced image artifacts in nonaxial echo-planar imaging”. In: Magnetic
Resonance in Medicine 48.3 (2002), pp. 509–515.
[70] Andrew D Scott et al. “Speech MRI: morphology and function”. In: Physica Medica
30.6 (2014), pp. 604–618.
[71] Sajan Goud Lingala et al. “Recommendations for real-time speech MRI”. In: Journal
of Magnetic Resonance Imaging 43.1 (2016), pp. 28–44.
[72] Sajan Goud Lingala et al. “A fast and flexible MRI system for the study of dynamic
vocal tract shaping”. In: Magnetic Resonance in Medicine 77.1 (2017), pp. 112–125.
[73] Yoon-Chul Kim et al. “Novel 16-channel receive coil array for accelerated upper airway
MRI at 3 Tesla”. In: Magnetic Resonance in Medicine 65.6 (2011), pp. 1711–1717.
[74] Sajan Goud Lingala et al. “Feasibility of through-time spiral generalized autocalibrating partial parallel acquisition for low latency accelerated real-time MRI of speech”.
In: Magnetic Resonance in Medicine 78.6 (2017), pp. 2275–2282.
[75] Yoon-Chul Kim, Shrikanth S Narayanan, and Krishna S Nayak. “Accelerated 3D
MRI of vocal tract shaping using compressed sensing and parallel imaging”. In: 2009
116

IEEE International Conference on Acoustics, Speech and Signal Processing. IEEE.
2009, pp. 389–392.
[76] Aaron Niebergall et al. “Real-time MRI of speaking at a resolution of 33 ms: undersampled radial FLASH with nonlinear inverse reconstruction”. In: Magnetic Resonance in Medicine 69.2 (2013), pp. 477–485.
[77] Shrikanth Narayanan et al. “An approach to real-time magnetic resonance imaging
for speech production”. In: The Journal of the Acoustical Society of America 113.4
(2003), pp. 2258–2258.
[78] Yoon-Chul Kim, Shrikanth S Narayanan, and Krishna S Nayak. “Accelerated threedimensional upper airway MRI using compressed sensing”. In: Magnetic Resonance
in Medicine 61.6 (2009), pp. 1434–1440.
[79] Michael I Proctor et al. “Rapid semi-automatic segmentation of real-time magnetic
resonance images for parametric vocal tract analysis.” In: Proceedings of the 11th Annual Conference on the International Speech Communication Association (INTERSPEECH). 2010, pp. 1576–1579.
[80] Samuel Silva and Ant´onio Teixeira. “Unsupervised segmentation of the vocal tract
from real-time MRI sequences”. In: Computer Speech & Language 33.1 (2015), pp. 25–
46.
[81] Vikram Ramanarayanan et al. “Analysis of speech production real-time MRI”. In:
Computer Speech & Language 52 (2018), pp. 1–22.
[82] Maojing Fu et al. “High-frame-rate full-vocal-tract 3D dynamic speech imaging”. In:
Magnetic Resonance in Medicine 77.4 (2017), pp. 1619–1629.
[83] Michael Burdumy et al. “One-second MRI of a three-dimensional vocal tract to measure dynamic articulator modifications”. In: Journal of Magnetic Resonance Imaging
46.1 (2017), pp. 94–101.
[84] Yongwan Lim et al. “3D dynamic MRI of the vocal tract during natural speech”. In:
Magnetic Resonance in Medicine 81.3 (2019), pp. 1511–1520.
[85] Xue Feng et al. “Assessment of velopharyngeal function with dual-planar high-resolution
real-time spiral dynamic MRI”. In: Magnetic Resonance in Medicine 80.4 (2018),
pp. 1467–1474.
[86] Shrikanth Narayanan et al. “An approach to real-time magnetic resonance imaging
for speech production”. In: The Journal of the Acoustical Society of America 115.4
(2004), pp. 1771–1776.
117

[87] Bradley P Sutton et al. “Faster dynamic imaging of speech with field inhomogeneity corrected spiral fast low angle shot (FLASH) at 3 T”. In: Journal of Magnetic
Resonance Imaging 32.5 (2010), pp. 1228–1237.
[88] Andreia C Freitas et al. “Real-time speech MRI: Commercial Cartesian and nonCartesian sequences at 3T and feasibility of offline TGV reconstruction to visualise
velopharyngeal motion”. In: Physica Medica 46 (2018), pp. 96–103.
[89] Adrienne E Campbell-Washburn et al. “Opportunities in interventional and diagnostic imaging by using high-performance low-field-strength MRI”. In: Radiology 293.2
(2019), pp. 384–393.
[90] Yongwan Lim, Prakash Kumar, and Krishna S Nayak. “Speech production real-time
MRI at 0.55 T”. In: Magnetic Resonance in Medicine 91.1 (2024), pp. 337–343.
[91] Yoon-Chul Kim et al. “Improved imaging of lingual articulation using real-time multislice MRI”. In: Journal of Magnetic Resonance Imaging 35.4 (2012), pp. 943–948.
[92] Maojing Fu et al. “High-resolution dynamic speech imaging with deformation estimation”. In: 2015 37th Annual International Conference of the IEEE Engineering in
Medicine and Biology Society (EMBC). IEEE. 2015, pp. 1568–1571.
[93] Shrikanth Narayanan, Dani Byrd, and Abigail Kaun. “Geometry, kinematics, and
acoustics of Tamil liquid consonants”. In: The Journal of the Acoustical Society of
America 106.4 (1999), pp. 1993–2007.
[94] Felipe Sanchez et al. “Detection of solid and subsolid pulmonary nodules with lung
MRI: performance of UTE, T1 gradient-echo, and single-shot T2 fast spin echo”. In:
Cancer Imaging 23.1 (2023), p. 17.
[95] Lilian Lonzetti et al. “Magnetic resonance imaging of interstitial lung diseases: A
state-of-the-art review”. In: Respiratory Medicine 155 (2019), pp. 79–85.
[96] Mark Wielp¨utz and Hans-Ulrich Kauczor. “MRI of the lung: state of the art”. In:
Diagnostic and Interventional Radiology 18.4 (2012), p. 344.
[97] Mark O Wielp¨utz, Monika Eichinger, and Michael Puderbach. “Magnetic resonance
imaging of cystic fibrosis lung disease”. In: Journal of Thoracic Imaging 28.3 (2013),
pp. 151–159.
[98] Yi-Xiang J Wang et al. “Magnetic resonance imaging for lung cancer screen”. In:
Journal of Thoracic Disease 6.9 (2014), p. 1340.
[99] G Wilson Miller et al. “Advances in functional and structural imaging of the human
lung using proton MRI”. In: NMR in Biomedicine 27.12 (2014), pp. 1542–1556.
118

[100] Andreas Voskrebenzev and Jens Vogel-Claussen. “Proton MRI of the lung: how to
tame scarce protons and fast signal decay”. In: Journal of Magnetic Resonance Imaging 53.5 (2021), pp. 1344–1357.
[101] Harald E M¨oller et al. “MRI of the lungs using hyperpolarized noble gases”. In:
Magnetic Resonance in Medicine 47.6 (2002), pp. 1029–1051.
[102] Michael Salerno et al. “Hyperpolarized noble gas MR imaging of the lung: potential
clinical applications”. In: European Journal of Radiology 40.1 (2001), pp. 33–44.
[103] Robert R Edelman et al. “Noninvasive assessment of regional ventilation in the human
lung using oxygen–enhanced magnetic resonance imaging”. In: Nature Medicine 2.11
(1996), pp. 1236–1239.
[104] Orso Pusterla, Grzegorz Bauman, and Oliver Bieri. “Three-dimensional oxygen-enhanced
MRI of the human lung at 1.5 T with ultra-fast balanced steady-state free precession”.
In: Magnetic Resonance in Medicine 79.1 (2018), pp. 246–255.
[105] Tanguy Boucneau et al. “3D magnetic resonance spirometry”. In: Scientific Reports
10.1 (2020), p. 9649.
[106] Julian Glandorf et al. “Comparison of phase-resolved functional lung (PREFUL) MRI
derived perfusion and ventilation parameters at 1.5 T and 3T in healthy volunteers”.
In: PloS One 15.12 (2020), e0244638.
[107] Ahsan Javed et al. “Self-gated 3D stack-of-spirals UTE pulmonary imaging at 0.55
T”. In: Magnetic Resonance in Medicine 87.4 (2022), pp. 1784–1798.
[108] Ahsan Javed et al. “Increasing the scan-efficiency of pulmonary imaging at 0. 55 T
using iterative concomitant field and motion-corrected reconstruction”. In: Magnetic
Resonance in Medicine 92.1 (2024), pp. 173–185.
[109] Nam G Lee et al. “Replication of the bSTAR sequence and open-source implementation”. In: Magnetic Resonance in Medicine 91.4 (2024), pp. 1464–1477.
[110] Ye Tian et al. “Single breath-hold volumetric lung imaging at 0. 55T using stack-ofspiral (SoS) out-in balanced SSFP”. In: Magnetic Resonance in Medicine (2024).
[111] Nara S Higano et al. “Retrospective respiratory self-gating and removal of bulk motion
in pulmonary UTE MRI of neonates and adults”. In: Magnetic Resonance in Medicine
77.3 (2017), pp. 1284–1295.
[112] Samal Munidasa et al. “Free-breathing MRI for monitoring ventilation changes following antibiotic treatment of pulmonary exacerbations in paediatric cystic fibrosis”.
In: European Respiratory Journal 57.4 (2021).
119

[113] Simon L´evy et al. “Free-breathing low-field MRI of the lungs detects functional alterations associated with persistent symptoms after COVID-19 infection”. In: Investigative Radiology 57.11 (2022), pp. 742–751.
[114] Justus E Roos et al. “Hyperpolarized gas MR imaging: technique and applications”.
In: Magnetic Resonance Imaging Clinics 23.2 (2015), pp. 217–229.
[115] Ipshita Bhattacharya et al. “Oxygen-enhanced functional lung imaging using a contemporary 0.55 T MRI system”. In: NMR in Biomedicine 34.8 (2021), e4562.
[116] Rafael Heiss et al. “Pulmonary dysfunction after pediatric COVID-19”. In: Radiology
306.3 (2022), e221250.
[117] Andrew D Hahn et al. “Pulmonary MRI of neonates in the intensive care unit using
3D ultrashort echo time and a small footprint MRI system”. In: Journal of Magnetic
Resonance Imaging 45.2 (2017), pp. 463–471.
[118] Kevin M Johnson et al. “Optimized 3D ultrashort echo time pulmonary MRI”. In:
Magnetic Resonance in Medicine 70.5 (2013), pp. 1241–1250.
[119] Julius F Heidenreich et al. “Three-dimensional ultrashort echo time MRI for functional lung imaging in cystic fibrosis”. In: Radiology 296.1 (2020), pp. 191–199.
[120] Fei Tan et al. “Motion-compensated low-rank reconstruction for simultaneous structural and functional UTE lung MRI”. In: Magnetic Resonance in Medicine 90.3
(2023), pp. 1101–1113.
[121] Corona Metz et al. “Pulmonary imaging of immunocompromised patients during
hematopoietic stem cell transplantation using non-contrast-enhanced three-dimensional
ultrashort echo time (3D-UTE) MRI”. In: R¨oFo-Fortschritte auf dem Gebiet der
R¨ontgenstrahlen und der bildgebenden Verfahren. Vol. 194. 01. Georg Thieme Verlag KG. 2022, pp. 39–48.
[122] Julia Geiger et al. “Clinical application of ultrashort echo-time MRI for lung pathologies in children”. In: Clinical Radiology 76.9 (2021), 708–e9.
[123] Hiroto Hatabu et al. “Expanding applications of pulmonary MRI in the clinical evaluation of lung disorders: Fleischner Society position paper”. In: Radiology 297.2 (2020),
pp. 286–301.
[124] Adrienne E Campbell-Washburn et al. “T2-weighted lung imaging using a 0.55-T
MRI system”. In: Radiology: Cardiothoracic Imaging 3.3 (2021), e200611.
120

[125] Ipshita Bhattacharya et al. “Assessment of lung structure and regional function using
0.55 T MRI in patients with lymphangioleiomyomatosis”. In: Investigative Radiology
57.3 (2022), pp. 178–186.
[126] John F Schenck. “The role of magnetic susceptibility in magnetic resonance imaging:
MRI magnetic compatibility of the first and second kinds”. In: Medical Physics 23.6
(1996), pp. 815–850.
[127] Bochao Li et al. “Lung parenchyma transverse relaxation rates at 0.55 T”. In: Magnetic Resonance in Medicine 89.4 (2023), pp. 1522–1530.
[128] Yicun Wang et al. “FMRI based on transition-band balanced SSFP in comparison
with EPI on a high-performance 0.55 T scanner”. In: Magnetic Resonance in Medicine
85.6 (2021), pp. 3196–3210.
[129] Sandra J England. “Current techniques for assessing pulmonary function in the newborn and infant: advantages and limitations”. In: Pediatric Pulmonology 4.1 (1988),
pp. 48–53.
[130] Grzegorz Bauman et al. “Non-contrast-enhanced perfusion and ventilation assessment
of the human lung by means of Fourier decomposition in proton MRI”. In: Magnetic
Resonance in Medicine 62.3 (2009), pp. 656–664.
[131] Grzegorz Bauman and Oliver Bieri. “Matrix pencil decomposition of time-resolved
proton MRI for robust and improved assessment of pulmonary ventilation and perfusion”. In: Magnetic Resonance in Medicine 77.1 (2017), pp. 336–342.
[132] Andreas Voskrebenzev et al. “Feasibility of quantitative regional ventilation and perfusion mapping with phase-resolved functional lung (PREFUL) MRI in healthy volunteers and COPD, CTEPH, and CF patients”. In: Magnetic Resonance in Medicine
79.4 (2018), pp. 2306–2314.
[133] F Klimeˇs et al. “Free-breathing quantification of regional ventilation derived by phaseresolved functional lung (PREFUL) MRI”. In: NMR in Biomedicine 32.6 (2019),
e4088.
[134] Gesa H P¨ohler et al. “Repeatability of phase-resolved functional lung (PREFUL)-
MRI ventilation and perfusion parameters in healthy subjects and COPD Patients”.
In: Journal of Magnetic Resonance Imaging 53.3 (2021), pp. 915–927.
[135] Filip Klimeˇs et al. “3D phase-resolved functional lung ventilation MR imaging in
healthy volunteers and patients with chronic pulmonary disease”. In: Magnetic Resonance in Medicine 85.2 (2021), pp. 912–925.
121

[136] Juan M Santos, Graham A Wright, and John M Pauly. “Flexible real-time magnetic
resonance imaging framework”. In: The 26th Annual International Conference of the
IEEE Engineering in Medicine and Biology Society. Vol. 1. IEEE. 2004, pp. 1048–
1051.
[137] Yoon-Chul Kim, Shrikanth S Narayanan, and Krishna S Nayak. “Flexible retrospective selection of temporal resolution in real-time speech MRI using a golden-ratio
spiral view order”. In: Magnetic Resonance in Medicine 65.5 (2011), pp. 1365–1371.
[138] Adam B Kerr et al. “Real-time interactive MRI on a conventional scanner”. In: Magnetic Resonance in Medicine 38.3 (1997), pp. 355–367.
[139] Stephen Boyd et al. “Distributed optimization and statistical learning via the alternating direction method of multipliers”. In: Foundations and Trends in Machine
Learning 3.1 (2011), pp. 1–122.
[140] Martin Uecker et al. “The BART toolbox for computational magnetic resonance imaging”. In: Proceedings 23rd Scientific Meeting, International Society for Magnetic Resonance in Medicine, Toronto, Canada. Vol. 23. 2015, p. 2486.
[141] Martin Uecker et al. “ESPIRiT—an eigenvalue approach to autocalibrating parallel
MRI: where SENSE meets GRAPPA”. In: Magnetic Resonance in Medicine 71.3
(2014), pp. 990–1001.
[142] Amanda L. Miller and Kenneth B. Finch. “Corrected high–frame rate anchored ultrasound with software alignment”. In: Journal of Speech, Language, and Hearing
Research 54.2 (2011), pp. 471–486.
[143] Adam C Lammert et al. “Vocal tract cross-distance estimation from real-time MRI using region-of-interest analysis.” In: Proceedings of the 14th Annual Conference on the
International Speech Communication Association (INTERSPEECH). 2013, pp. 959–
962.
[144] Miran Oh et al. “Enriching the understanding of glottalic consonant production:
Vertical larynx movement in Hausa ejectives and implosives”. In: The Journal of the
Acoustical Society of America 144.3 (2018), pp. 1940–1941.
[145] Nimisha Patil et al. “Comparison of Basic Beatboxing Articulations Between Expert
and Novice Artists Using Real-Time Magnetic Resonance Imaging.” In: Proceedings of
the 18th Annual Conference on the International Speech Communication Association
(INTERSPEECH). 2017, pp. 2277–2281.
[146] Yongwan Lim et al. “Dynamic off-resonance correction for spiral real-time MRI of
speech”. In: Magnetic Resonance in Medicine 81.1 (2019), pp. 234–246.
122

[147] David Y Zeng et al. “Deep residual network for off-resonance artifact correction
with application to pediatric body MRA with 3D cones”. In: Magnetic Resonance
in Medicine 82.4 (2019), pp. 1398–1411.
[148] Yongwan Lim et al. “Deblurring for spiral real-time MRI using convolutional neural
networks”. In: Magnetic Resonance in Medicine 84.6 (2020), pp. 3438–3452.
[149] Oliver Bieri and Klaus Scheffler. “Fundamentals of balanced steady state free precession MRI”. In: Journal of Magnetic Resonance Imaging 38.1 (2013), pp. 2–11.
[150] Wenwen Jiang et al. “Motion robust high resolution 3D free-breathing pulmonary
MRI using dynamic 3D image self-navigator”. In: Magnetic Resonance in Medicine
79.6 (2018), pp. 2954–2967.
[151] Lalit Maurya et al. “Non-contact breathing rate monitoring in newborns: A review”.
In: Computers in Biology and Medicine 132 (2021), p. 104321.
[152] Oliver Bieri, Michael Markl, and Klaus Scheffler. “Analysis and compensation of
eddy currents in balanced SSFP”. In: Magnetic Resonance in Medicine 54.1 (2005),
pp. 129–137.
[153] K Scheffler and J Hennig. “Eddy current optimized phase encoding schemes to reduce
artifacts in balanced SSFP imaging”. In: Proceedings of the 11th Annual Meeting of
International Society for Magnetic Resonance in Medicine, Toronto, Canada. 2003,
p. 294.
[154] Stefan Wundrak et al. “Golden ratio sparse MRI using tiny golden angles”. In: Magnetic Resonance in Medicine 75.6 (2016), pp. 2372–2378.
[155] Patrick Le Roux. “Simplified model and stabilization of SSFP sequences”. In: Journal
of Magnetic Resonance 163.1 (2003), pp. 23–37.
[156] Kelvin J Layton et al. “Pulseq: a rapid and hardware-independent pulse sequence
prototyping framework”. In: Magnetic Resonance in Medicine 77.4 (2017), pp. 1544–
1552.
[157] Lorenzo Di Sopra et al. “An automated approach to fully self-gated free-running cardiac and respiratory motion-resolved 5D whole-heart MRI”. In: Magnetic Resonance
in Medicine 82.6 (2019), pp. 2118–2132.
[158] Souheil J Inati et al. “ISMRM Raw data format: A proposed standard for MRI raw
datasets”. In: Magnetic Resonance in Medicine 77.1 (2017), pp. 411–421.
[159] Xia Zhao et al. “Impact of gradient imperfections on bone water quantification with
UTE MRI”. In: Magnetic Resonance in Medicine 84.4 (2020), pp. 2034–2047.
123

[160] Ye Tian et al. “Aliasing artifact reduction in spiral real-time MRI”. In: Magnetic
Resonance in Medicine 86.2 (2021), pp. 916–925.
[161] Sedi Bartz, Rub´en Campoy, and Hung M Phan. “An adaptive alternating direction
method of multipliers”. In: Journal of Optimization Theory and Applications 195.3
(2022), pp. 1019–1055.
[162] Jakob Wasserthal et al. “TotalSegmentator: robust segmentation of 104 anatomic
structures in CT images”. In: Radiology: Artificial Intelligence 5.5 (2023).
[163] J-P Thirion. “Image matching as a diffusion process: an analogy with Maxwell’s
demons”. In: Medical Image Analysis 2.3 (1998), pp. 243–260.
[164] Thomas Goodhart et al. “Dynamic changes in lung water density and volume following supine body positioning”. In: Magnetic Resonance in Medicine 91.6 (2024),
pp. 2612–2620.
[165] Felicia Seemann et al. “Imaging gravity-induced lung water redistribution with automated inline processing at 0.55 T cardiovascular magnetic resonance”. In: Journal of
Cardiovascular Magnetic Resonance 24.1 (2022), p. 35.
[166] Ashley G Anderson III et al. “Intra-Scan Center Frequency Drift Correction for 3D
Spiral Exams”. In: Proceedings 25th Scientific Meeting, International Society for Magnetic Resonance in Medicine, Honolulu. Vol. 25. 2017, p. 1482.
[167] Johan Petersson et al. “Posture primarily affects lung tissue distribution with minor
effect on blood flow and ventilation”. In: Respiratory physiology & neurobiology 156.3
(2007), pp. 293–303.
[168] Guido Musch et al. “Topographical distribution of pulmonary perfusion and ventilation, assessed by PET in supine and prone humans”. In: Journal of applied physiology
93.5 (2002), pp. 1841–1851.
[169] Samal Munidasa et al. “Comparison of 3D UTE free-breathing lung MRI with hyperpolarized 129Xe MRI in pediatric cystic fibrosis”. In: Magnetic Resonance in Medicine
(2024).
[170] J. Pauly. “A k-space analysis of small-tip-angle excitation”. In: Journal of Magnetic
Resonance 213.2 (2011), pp. 544–557.
[171] Kay Nehrke et al. “On the performance and accuracy of 2D navigator pulses”. In:
Magnetic Resonance Imaging 17.8 (1999), pp. 1173–1181.
[172] Marcus T Alley et al. “Angiographic imaging with 2D RF pulses”. In: Magnetic
Resonance in Medicine 37.2 (1997), pp. 260–267.
124

[173] Susanne Rieseberg, Jens Frahm, and J¨urgen Finsterbusch. “Two-dimensional spatiallyselective RF excitation pulses in echo-planar imaging”. In: Magnetic Resonance in
Medicine 47.6 (2002), pp. 1186–1193.
[174] Chun-yu Yip, Jeffrey A Fessler, and Douglas C Noll. “Iterative RF pulse design for
multidimensional, small-tip-angle selective excitation”. In: Magnetic Resonance in
Medicine 54.4 (2005), pp. 908–917.
[175] Ulrich Katscher et al. “Transmit sense”. In: Magnetic Resonance in Medicine 49.1
(2003), pp. 144–150.
[176] Yudong Zhu. “Parallel excitation with an array of transmit coils”. In: Magnetic Resonance in Medicine 51.4 (2004), pp. 775–784.
[177] Suwit Saekho et al. “Small tip angle three-dimensional tailored radiofrequency slabselect pulse for reduced B1 inhomogeneity at 3 T”. In: Magnetic Resonance in Medicine
53.2 (2005), pp. 479–484.
[178] Francesco Padormo et al. “Parallel transmission for ultrahigh-field imaging”. In: NMR
in Biomedicine 29.9 (2016), pp. 1145–1161.
[179] William A Grissom et al. “Fast large-tip-angle multidimensional and parallel RF pulse
design in MRI”. In: IEEE Transactions on Medical Imaging 28.10 (2009), pp. 1548–
1559.
[180] William Grissom et al. “Spatial domain method for the design of RF pulses in multicoil
parallel excitation”. In: Magnetic Resonance in Medicine 56.3 (2006), pp. 620–629.
[181] Steven Conolly et al. “Variable-rate selective excitation”. In: Journal of Magnetic
Resonance 78.3 (1988), pp. 440–458.
[182] Neville D Gai and Yuval Zur. “Design and optimization for variable rate selective
excitation using an analytic RF scaling function”. In: Journal of Magnetic Resonance
189.1 (2007), pp. 78–89.
[183] Charles H Cunningham, Graham A Wright, and Michael L Wood. “High-order multiband encoding in the heart”. In: Magnetic Resonance in Medicine 48.4 (2002), pp. 689–
698.
[184] Dan Xu, Kevin F King, and Zhi-Pei Liang. “Variable slew-rate spiral design: theory
and application to peak B1 amplitude reduction in 2D RF pulse design”. In: Magnetic
Resonance in Medicine 58.4 (2007), pp. 835–842.
[185] Brian A Hargreaves et al. “Variable-rate selective excitation for rapid MRI sequences”.
In: Magnetic Resonance in Medicine 52.3 (2004), pp. 590–597.
125

[186] Daeho Lee et al. “Time-optimal design for multidimensional and parallel transmit
variable-rate selective excitation”. In: Magnetic Resonance in Medicine 61.6 (2009),
pp. 1471–1479.
[187] Michael Lustig, Seung-Jean Kim, and John M Pauly. “A fast method for designing
time-optimal gradient waveforms for arbitrary k-space trajectories”. In: IEEE Transactions on Medical Imaging 27.6 (2008), pp. 866–873.
[188] Shaihan J Malik and Joseph V Hajnal. “Phase relaxed localized excitation pulses
for inner volume fast spin echo imaging”. In: Magnetic Resonance in Medicine 76.3
(2016), pp. 848–861.
[189] CH Boesch, R Gruetter, and E Martin. “Temporal and spatial analysis of fields generated by eddy currents in superconducting magnets: optimization of corrections and
quantitative characterization of magnet/gradient systems”. In: Magnetic Resonance
in Medicine 20.2 (1991), pp. 268–284.
[190] Q Liu, DG Hughes, and PS Allen. “Quantitative characterization of the eddy current
fields in a 40-cm bore superconducting magnet”. In: Magnetic Resonance in Medicine
31.1 (1994), pp. 73–76.
[191] Yuhua Wu et al. “Gradient-induced acoustic and magnetic field fluctuations in a 4T
whole-body MR imager”. In: Magnetic Resonance in Medicine 44.4 (2000), pp. 532–
536.
[192] Bernd U Foerster, Dardo Tomasi, and Elisabeth C Caparelli. “Magnetic field shift
due to mechanical vibration in functional magnetic resonance imaging”. In: Magnetic
Resonance in Medicine 54.5 (2005), pp. 1261–1267.
[193] Jeff H Duyn et al. “Simple correction method for k-space trajectory deviations in
MRI”. In: Journal of Magnetic Resonance 132.1 (1998), pp. 150–153.
[194] David J Goodyear et al. “Single point measurements of magnetic field gradient waveform”. In: Journal of Magnetic Resonance 163.1 (2003), pp. 1–7.
[195] Marcus T Alley, Gary H Glover, and Norbert J Pelc. “Gradient characterization using
a Fourier-transform technique”. In: Magnetic Resonance in Medicine 39.4 (1998),
pp. 581–587.
[196] Signe J Vannesjo et al. “Gradient system characterization by impulse response measurements with a dynamic field camera”. In: Magnetic Resonance in Medicine 69.2
(2013), pp. 583–593.
[197] Mustafa C¸ avu¸so˘glu et al. “VERSE-guided parallel RF excitations using dynamic field
correction”. In: NMR in Biomedicine 30.6 (2017), e3697.
126

[198] Hai Zheng et al. “Parallel transmission RF pulse design for eddy current correction
at ultra high field”. In: Journal of Magnetic Resonance 221 (2012), pp. 139–146.
[199] Norman F Ramsey. “Resonance transitions induced by perturbations at two or more
different frequencies”. In: Physical Review 100.4 (1955), p. 1191.
[200] F Bloch and A Siegert. “Magnetic resonance for nonrotating fields”. In: Physical
Review 57.6 (1940), p. 522.
[201] Laura I Sacolick et al. “B1 mapping by Bloch-Siegert shift”. In: Magnetic Resonance
in Medicine 63.5 (2010), pp. 1315–1322.
[202] Peter B¨ornert and Bernd Aldefeld. “On spatially selective RF excitation and its analogy with spiral MR image acquisition”. In: Magnetic Resonance Materials in Physics,
Biology and Medicine 7 (1998), pp. 166–178.
[203] John David Jackson and Ronald F Fox. Classical electrodynamics. 1999.
[204] Peter M Joseph, Leon Axel, and Matthew O’Donnell. “Potential problems with selective pulses in NMR imaging systems”. In: Medical Physics 11.6 (1984), pp. 772–
777.
[205] Adrienne E Campbell-Washburn et al. “Real-time distortion correction of spiral and
echo planar images using the gradient system impulse response function”. In: Magnetic
Resonance in Medicine 75.6 (2016), pp. 2278–2285.
[206] Sohae Chung et al. “Rapid B1+ mapping using a preconditioning RF pulse with
TurboFLASH readout”. In: Magnetic Resonance in Medicine 64.2 (2010), pp. 439–
446.
[207] W Patricia Bandettini et al. “A comparison of cine CMR imaging at 0.55 T and 1.5
T”. In: Journal of Cardiovascular Magnetic Resonance 22.1 (2020), p. 37.
[208] Matthew C Restivo et al. “Efficient spiral in-out and EPI balanced steady-state free
precession cine imaging using a high-performance 0.55 T MRI”. In: Magnetic Resonance in Medicine 84.5 (2020), pp. 2364–2375.
[209] Christopher M Collins, Shizhe Li, and Michael B Smith. “SAR and B1 field distributions in a heterogeneous human head model within a birdcage coil”. In: Magnetic
Resonance in Medicine 40.6 (1998), pp. 847–856.
[210] Jianming Jin and Ji Chen. “On the SAR and field inhomogeneity of birdcage coils
loaded with the human head”. In: Magnetic Resonance in Medicine 38.6 (1997),
pp. 953–963.
127

[211] Martin Uecker et al. “Real-time MRI at a resolution of 20 ms”. In: NMR in Biomedicine
23.8 (2010), pp. 986–994.
[212] Manuel Freiberger et al. “The agile library for biomedical image reconstruction using
GPU acceleration”. In: Computing in Science & Engineering 15.1 (2012), pp. 34–44.
[213] Xiao-Long Wu et al. “Impatient MRI: Illinois Massively Parallel Acceleration Toolkit
for image reconstruction with enhanced throughput in MRI”. In: 2011 IEEE International Symposium on Biomedical Imaging: from nano to macro. 2011, pp. 69–72.
[214] Neil B Mascarenhas et al. “Fast 3D cine steady-state free precession imaging with
sensitivity encoding for assessment of left ventricular function in a single breathhold”. In: American Journal of Roentgenology 187.5 (2006), pp. 1235–1239.
[215] Kestutis J Barkauskas et al. “Quantification of left ventricular functional parameter values using 3D spiral bSSFP and through-time Non-Cartesian GRAPPA”. In:
Journal of Cardiovascular Magnetic Resonance 16.1 (2014), p. 65.
[216] Javier Montalt-Tordera et al. “Rapid 3D whole-heart cine imaging using golden ratio
stack of spirals”. In: Magnetic Resonance Imaging 72 (2020), pp. 1–7.
[217] Steven C Persak et al. “Noninvasive estimation of pharyngeal airway resistance and
compliance in children based on volume-gated dynamic MRI and computational fluid
dynamics”. In: Journal of Applied Physiology 111.6 (2011), pp. 1819–1827.
[218] Yoon-Chul Kim et al. “Real-time 3D magnetic resonance imaging of the pharyngeal
airway in sleep apnea”. In: Magnetic Resonance in Medicine 71.4 (2014), pp. 1501–
1510.
[219] Stephen S Henrichon et al. “Dynamic MRI of the wrist in less than 20 seconds: normal
midcarpal motion and reader reliability”. In: Skeletal Radiology 49.2 (2020), pp. 241–
248.
[220] Calvin B Shaw et al. “Real-time three-dimensional MRI for the assessment of dynamic
carpal instability”. In: PLoS One 14.9 (2019), e0222704.
[221] Ziwu Zhou et al. “Golden-ratio rotated stack-of-stars acquisition for improved volumetric MRI”. In: Magnetic Resonance in Medicine 78.6 (2017), pp. 2290–2298.
[222] Jason Kraig Mendes et al. “Quantitative 3D myocardial perfusion with an efficient
arterial input function”. In: Magnetic Resonance in Medicine 83.6 (2020), pp. 1949–
1963.
128

[223] Kai Tobias Block et al. “Towards routine clinical use of radial stack-of-stars 3D
gradient-echo sequences for reducing motion sensitivity”. In: Investigative Magnetic
Resonance Imaging 18.2 (2014), pp. 87–106.
[224] Ziwei Zhao et al. “Lung perfusion at 0.55 T using ASL: Feasibility and Initial Results”.
In: Proceedings 30th Scientific Meeting, International Society for Magnetic Resonance
in Medicine, London, England, UK. 2022, p. 3365.
129

Abstract (if available)

Abstract Magnetic resonance imaging (MRI), without ionizing radiation, is a non-invasive medical imaging modality that can visualize structural anatomy as well as functional responses. three-dimensional (3D) dynamic MRI provides an efficient solution to visualize dynamic movements and to reveal functional responses volumetrically. This dissertation is focused on 3D dynamic MRI methods on vocal tract during real-time, and 3D pulmonary ventilation mapping during free-breathing. It targets on improving three limitations faced by these applications: limited spatial-temporal resolution, low parenchyma signals, and system imperfections.
First, I describe a 3D real-time MRI of speech production at 1.5T with improved spatio- temporal sharpness using randomized, variable-density, stack-of-spiral (SOS) sampling combined with a 3D spatiotemporally constrained reconstruction. Real-time MRI has emerged as one of the most powerful tools for studying human speech production. In this work, I evaluated five candidate (k, t) sampling strategies using a previously proposed gradient- echo stack-of-spiral sequence and a 3D constrained reconstruction with spatial and temporal penalties. The strategy yielding highest image quality was chosen as the proposed method. The proposed method provides the best qualities with the improved spatiotemporal resolution for better visualization of fast lip and tongue movements. The significant improvements are demonstrated in tongue boundary sharpness (p<0.001) during normal and 1.5× speeded speech. This allows researchers to demonstrate both functional and morphological aspects of the vocal tract during speech and to understand the complex spatiotemporal coordination of upper airway structures in motion.
Next, I describe an efficient free-breathing 3D functional lung imaging using SOS Out- in balanced steady-state free-precession (bSSFP). MRI, as an alternative for computed tomography, also provides functional response in screening and diagnostic assessment for many lung diseases, such as cancer, emphysema, and chronic obstructive pulmonary disease, without ionizing radiation. The proposed method, for the first time, demonstrates 3D regional pulmonary ventilation mapping at 0.55T within 5min scan using 3D bSSFP acquisition, constrained reconstruction and non-rigid registration. It overcome the challenges of low parenchyma signals of lung MRI in mid- and high-field strengths (≥ 1.5T). In healthy volunteers, the SOS out-in lung images provided sufficient vessel-parenchyma contrast and boundary sharpness to support accurate ventilation estimation. Regional ventilation measurements from 3D SOS out-in demonstrated good repeatability (relative differences<10%). Ventilation maps from 3D SOS out-in strongly correlated with Phase-resolved Functional Lung (PREFUL) MRI on a slice-matched basis, as well as with global tidal volume (R2 >0.7, p<0.001). Therefore, ventilation measurements are sensitive, consistent, and in good agreement with PREFUL and spirometry.
Lastly, I describe my efforts in developing a small-tip multidimensional radiofrequency (RF) pulse design procedure that incorporates linear time-invariant gradient imperfections and concomitant field effects. I developed an extension of the small-tip excitation k-space formalism, where concomitant fields were approximated as a Bloch-Siegert shift in the rotating frame. This was evaluated using realistic simulations of 2D selective excitation at various field strengths (0.2T, 0.55T, 1.5T,3T, and 7T) with single and parallel transmit. The extended formalism provides improved 2D excitation profiles in all scenarios simulated, compared against the original formalism. The proposed method corrects the concomitant field effects on 2D selective excitations for B0 ≥ 0.2T when the magnitude of the B0 is far larger than that of nonrotating concomitant fields. The proposed method providing sharper and more accurate excitation profiles at off-isocenter distances up to 15 cm. The impact of the proposed method is greatest in scenarios where concomitant fields are substantial, such as low field strengths and off-isocenter. This could be particularly important for contemporary low-field MRI systems with high-performance gradients, with potential applications in off-isocenter navigators during real-time volumetric imaging.
A classic description of the principles of MRI, as well as the summary and future remarks about this dissertation are also presented.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Fast flexible dynamic three-dimensional magnetic resonance imaging

PDF

New methods for carotid MRI

PDF

Correction, coregistration and connectivity analysis of multi-contrast brain MRI

PDF

Improved brain dynamic contrast enhanced MRI using model-based reconstruction

PDF

Shift-invariant autoregressive reconstruction for MRI

PDF

Estimating liver iron non-invasively with high-field MRI

PDF

Efficient and accurate 3D FISP-MRF at 0.55 Tesla

PDF

Seeing sleep: real-time MRI methods for the evaluation of sleep apnea

PDF

New theory and methods for accelerated MRI reconstruction

PDF

Functional real-time MRI of the upper airway

PDF

Optimization methods and algorithms for constrained magnetic resonance imaging

PDF

Methods for improving reliability and consistency in diffusion MRI analysis

PDF

Fast upper airway MRI of speech

PDF

Novel computational techniques for connectome analysis based on diffusion MRI

PDF

Dynamic cardiovascular magnetic resonance imaging for improved assessment of ischemic heart disease

PDF

Dynamic functional magnetic resonance imaging to localize activated neurons

PDF

3D object detection in industrial site point clouds

PDF

Articulated human body deformation from in-vivo 3d image scans

PDF

High-dimensional magnetic resonance imaging of microstructure

PDF

Novel theoretical characterization and optimization of experimental efficiency for diffusion MRI

Asset Metadata

Creator Zhao, Ziwei (author)

Core Title Technology for improved 3D dynamic MRI

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Electrical and Computer Engineering

Degree Conferral Date 2024-12

Publication Date 01/15/2025

Defense Date 12/14/2025

Publisher Los Angeles, California (original), University of Southern California (original), University of Southern California. Libraries (digital)

Tag concomitant fields,low field MRI,Lung,OAI-PMH Harvest,speech,three-dimensional dynamic iImaging

Format theses (aat)

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Nayak, Krishna (committee chair), Leahy, Richard (committee member), Wood, John (committee member)

Creator Email ziweiz@usc.edu,ziweizhao626@gmail.com

Unique identifier UC11399FAXB

Identifier etd-ZhaoZiwei-13760.pdf (filename)

Legacy Identifier etd-ZhaoZiwei-13760

Document Type Dissertation

Format theses (aat)

Rights Zhao, Ziwei

Internet Media Type application/pdf

Type texts

Source 20250115-usctheses-batch-1235 (batch), 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 author, as the original true and official version of the work, but does not grant the reader permission to use the work if the desired use is covered by copyright. It is the author, as rights holder, who must provide use permission if such use is covered by copyright.

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 90089-2810, USA

Repository Email cisadmin@lib.usc.edu