Close
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected
Invert selection
Deselect all
Deselect all
Click here to refresh results
Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Dynamic functional magnetic resonance imaging to localize activated neurons
(USC Thesis Other)
Dynamic functional magnetic resonance imaging to localize activated neurons
PDF
Download
Share
Open document
Flip pages
Contact Us
Contact Us
Copy asset link
Request this asset
Transcript (if available)
Content
DYNAMIC FUNCTIONAL MAGNETIC RESONANCE IMAGING
TO LOCALIZE ACTIVATED NEURONS
by
Witaya Sungkarat
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
(BIOMEDICAL ENGINEERING)
May 2007
Copyright 2007 Witaya Sungkarat
ii
Acknowledgements
I would like to thank my supervisor, Prof. Manbir Singh, for the warm support
and wonderful advice that he has provided throughout the duration of my Ph.D. I am
very grateful to my Ph.D. guidance committee, Professors Walter Wolf, Jesse T. Yen,
Norberto M. Grzywacz, and K. Kirk Shung, for their time, questions and advice. I
also thank Prof. Richard M. Leahy for his kindness and advice, especially, to study
useful electrical engineering classes beyond the master degree level to prepare for
Ph.D. work. A special thanks to Nanthana Tinroongroj for her help to make a nice
graphic user interface for our software. I would also like to thank all the people with
whom I have worked over the past years, and especially those who have patiently
undergone MRI scans for my experiments. Without them, none of this work could
have been done. Finally I would like to thank the Department of Biomedical
Engineering and the University of Southern California.
iii
Table of Contents
Acknowledgements ........................................................................................................ii
List of Tables.................................................................................................................vi
List of Figures...............................................................................................................vii
Abbreviations ................................................................................................................xi
Abstract........................................................................................................................xiii
Chapter 1 Introduction.................................................................................................1
Chapter 2 Basics..........................................................................................................3
2.1 Magnetic Resonance Imaging (MRI) ................................................................3
2.1.1 Nuclear Magnetism, Spin, and Resonance..................................................4
2.1.2 Nuclear Magnetic Resonance Excitation.....................................................5
2.1.3 T
1
and T
2
Relaxation Time ..........................................................................5
2.1.4 Bloch Equations...........................................................................................6
2.1.5 Spectroscopy and Imaging ..........................................................................6
2.1.6 Decay Time constant, T
2
* ...........................................................................6
2.1.7 Echo Planar Imaging (EPI)..........................................................................7
2.2 Functional Magnetic Resonance Imaging (fMRI).............................................8
2.2.1 Brain Physiology: Regional Vascular Responses........................................9
2.2.2 Blood Oxygen Level Dependent (BOLD) Signal .......................................9
2.2.3 Arterial Spin Labeling (ASL)....................................................................11
2.3 Detection Methods in fMRI.............................................................................12
2.3.1 Hemodynamic Response Function (HRF) Models....................................13
2.3.2 Correlation Coefficients ............................................................................14
2.3.3 General Linear Model (GLM)...................................................................15
2.3.4 Statistical Parametric Map (SPM) .............................................................16
Chapter 3 Dynamics of fMRI Signals .......................................................................17
3.1 Neuro–vascular Anatomical Background........................................................17
3.2 Big Vein BOLD Contributions........................................................................18
3.3 Multi-focal Brain Activations..........................................................................19
3.4 Regional Vasculature Variation and BOLD Specificity .................................20
3.5 Spatiotemporal Resolution of fMRI ................................................................20
3.5.1 Spatial Resolution......................................................................................21
3.5.2 Temporal Resolution .................................................................................21
3.6 BOLD fMRI Signal Flow................................................................................22
3.7 Advantages of Dynamic fMRI (dfMRI)..........................................................22
3.7.1 Space: Capillaries, Venules, Veins, or a Sinus..........................................22
iv
3.7.2 Time: Fast or Slow Responses ..................................................................22
3.7.3 Spatiotemporal Network............................................................................23
3.8 Potential Problems and Solutions of dfMRI....................................................23
3.8.1 Subject Specific BOLD HRF ....................................................................23
3.8.2 Regional Brain Differences .......................................................................24
3.8.3 Sub-threshold activations ..........................................................................24
3.8.4 Sub-regional (Only Capillaries and Venules ) Deconvolution..................24
Chapter 4 Methods of Dynamic fMRI.......................................................................25
4.1 Model-based Detection....................................................................................25
4.2 Independent Assumption in the Detection Process .........................................26
4.3 Fine Resolution Modeling of the BOLD Responses .......................................26
4.4 Down Sampling: Right Place and Right Time ...............................................26
4.5 Multi-Models for Timing Differences in Fine Resolution ..............................27
4.6 Multiple Detection per Voxel..........................................................................28
4.6.1 Multiple Detection Method .......................................................................28
4.6.2 Matrix Formulation for Multiple Detection Method.................................28
4.6.3 Logic Test if Needed .................................................................................30
4.7 Best to Display Results in a Movie or Time Frames.......................................30
4.8 Fine Tuning by Deconvolution at Each ROI for Regional HRF .....................30
4.9 Flow Chart of dfMRI Methodology ................................................................32
Chapter 5 Simulation: Limit of the Temporal Resolution.........................................33
5.1 Monte Carlo Simulation ..................................................................................34
5.1.1 Simulation Diagram...................................................................................36
5.2 Simulation Results...........................................................................................37
5.3 The Temporal Limit ........................................................................................39
Chapter 6 Experimental Validation: Temporal Resolution.......................................40
6.1 A Real Human Experiment to Validate the Simulation ..................................42
6.2 Consistent Result (vs. Simulation) ..................................................................43
6.3 The Case Studies of dfMRI .............................................................................45
6.4 Movie Frames: The dfMRI Results.................................................................46
6.5 Analyzing the dfMRI Results..........................................................................48
6.6 Directly Estimating Temporal Resolution from 2-Event Studies....................51
6.6.1 Visio-motor Tasks: The Method & Design ...............................................51
6.6.2 Data Processing: dfMRI Processing..........................................................52
6.6.3 Results and Discussion..............................................................................53
Chapter 7 Applications of dfMRI..............................................................................55
7.1 Motor Task in Patients.....................................................................................55
7.2 Picture Naming Task to Identify Language Areas ..........................................56
7.3 Artifact Identification and Reduction..............................................................57
7.4 Working Memory ............................................................................................57
7.5 The Software at Work: Step by Step ...............................................................59
v
Chapter 8 Neuromagnetic Signal Detection with MRI .............................................67
8.1 Not a Blood Related Functional Study............................................................68
8.2 B
0
and Working Neurons.................................................................................68
8.3 Fine Resolution Modeling and Aliasing..........................................................71
8.4 Experiment Design and Formulae: Exploiting Aliasing .................................74
8.5 Slice Timing Correction: A Special Case........................................................75
8.6 The Experiment ...............................................................................................76
8.6.1 Results .......................................................................................................78
8.6.2 Discussion..................................................................................................82
8.7 Spatiotemporal Dynamics Once Again ...........................................................83
Chapter 9 Neuromagnetic Signal Detection Difficulties...........................................84
9.1 Noise I: Low Frequency ..................................................................................84
9.2 Noise II: High Frequency ................................................................................85
9.3 Sampling and Task/Stimulus Rates.................................................................85
9.4 The Possible Solutions ....................................................................................85
9.4.1 Design........................................................................................................85
9.4.2 Filters.........................................................................................................86
Chapter 10 Neuromagnetic fMRI: An Optimal Design ............................................87
10.1 Choosing Aliasing Frequency of the Signal: The Design .............................88
10.2 Designing A Digital Zero-Phase Bandpass Filter .........................................89
10.3 The Experiment .............................................................................................89
10.4 Data Processing Methods ..............................................................................90
10.5 Results and Discussion..................................................................................90
Chapter 11 Future Work and Contributions ..............................................................93
11.1 Original Contributions...................................................................................94
11.1.1 Accurate Slice Timing Correction...........................................................94
11.1.2 Multiple Reference Functions to Address Variability.............................94
11.1.3 Making Complicated Computation Simple .............................................95
11.1.4 Spatio-temporal fMRI Visualization .......................................................95
11.1.5 Addressing Difficulties in Certain fMRI Experiments............................95
11.1.6 Detecting Neuro-magnetic Signals Directly with MRI ...........................95
References ....................................................................................................................96
vi
List of Tables
Table 6.1. Comparison of experimentally determined standard
deviation σ(exp) in ms, to that obtained from the
Monte-Carlo simulation σ(sim). …………………………………..… 45
Table 6.2. The difference of designed differential timing and
detected differential timing in the visio-motor
event-related experiment. ……………………………………………. 54
Table 8.1. Sum of t-score for activated voxels lying within
the visual cortex at p < .05 …………………………... ……………... 81
Table 10.1. Examples of design parameters for neuromagnetic
signal detection with MRI. …............................................................... 88
vii
List of Figures
Figure 2.1. Pulse sequence diagram for Echo Planar Imaging (EPI). ..………... 7
Figure 2.2. A model of the Hemodynamic Resonse Function (HRF)
With its peak around 4-5 seconds. ………………………..……….. 11
Figure 2.3. The HRF model plot according to Eq. 2.5. ……………………..…. 13
Figure 2.4. Top row shows the convolution between a stimulus train and
HRF to generate a model to be used before sampling at slice
timing to obtain a discrete reference function. …………………….. 15
Figure. 3.1. Neurons are located in a capillary bed between
an artery and a vein. …………………………………………..…… 18
Figure 4.1. Top two rows are the two slices of dfMRI movie frames.
The two bottom rows are the continuation from the top 2 rows. ….. 31
Figure 4.2. The block diagram of the dfMRI methodology. ………..…………. 32
Figure 5.1. (left) Temporal sampling of the hemodynamic response
function by combining samples from multiple trials.
As an example, the time-points when slice 1 is acquired
are shown by open circles. As shown at the right, the
interval between successive temporal samples is under 200 ms
for almost all samples and under 100 ms for a majority. .…………. 33
Figure 5.2. The block diagram of the simulation procedure. ………………… 36
Figure 5.3. (top 3 rows) Experimental time-courses from three pixels
activated in the finger movement protocol. The SNR in
these time courses are (top-to-bottom) 4.26, 1.84 and 1.25
respectively. The bottom row shows a simulated time-course
at SNR = 4.0. ……………………………………………………… 38
Figure 5.4. (Left) A plot of the standard deviation (std) in ms from
10,000 trials of a Monte-Carlo study, indicating the
uncertainty in determining the temporal delay of a pixel
as a function of SNR in its time-course.
(Right) Plots from simulations with added physiological
noise (1, 3, 5, and 7 times of the estimated noise power
from the real data) to the simulated signal. …………..……………. 38
viii
Figure 6.1. An example of four reference functions with Δt = 400 ms.
In the actual studies, Δt = 200 and 100 ms were used.
Only the first four lobes of the 19-lobe reference functions
are shown here for clarity. Notice that in addition to the
shift in the peak position, the shape of the reference
function also changes with Δt. …………………………………….. 42
Figure 6.2. The block diagram of the approach to validate the temporal
resolution. …………………………………………………………..44
Figure 6.3. Selected frames 400 ms apart showing the spatiotemporal
relationships among pixels activated within the SMA, PMA
and SA during the finger movement experiment for subject2.
The subject moved her left finger. ……….……………………….. 47
Figure 6.4. (a). Clusters in the SMA and plots of corresponding
activation metric 'A' defined in the text. ……………….………….. 49
(b). Clusters in the PMA and plots of corresponding
activation metric ‘A’ defined in the text. …...…………...………… 50
(c). Clusters in the SA and plots of corresponding
activation metric 'A' defined in the text. …………….……………. 50
Figure 6.5. The design diagram of the experiment with visual and motor
Tasks to estimate fMRI temporal resolution by using dfMRI. ……. 52
Figure 6.6. Selected frames 400 ms apart (subject # 2) showing the
spatiotemporal relationships among voxels in the visio-motor
experiment. ……………………………………………………….. 53
.
Figure 7.1. A motor study in a brain tumor patient.
Top row - activations associated with right hand movement.
Bottom row - activations associated with left hand movement. .…. 56
Figure 7.2. A study to identify speech areas by naming pictures. …………….. 57
Figure 7.3. In a working memory study, artifacts due to subject’s verbal
responses were detected at 7 seconds by a reference function
designed to detect prompt effects, while the working memory
activations were detected 4 seconds later. ……………………….... 58
ix
Figure 7.4. (a). Running the software:
making multiple reference functions. …………………...………… 61
(b). Running the software:
importing EPI data time-series. ……………...…………………….. 61
(c). Running the software:
importing T1 anatomical images. ……………………………...….. 62
(d). Running the software:
registering T1 images to EPI time-series. ……………...………….. 62
(e). Running the software:
making a brain mask. ……………………………………...………. 63
(f). Running the software:
running dfMRI engine. ……………………………...……………...63
(g). Running the software:
resultant timeframes 1-4. …………………………………...……... 64
(h). Running the software:
resultant timeframes 5-8. ………………………...………………... 64
(i). Running the software:
resultant timeframes 9-12. ……………………………...…………. 65
(j). Running the software:
resultant timeframes 13-16. …………………...…………………... 65
(k). Running the software:
resultant timeframes 17-20. ……………………………...………... 66
(l). Running the software:
optimal views of sensorimotor activation for
(top 2 rows) right hand movement and
(bottom 2 rows) left hand movement on 10 slices. ……...………… 66
Figure 8.1. Representative EEG using the same 4 Hz visual stimulus.
Two main peaks can be observed, a sharp peak at 50 ms
and a broad peak at 125-175 ms. …...……………………………... 71
Figure 8.2 (a). A conceptual depiction of how a gradual shift is
obtained between the stimulus presentation and slice
acquisition (or slice sampling) time. ……………………………..... 73
(b). Formation of reference functions by a model of the
expected neural current or neuromagnetic signal (black line)
and its sampling by the two slices (o and *). ……………………… 73
Figure 8.3. An example of conventional event-related BOLD imaging
for the flash checker board stimulus with TI = 20 s at
p <0 .001 for the same two subject shown in Fig. 8.4.
Strong activation is in the macrovasculature as expected. ….…….. 79
x
Figure 8.4. Images representing activity detected as a function of
latency in equivalent steps of 25 ms (mapped to 275 ms). ………. 80
Figure 8.5. Plot of t-score averaged over all scans. They are consistent
with EEG latencies. ……………………………………………….. 81
Figure 10.1. Power spectrum of a voxel time-course showing the aliasing
signal of interest peak and physiological signal peaks. …………… 87
Figure 10.2. Detected neuromagnetic signals at phases: 20, 40, …, 180 ms
on 2 selected subjects (S1-S2), corrected p <0.001.
S1 and S2 show activity superimposed echo planar images. ...……. 91
Figure 10.3. Plot of t-score in the visual cortex ROI of 5 subjects. ……..……… 91
xi
Abbreviations
3T 3 Tesla
ASL arterial spin labeling
BOLD blood oxygenation level dependent
CBF cerebral blood flow
CBV cerebral blood volume
CSF cerebrospinal fluid
EEG electroencephalography
EPI echo-planar imaging
ERP event-related potential
s
f sampling frequency
GLM general linear model
GWN Gaussian white noise
HRF hemodynamic response function
IRF impulse response function
ISI inter stimulus interval
LUT look up table
MEG megnetoencephalography
MRI magnetic resonance imaging
fMRI functional Magnetic Resonance Imaging
dfMRI dynamic fMRI
nfMRI neuromagnetic fMRI
xii
MRS magnetic resonance spectroscopy
NMR nuclear magnetic resonance
PET positron emission topography
PMA primary motor area
psd power spectral density
RF radiofrequency
SA sensory area
SMA supplementary motor area
SNR signal to noise ratio
SPM Statistical Parametric Map
std standard deviation
T signal period
1
T longitudinal (spin-lattice) relaxation time constant
2
T spin-spin relaxation time constant
2
* T spin-spin relaxation time constant including magnetic field
inhomogeneity
TA acquisition time
TE echo time
TI stimulus presentation interval (= ISI)
T
R
repetition time
s
T sampling interval
VOI volume of interest
xiii
Abstract
An important challenge in functional magnetic resonance imaging (fMRI) is to
achieve the most spatially accurate results, i.e., to localize activation as close as
possible to the actual site of working neurons. Because fMRI detection methods are
based on the blood oxygen level dependent (BOLD) property of brain microvascular
and venous systems, any model-based fMRI detection method will provide inaccurate
results if the model is biased to the BOLD response of a draining venous system, e.g.,
a vein. Here, a novel detection method is proposed to achieve accurate localization by
measuring the spatiotemporal dynamics of brain activations in fMRI to separate
microvasculature from venous activations. The method utilizes precisely down-
sampled multiple models (with positive and negative delays). Thus conventional 3D
fMRI results are expanded into 4D (with an additional temporal dimension). Example
fMRI studies underline the usefulness of the approach.
Even though early BOLD activations may be very close to the working
neurons, they are not guaranteed to be at the actual site of working neurons. Because
the BOLD effect is a slow process (a few seconds to its peak), one cannot use it to
detect fast electrical signals produced by working neurons. The proposed method is
extended to directly detect fast electrical signals of neurons by a design that converts
millisecond latencies of the electric fields to relatively slow latencies through aliasing.
The dephasing of the MRI magnetic field caused by the neural currents can then be
detected at the relatively slow aliased frequencies by using the proposed method to
measure spatiotemporal dynamics.
xiv
Keywords: Spatiotemporal Dynamics, fMRI, dfMRI, MRI, Direct Neuromagnetic
Signal Detection, nfMRI, Neuroelectric Signal, Magnetic Source, msMRI, Aliasing.
1
Chapter 1 Introduction
After the discovery of functional Magnetic Resonance Imaging (fMRI) in
1991, the studies of brain functions have never been the same. Non-invasively, fMRI
can be conducted without any contrast agent. Moreover, there is no known hazard for
repeating fMRI experiments on a human. It has been adopted very fast not only in
pre-neurosurgical planning but also in psychological and neuroscientific research.
FMRI is generally acknowledged as a tool to obtain information-rich data. Standard
fMRI is based on detecting Blood Oxygenation Level Dependent (BOLD) signals.
Even though fMRI seems to be a mature field of research, there are still many
challenges to be solved. A major one is that detected fMRI activations may not be at
the actual sites of active neurons. Another is that the sampling rate of fMRI is
relatively slow, resulting in aliasing artifacts. In addition, standard fMRI does not
demonstrate spatiotemporal relationship among activated brain regions.
In this dissertation, a new methodology is proposed to achieve better fMRI
spatiotemporal resolution and specificity, and to detect and display spatiotemporal
relationship of brain activations beyond standard fMRI. Also, the method is extended
to directly detect neuromagnetic signals with MRI.
After a brief overview on MRI and fMRI basics in Chapter 2, details of the
new purposed methodology including a Monte Carlo simulation, its validation, and
several applications are presented in Chapter 3-7. In chapters 8, 9 and 10, the method
2
is proposed to directly detect neuronal activity by a special design that relies on
exploiting aliasing signal properties. The final chapter is the future work.
3
Chapter 2 Basics
The essential basics of magnetic resonance imaging, physiology, mathematics,
and statistics are presented in this chapter to build up the understanding of functional
magnetic resonance imaging. Subsequently, additional fundamentals of
spatiotemporal dynamic functional magnetic resonance imaging are presented in
chapter 3, and the key concepts of direct neuromagnetic magnetic resonance imaging
are presented in chapter 8.
2.1 Magnetic Resonance Imaging (MRI)
In 1946, Bloch et al. at Stanford and Purcell et al. at Harvard independently
demonstrated nuclear magnetic resonance (NMR) phenomenon. Thereafter, NMR
spectroscopy has become an important analytical technique in chemistry. Because the
NMR frequency of a nucleus is sensitive and specific to its local environment, specific
spectra and images of a bulk matter can be reconstructed.
Lauterbur could reconstruct the first NMR images of a live animal in 1974, one
year after his group as well as Mansfield and Grannell independently published
procedures for NMR imaging [33, 38]. Because internal organs could be non-
invasively imaged, the medical profession rapidly adopted NMR imaging as a medical
investigation tool. Around 1980, the ‘nuclear’ term was eliminated from nuclear
magnetic resonance imaging in order to not alarm patients. (Actually, there is no
radioactive material in NMR imaging.) Consequently, nuclear magnetic resonance
imaging became magnetic resonance imaging (MRI).
4
2.1.1 Nuclear Magnetism, Spin, and Resonance
An atom consists of a nucleus and an electron cloud as its external shell. The
nucleus has a positive electric charge and is spinning. Therefore, an electric current is
produced. From Faraday’s law, the spinning nucleus (thus also called a ‘spin’)
generates a magnetic field because of the electric current. However, it is necessary
that there be an odd number of nucleons (protons or neutrons, but preferably protons)
to produce a net magnetic moment as even number of nucleons spin in opposite
directions and cancel the magnetization vector.
In MRI, the nucleus of a hydrogen atom (a proton) is the spin of interest,
because it has the highest gyromagnetic ratio and sensitivity of all diagnostically
relevant nuclei, and there are abundant water molecules and organic compounds that
contain hydrogen in the body.
If spins are placed in a magnetic field, they will align their magnetic fields with
the applied magnetic field in parallel (spin-up) or anti-parallel (spin-down). The spin-
down state has higher energy than the spin-up state. However, the spins still rotate
about the axis parallel to the external magnetic field at an angular frequency,
00
B ω γ = (2.1)
where γ is a constant (i.e., a gyromagnetic ratio) and
0
B is the applied magnetic field.
The angular frequency
0
ω (i.e., Larmor frequency) is the particular (resonance)
frequency of a radiofrequency (RF) irradiation causing the spins to jump between the
two energy states. (Hence, the ‘nuclear magnetic resonance’ term was used.)
5
2.1.2 Nuclear Magnetic Resonance Excitation
The local magnetization of spins within a voxel can be represented as a net
resultant vector of spin magnetization, M. Therefore, at equilibrium, M is parallel to
0
B (an applied magnetic field). To simplify the system, both M and
0
B can be
considered as if they are on a rotating Cartesian coordinate system which rotates at
0
ω , and where
0
B is parallel to its z-axis (longitudinal).
If the spins are perturbed by an RF irradiation (tuned to the resonance
frequency of the spins by an RF pulse duration and intensity) to generate a perturbing
magnetic field
1
B in the x-y (transverse) plane of the rotating coordinate, the net
magnetization vector M will flip from the z-axis toward the x-y plane to a certain
degree, e.g. 45°, 90°, 180°, depending on the strength and duration of the RF pulse.
2.1.3 T
1
and T
2
Relaxation Time
After the spins are excited, they will relax and return to their original
equilibrium states. The spins exchange their energy with the immediate environment
(lattice) during the relaxation, governed by
1
T, the longitudinal (spin-lattice)
relaxation time. The spins and nearby spins also begin to dephase because of nuclear
interactions. The spin-spin relaxation time (
2
T ) governs these processes in the x-y
(transverse) plane magnetization (so does
2
* T , the decay time constant, to be
explained later in 2.1.6).
6
2.1.4 Bloch Equations
The relationship between spin magnetization and magnetic field after an
excitation can be described with a set of equations (i.e., Bloch equations) as below:
2
() /
xxx
M MB MT γ=× −
&
(2.2)
2
() /
yyy
M MB M T γ=× −
&
(2.3)
01
()( )/
xzz
M MB M M T γ=× − −
&
(2.4)
where M is the magnetization vector; M
x
, M
y
, and M
z
are scalar (length) components
of the magnetization along the x-, y-, and z-axis, respectively; M
0
is the original
length of the magnetization (along the z-axis) before the excitation.
2.1.5 Spectroscopy and Imaging
Because of the specificity of the Larmor frequency (Eq. 2.1) and the physics of
the spins (Eqs. 2.2, 2.3, and 2.4), many protocols were designed to detect distinct
responses from an object. The protocols include the applications of magnetic
gradients and RF pulses (i.e., a slice selection, phase encoding, and frequency
encoding). There are 2 major NMR techniques: magnetic resonance spectroscopy
(MRS) for detecting frequency responses from a location in the object, and magnetic
resonance imaging (MRI) for imaging inside the object.
2.1.6 Decay Time constant, T
2
*
Macroscopically (in a voxel level), spin dephasing is also affected by nearby
spins that experience slightly different magnetic fields due to magnetic inhomogeneity
7
G
z
G
y
G
x
RF
G
z
G
y
G
x
RF
Fig. 2.1. Pulse sequence diagram for Echo Planar Imaging (EPI).
(because of imperfect instrumentations or certain physiologic changes). Thus, the
NMR signals diminish (due to dephasing) with another decay time constant
2
* T , in
addition to the standard spin-spin relaxation time
2
T . A
2
* T image can be achieved
by MRI techniques such as a gradient-echo sequencing and Echo-Planar Imaging
(EPI). Similar to
2
T images,
2
* T images provide a possibility to image functioning
brain (2.2.1) and to directly image working neuron clusters (chapter 8).
2.1.7 Echo Planar Imaging (EPI)
Echo Planar Imaging is a very fast MRI technique. Therefore, it is very useful
in imaging moving organs such as the heart, and brain activation. The technique
demands MRI hardware that can rapidly flip magnetic gradients. Fig. 2.1 shows the
pulse sequence of EPI.
8
The gradient along the z-axis (G
z
) is applied to select an anatomical slice
location. The applied RF pulse flips the spin magnetization onto the xy-plane so that
NMR signals are generated and then can be recorded. The flipping gradient along the
x-axis (G
x
) causes the regrowth of the NMR signal (called Gradient Echo). This
technique enables the encoding and recording the signal along the x-axis before the
spin-spin relaxation (
2
* T ) attenuates the transverse (xy-plane) magnetization. In
applying the very brief gradient (blip) along the y-axis at every cycle of gradient
flipping along the x-axis, the encoding for the y-direction is achievable.
2.2 Functional Magnetic Resonance Imaging (fMRI)
Being able to scan a subject repeatedly and rapidly, MRI has been applied to
visualize and analyze a functioning brain (i.e., functional magnetic resonance imaging,
fMRI). In general, an fMRI result is in the form of a statistical brain activation map
superimposed on a standard structural MRI image.
An fMRI requires a subject (animal or human) to perform a task or to sense a
stimulus while the subject’s brain is scanned many times. The resultant set of images
is then statistically analyzed to detect changes correlated to the task or the stimulus.
During the initial development of fMRI [8], MRI contrast agents such as gadolinium-
based compounds were used as injected exogenous paramagnetic MRI contrast agents
to detect changes in regional cerebral blood flow and/or cerebral blood volume caused
by active neurons.
9
2.2.1 Brain Physiology: Regional Vascular Responses
Working neurons require extra oxygen, which is transported by hemoglobin
molecules in red blood cells via the vascular system. To meet the increased oxygen
demand, within 1 second [58], the nearby arteriole smooth muscles relax, causing
arteriole dilatation, i.e., reducing intra-arteriole resistance. Consequently, regional
cerebral blood flow (CBF) is increased, bringing more oxyhemoglobin to the region.
At the same time, the local cerebral blood volume (CBV) is also increased in the
microvasculatures (arterioles, capillaries, and venules) but at a slower rate [36, 46]
compared to that of CBF. For prolonged activation, the increased CBF and CBV will
reach a steady state.
After the neurons stop being active, there is no demand for extra oxygen for
that region. The arterioles stop being dilated, causing the reduction of regional CBF
and CBV. Again, the CBV response is slower than that of CBF.
Approximately, arteriole-capillary-venule blood transit time is 2 seconds [35].
In general, this time is the shortest delay between the onset of the stimulus and the
detected fMRI signal from venules due to the stimulus.
2.2.2 Blood Oxygen Level Dependent (BOLD) Signal
After releasing 2 oxygen molecules, a diamagnetic oxyhemoglobin converts to
a paramagnetic deoxyhemoglobin, which can be used as an endogenous MRI contrast
agent. The action occurs in capillaries and it is very fast, causing an increase in
deoxyhemoglobin concentration in venules. Therefore, the MRI
2
* T signal from
microvasculatures adjacent to functioning neurons becomes smaller than that in a
10
normal (inactive) state, causing the initial dip of the
2
* T signal [3, 26, 61, 62, 63, 64].
Later, regional CBF and CBV are increased for the needed oxygen as mentioned
earlier. However, they overcompensate for the need for oxygen. Thus, in contrast to
the initial dip phenomenon, the deoxyhemoglobin in the region is reduced by the over
supply of oxygenated blood, causing an increase in the
2
* T signals from activated
brain regions.
In general, the positive
2
* T signal can be detected at 2-2.5 seconds after the
stimulus, and then it reaches its peak at around 4-5 seconds. However, if the stimulus
or the task is continued, the relatively diluted deoxyhemoglobin state (positive signal)
will reach a plateau at 5-6 seconds.
After the stimulus stops, there is no need for extra oxygen. As CBF and CBV
decrease, the
2
* T signal decreases according to the deoxyhemoglobin concentration
turning back to normal from the relatively diluted state. However, CBF decreases
relatively faster than CBV (mentioned in 2.2.1). The incoming oxyhemoglobin via
CBF is therefore relatively diluted, i.e., the deoxyhemoglobin is relatively
concentrated. Hence, before returning to the normal state, the
2
* T signal is
temporarily reduced than the normal state, causing a signal undershoot [31].
In summary, regional brain activity causes
2
* T signal changes. The signals
depend on the level of blood oxygen. Therefore they are called blood oxygen level
dependent (BOLD) signals [42]. In other words, the BOLD signals depend on a
hemodynamic response to regional working neurons. If a stimulus is brief (less than
1-2 seconds), the response can be considered as an impulse response function (IRF).
11
In fMRI, the IRF is called a hemodynamic response function (HRF). In general, fMRI
is based on detecting the BOLD signals that are detectable in capillaries, venules, and
vein systems. In practice, HRF is estimated by a deconvolution process by
considering the brain hemodynamic responses as a linear time invariant (LTI) system
[10].
Fig. 2.2 shows a model of the BOLD signal (HRF) with the initial dip, the
main signal, and the undershoot.
For further details, Buxton et al. purposed a complex physical and physiologic
model (i.e., the Balloon model) for the BOLD phenomenon [9]. In addition, Friston et
al. estimated physiologically plausible Volterra kernels for the Balloon model [17].
2.2.3 Arterial Spin Labeling (ASL)
Arterial spin labeling (ASL) is an MRI technique, which involves applying
RF pulses at a carotid artery to label spins in the blood water, before flowing into the
Figure 2.2. A model of the Hemodynamic Resonse Function (HRF) with its peak around 4-5
seconds.
12
brain [12]. Thus, ASL images can be processed to estimate and display CBF. Even
though ASL fMRI can estimate the relevant CBF at arterioles and capillaries adjacent
to working neurons, it is not performed as routinely as BOLD fMRI, because not only
the ASL signal is smaller than the BOLD signal but also the experiment setup,
instruments, and experimental design are more complicated. Moreover, the method is
time consuming and insensitive if the acquired volume covers many slices.
2.3 Detection Methods in fMRI
To detect fMRI signals (in order to obtain an activation map), many methods
can be applied. A simple one is an image subtraction, by which a baseline (no
stimulus) MRI image is subtracted from an active (functioning) brain image. The
resultant image is an activation map of positive BOLD signals.
However, fMRI signals are small and noisy. (In general, the signal is in 1 – 5
% range above baseline.) To achieve sufficient statistical power, a subject can be
scanned many times or the experiment can be repeated for a number of times. Then,
averaging and hypothesis testing (e.g., t-test) are applied to obtain a resultant
activation map.
As mentioned previously (2.2.1, 2.2.2), the brain response is a complex
process. However, if the system is considered as a linear time invariant (LTI) or a
non-linear system, a system black box can be represented by an IRF (HRF in BOLD
fMRI). Then signal processing techniques can be applied to the fMRI detection.
13
Figure 2.3. The HRF model plot according to Eq. 2.5.
2.3.1 Hemodynamic Response Function (HRF) Models
In signal processing, if the IRF of a system is known, the signal detection
becomes simple. Based on fMRI experiments and deconvolution, many parametric
BOLD HRF models were proposed (e.g., Cohen [10]; Friston et al. [14, 16]; Goutte et
al. [21]; Glover [20]; Lange and Zeger [32]). The simple Cohen’s model [10] is
8.6 / 0.547
()
t
ht t e
−
= (2.5)
where h(t) is the HRF, a function of time t. In Fig. 2.3, the relation of HRF BOLD
signal vs. time is plotted. The more sophisticated Friston et al. model [14, 16], which
can also simulate an undershoot, is
515
()
(6) (16)
tt
te t e
ht
− −
=−
ΓΓ
(2.6)
where
0
() ( 1)
t
tedt ττ
∞
−
Γ= −
∫
(2.7)
14
2.3.2 Correlation Coefficients
The fMRI intensity sequence (in time) of a particular voxel can be considered
as a signal vector (i.e., either a time-series or a time-course). To detect BOLD fMRI
signals, Bandettini et al. [5] computed the correlation coefficient of every time-course
versus a reference function, and then searched for the pixel correlation coefficient that
was higher than a certain threshold (which can be converted to a p-value for a
statistical significance). An activation map consists of all voxels whose correlation
coefficients surpass the threshold. The reference function can simply be the sequence
coding of fMRI experimental conditions and control conditions or a realistic model by
convolving the sequence coding with a HRF. Fig. 2.4 shows examples of reference
function. The correlation coefficient method is sensitive but can be applied only to a
2-condition (or condition-control) experiment. The correlation coefficient (cc)
formula is
1
22
11
()( )
() ( )
N
iy i x
n
NN
iy i x
nn
yx
cc
yx
μμ
μμ
=
==
−−
=
−−
∑
∑∑
(2.8)
where N is the total acquired fMRI images or the length of a pixel intensity time-
series, y
i
is the observed pixel intensity at i-th image (time), x
i
is the element of the
reference function at time i, μ
y
and μ
x
are the means of all y
i
and x
i
, respectively.
15
2.3.3 General Linear Model (GLM)
For a multi-condition fMRI experiment, a general linear model, GLM, is
suitable for the fMRI detection. Because the conditions can be encoded as regressors,
a GLM for a pixel intensity time-series can be expressed as
y β ε = + X (2.9)
where y is the observed pixel intensity vector, X is a design matrix (the collection of
regressors or reference functions), β is a to-be-estimated parameter vector, ε is an error
vector such that its elements are independent and identically distributed normal
random variables with zero mean and variance σ
2
.
In a null hypothesis testing, whether a voxel time-series is activated, a t-score
can be computed as
21
()
T
TT
c
t
cX X c
β
σ
−
= (2.10)
where t is a t-score to be compared with a Student’s t-distribution, c is a contrast
vector for the null hypothesis testing.
⊗
Figure 2.4. Top row shows the convolution between a stimulus train and HRF to generate a
model to be used before sampling at slice timing to obtain a discrete reference function.
16
2.3.4 Statistical Parametric Map (SPM)
The Wellcome Department of Cognitive Neurology Functional Imaging
Laboratory (http://www.fil.ion.ucl.ac.uk) implemented the GLM in a free software
(i.e., Statistical Parametric Map, SPM) to process positron emission topography (PET)
and fMRI data [15]. The SPM software has become a standard fMRI processing
software, and it also includes various statistical signal and image processing tools
(such as a spatial transformation and filters). In late 2005, the latest version (SPM5)
was released. It adds a capability to process and model electroencephalography (EEG)
data and an improved spatial image data pre-processing.
17
Chapter 3 Dynamics of fMRI Signals
Standard fMRI can be compared with taking a photographic snap shot of the
working brain to obtain an activation map. Thus, only static or time-averaged spatial
information and relative fMRI signal strength or statistical map with a certain
confidence can be observed. Thomas Edison invented a movie technology by
projecting a series of snap shot photographs. Unlike a photograph, a movie can
present the sense of spatiotemporal dynamics. In these few chapters, a new method is
presented to demonstrate how to achieve spatiotemporal dynamic fMRI (dfMRI)
beyond standard fMRI. As a movie does to a photograph, dfMRI adds spatiotemporal
dynamics to standard fMRI.
3.1 Neuro–vascular Anatomical Background
Arranged in layers, neurons are in the gray matter of a brain inside a skull.
Neuronal dendrite and axon carry electrical input and output signals, respectively.
Bundles of axons comprise most of the brain white matter.
Even though there is a special blood supply auto-regulation, brain vascular
organization is quite similar to other organs. Sequentially, blood flows through a big
artery, a small artery, an arteriole, a capillary, a venule, a small vein, and a vein.
However, there are venous sinuses collecting venous blood before draining it out of
the skull.
18
Figure. 3.1. Neurons are located in a capillary bed between an artery and a vein.
Located closer to capillaries than other vascular structures, the functioning
neuron receives nutrients and extra oxygen from the capillary blood, consequently,
giving the BOLD signal mentioned in 2.2.1 and 2.2.2.
Fig. 3.1 shows the spatial relationship of brain vascular system where neurons
are located in the gray matter, in a capillary bed between an artery and a vein.
3.2 Big Vein BOLD Contributions
In the venous drainage system, many capillaries drain into a venule, venules
drain into a vein, and then veins drain into a bigger vein. Because an MRI signal
depends on a measurement per a unit volume (voxel), therefore, down the drainage
stream and later in time, the BOLD signal becomes stronger [22, 54]. Consequently,
19
most fMRI detection methods, which are mainly based on signal changes, are more
sensitive to the BOLD signals at big veins than those at small veins, venules, or
capillaries.
However, for spatial specificity, the BOLD signals from bigger veins are
undesired, because the locations can be far from the working neurons. There are some
techniques to mitigate the big vein BOLD signals or to enhance the microvascular
(i.e., arteriole, capillary, and venule) BOLD signals, e.g. applying stronger magnetic
fields [13], masking fMRI by MRI venography [22], increasing fMRI detectability by
acquiring a lot of data [44], detecting the BOLD initial dip [61, 62, 63], combining
fMRI and perfusion MRI [47].
Because it is able to detect the BOLD spatiotemporal characteristics, the new
dfMRI method developed here enables us to visualize and detect not only the big vein
BOLD maps but also, perhaps, the early microvascular BOLD maps.
3.3 Multi-focal Brain Activations
Even in a simple task such as tapping a finger, many parts of the brain work in
concert. Typically, there are fMRI activations not only in the contralateral pre-central
and post-central gyri (predominantly in the central sulci due to the big vein BOLD
contributions) but also in the bilateral pre-motor cortex and the supplementary motor
area. From electroencephalography (EEG) studies, brain activities are also the
functions of space and time.
Better than the fMRI which displays only the spatial information, a dfMRI
movie can display the spatiotemporal relationships of brain activation.
20
3.4 Regional Vasculature Variation and BOLD Specificity
Macroscopically, the brain vascular organization seems to be the same in all
areas of the brain. Microscopically, it is nearly impossible that the blood vessel at a
certain level always has the same length, size, or amount. Consequently, the HRF or
BOLD response should be distinct in each brain region. In a simple motor study, the
BOLD responses of different motor areas were not the same [41]. However, based on
averaging, the estimated HRF (by a deconvolution process) of the motor area and the
visual cortex were not significantly different [10]. Actually, to be very accurate in
finding a BOLD HRF, a deconvolution should be performed with a limited spatial
extent, because many levels of the brain vasculature are in a large volume of interest
(VOI) and each level can have a different HRF. In addition, the average across several
VOIs from different brain regions will definitely provide an inaccurate HRF.
One of the advantages of dfMRI is the ability to display spatiotemporal
activation maps. Thus, VOIs can be selected as a function of time, then, HRFs can be
estimated accurately.
3.5 Spatiotemporal Resolution of fMRI
Requiring an MRI scanner, an experiment protocol, and a biological subject,
fMRI is the final product of multi-factorial considerations and processes. A standard
fMRI shows only a spatial result (i.e., an activation map), computed from the temporal
information (time-series or time-courses) for all scanned regions.
21
3.5.1 Spatial Resolution
Mainly, the MRI spatial resolution depends on the physical properties of the
MRI scanner (e.g., magnetic field strength, an RF coil), an acquisition time (TA), and
the field of view (FOV) [37]. Due to being computed from a sequence of images,
fMRI spatial resolution depends on a sampling period (repetition time, T
R
) and
intrinsic (anatomic and physiologic) BOLD hemodynamic response.
At the standard magnetic field (1.5 Tesla) for clinical uses, the fMRI spatial
resolution is in 1.5 – 3 mm range. Practically, the resolution of 3 mm is preferred in
order to achieve a high signal to noise ratio (SNR) of MRI and a high sampling
frequency (a short T
R
) for efficient fMRI data processing.
3.5.2 Temporal Resolution
Whereas a sampling period (T
R
) can be shorter than 100 ms, the fMRI
temporal resolution is in 100-400 ms range because of the SNR and the intrinsic
(anatomic and physiologic) BOLD hemodynamic response [54]. Unlike the resolution
in standard signal processing, the fMRI temporal resolution depends on the slow
BOLD signal change. Attempts have been made to improve the fMRI resolution to be
finer than the twice of sampling period by controlling stimulus timing of repeated
experiments [24, 28, 40, 43, 54].
22
3.6 BOLD fMRI Signal Flow
Because fMRI is the result of the relevant BOLD detection, dfMRI displays the
BOLD activation map as a function of time. Therefore, theoretically, in a certain brain
region, the detected BOLD signal begins at a capillary, then, flows to a venule, a vein,
a bigger vein, and a venous sinus finally. In addition, if the task involves many
regions of the brain, dfMRI will show temporal relations and evolutions among the
regions.
3.7 Advantages of Dynamic fMRI (dfMRI)
Adding another dimension (time) to fMRI, dfMRI benefits not only temporal
but also spatial information. Moreover, spatiotemporal brain networks can be
visualized and further analyzed.
3.7.1 Space: Capillaries, Venules, Veins, or a Sinus
Standard fMRI detects BOLD signals with a HRF model or a significant signal
chance. However, in dfMRI, early or late fMRI signals in the same brain region are
associated with smaller or bigger draining veins respectively.
3.7.2 Time: Fast or Slow Responses
It is possible that different brain regions response to a task at different rates
because of brain control mechanisms such as excitatory, inhibitory, and networks.
The variability of supporting vasculatures in different brain regions also contribute to
the different responses. Moreover, a complex task may involve sequential brain
23
activities. The dfMRI allows us not only to visualize but also to analyze
spatiotemporal responses.
3.7.3 Spatiotemporal Network
Even in performing a simple task such as moving a finger, the brain works in
concert involving timing and many regions. Thus, brain activation map will be more
informative if displayed as a function of time, not just in 3-D or 2-D as a photograph.
A spatiotemporal network of the brain for a task can be visualized by dfMRI.
3.8 Potential Problems and Solutions of dfMRI
In order to have dfMRI be very specific and accurate, many problems are still
needed to be solved. Even in standard fMRI, there are some controversies and
uncertainties [4, 23].
3.8.1 Subject Specific BOLD HRF
Even though there was little variability of the HRF within a subject, substantial
variability across subjects was observed [1]. Practically and generally, fMRI detection
is based on a HRF model (2.3.1). In the analysis of brain activity relationship,
hemodynamic deconvolution should be performed on every subject [19].
24
3.8.2 Regional Brain Differences
Because the HRF is also varied in different brain regions, one should always
be cautious about any interpretation of brain networking or dynamics. Moreover,
there might be feedback, sustained, and transient activities in the same brain region to
complicate brain activation interpretations.
3.8.3 Sub-threshold activations
Most of detection processes rely on thresholds. It is likely that there are some
sub-threshold activations in fMRI. Perhaps, the active regions are too small or too
transient to be detected. Noise is also the major cause of undetectable activations.
3.8.4 Sub-regional (Only Capillaries and Venules ) Deconvolution
To address the uncertainties, selective deconvolution should be performed on
fMRI signals only from capillaries or venules. Then the resultant HRF should be used
in a second round of dfMRI detection. An average of multiple resultant HRFs might
be useful in dfMRI or detecting brain networking in order not to be biased against a
certain region.
25
Chapter 4 Methods of Dynamic fMRI
In chapter 2.3, a couple of fMRI detection methods were explained. The
methodology of dfMRI is present here. The main concepts are: modeling the BOLD
responses in continuous time domain with many negative and positive latencies
(delays); sampling the models at actual MRI acquisition time for each slice (e.g., along
the z-axis); detecting the BOLD signals in every voxel; arranging the detected maps
according to latencies.
4.1 Model-based Detection
Generally, during an fMRI experiment, there is a task for a subject to execute
or a stimulus for the subject’s brain to perceive. To improve reliability and
detectability, the task or the stimulus is repeated many times. The same experiment
can be repeated for many sessions. In dfMRI, only one kind of task or stimulus
(condition) in an experiment is recommended in order to avoid the interaction or
relation between different conditions. (A multi-condition experiment will be
presented as a further application in Ch. 7.) Furthermore, event-related experimental
design is preferred as it has many repeated responses to be detected, even though the
detectability of event-related design is poorer than a conventional (box-car, i.e.,
stimulus / no stimulus) design. Finally, to avoid non-linear superimposition of the
BOLD response, an inter stimulus interval (ISI) should be wide enough, e.g., 18-20
seconds, in order to avoid the superimposition effects of hemodynamic responses
(detected as BOLD signals).
26
4.2 Independent Assumption in the Detection Process
Similar to a standard fMRI detection method, each voxel is assumed to be
independent (voxel-wise). However, in general, a resultant fMRI activation map
comprises clusters of voxels. Because dfMRI results are separated into many
activation maps according to different time latencies, a cluster might not be large at a
latency, but it might be much larger than that in a standard fMRI activation map if the
resultant maps of a few latencies are superimposed. A dfMRI result will demonstrate
that each voxel is actually not independent in time, even though a spatial
independence is assumed in the detection processes.
4.3 Fine Resolution Modeling of the BOLD Responses
Because the brain works in a continuous time domain, the BOLD responses are
modeled at a fine resolution (10-200 ms scale), even though a chosen MRI sampling
interval (T
R
) is in the level of seconds in order to achieve a high signal to noise ratio,
SNR. To generate a model, an impulse train corresponding to the experiment stimulus
timings is convolved with a HRF (Eq. 2.5 or 2.6).
4.4 Down Sampling: Right Place and Right Time
Practically, in fMRI, when an MRI scanner acquires a volume data, it is
programmed to scan slice by slice (e.g., along z-axis). Thus, each slice is scanned at a
different time. The sampling rate of the sequence of volume images is 1/
R
T .
Approximately, the interval between consecutive scanned slices is equal to
R
T divided
by the number of slices in the scanned volume. To avoid interferences between
27
spatially consecutive slices, the order of scanned slices is usually interleaved, e.g., 1 3
5 7 9 2 4 6 8 10. If
R
T is large, the first and the last acquired slices in a MRI image
volume are sampled at nearly one
R
T apart.
To solve the difference in sampled slice timing problem, most fMRI detection
programs shift the data time-series of a pixel from a particular slice to a reference slice
by transforming the time-series to Fourier frequency domain, shifting in frequency
domain, then inversing the data back to a time-series. The process can be
implemented as a pre-processing step, i.e., slice timing correction [55, 56]. The
advantage of this method is that only one model is required for the fMRI detection
process.
For accuracy, instead of shifting the data, the HRF model may be sampled at
actual acquisition time for each slice [54]. There are as many sampled models as the
number of slices in the acquisition volume. Thus, there is no need to perform slice-
timing correction as a pre-processing step. Equivalently, shifting the fine-resolution
model by the time difference before down-sampling will achieve similar results [59].
4.5 Multi-Models for Timing Differences in Fine Resolution
To detect the latency (delay) differences of voxel time-series data, a BOLD
response is modeled in a fine resolution, shifted negatively and positively in time, and
then down sampled at the acquisition timing of a slice. In the previous section, there
are many HRF models due to slice timing. Now there are many models due to the
BOLD response delay for each acquired slice.
28
4.6 Multiple Detection per Voxel
Typically, there are many acquired slices in an fMRI study. For the ease and
convenience of notation, the proposed detection method is explained for each acquired
slice. The detection must be performed for every slice in the scanned volume.
4.6.1 Multiple Detection Method
For every voxel time-series data, a correlation coefficient between the data
time-series and every (sampled) BOLD response model, i.e., a reference function, is
computed using Eq. 2.8. Therefore, there are as many resultant correlation
coefficients as the number of latencies (delays) to be detected.
4.6.2 Matrix Formulation for Multiple Detection Method
For compactness and ease in programming, the detection algorithm can be
represented in a simple vector-matrix formulation of the form
()( )
22
T
T
XY
C
X Y
=
∑ ∑
(4.1)
where the every latency (total = L) reference function (zero mean) collection matrix
12
...
L
X xx x
⎡ ⎤
=
⎣ ⎦
uruuruu r
(4.2)
and the every pixel (total = M) time-series (zero mean) collection data matrix
12
...
M
Yyy y
⎡ ⎤
=
⎣ ⎦
u u ruuruuu r
(4.3)
Here, on matrices, the
∑
operation is a column-wise summation, while the square, the
square root, and the division operations are element-wise. Hence, the resultant matrix
29
C is a collection of correlation coefficients of every pixel time-series and every
modeled reference functions of all latencies.
The content of C matrix can be further masked (with a threshold) and/or
logical tested (Ch. 4.7), then displayed as movie frames to visualize spatiotemporal
dynamics of a brain performing a task.
With such a compact notation, implementation of the algorithm is
straightforward and not prone to programming errors. Moreover, data processing time
is drastically reduced if an optimized-matrix-operation program (such as MATLAB
[39]) is utilized.
In the slice-by-slice processing, Eq. 4.1 operation reduces multi-loop
operations in 4 dimensions (2 dimensions from a voxel location in a slice, 1 dimension
from the voxel time-series, and 1 dimension from latency/delay) into simple 2-
dimension matrix operations.
30
4.6.3 Logic Test if Needed
At a voxel, there are a series of resultant correlation coefficients corresponding
to the set of modeled reference functions. The reference functions with too early or
too late latencies will present a lower correlation coefficient value than those with the
most correct latency. The reference function associated with the maximum correlation
coefficient presents a relative activation latency,
argmax
ll
tc = (4.4)
where
l
c is a resultant correlation coefficient associated with a modeled reference
function in the C matrix (Eq. 4.1). Relative to the same reference function, each
voxel activation may be different or have the same value as
l
t . Thus
l
t can be used as
relative timing or latencies for spatially different voxel activations.
4.7 Best to Display Results in a Movie or Time Frames
As mentioned in Ch. 3.7 that dfMRI are 4D maps. Therefore, to appreciate the
results, the maps should be presented in a movie format or frames to demonstrate the
spatiotemporal dynamics of the brain activations. Fig. 4.1 shows the brain activations
of a simple motor task in movie frames [54].
4.8 Fine Tuning by Deconvolution at Each ROI for Regional HRF
The brain venous vasculature may be different in different brain regions.
Therefore, in comparing activation latencies (temporal information) between far apart
regions, only microvascular and venule latencies should be utilized to avoid the
venous structure uncertainties. In Ch. 3.8.4, the regional deconvolution of
31
microvascular and venule early responses is proposed. After early-response
deconvolution of different regions are performed, the resultant HRFs should be
averaged before conducting dfMRI detection again in order not to be biased to any
particular brain region. However, in practice, the fMRI time-series data are noisy.
Thus it is difficult to perform deconvolution except by estimating model parameters.
In general, the standard HRF models are sufficient for the detection processes and
dfMRI visulization.
Figure 4.1. Top two rows are the two slices of dfMRI movie frames. The two bottom rows are
the continuation from the top 2 rows [54].
32
4.9 Flow Chart of dfMRI Methodology
Fig. 4.2 shows the block diagram of the dfMRI methodology. Standard fMRI
processing steps such as image registration, spatial smoothing, and activation map
superimposition on anatomical images are not included in the diagram.
stimulus / task sequence (high resolution)
matrix of time courses
for a slice
HRF (high resolution)
⊗
high resolution model
multiple high resolution models
shift + nT; n = 1, 2 …, L
total (2L+1)(# of slices) reference
functions at effective TR resolution
matrix of reference functions
for a slice
volume image
sequence
multiple detection process
cc or GLM
cc or t-score matrix
for a slice
optional max cc or t-score
optional assigning latencies optional filtering out low cc
activation map
time frames
thresholding
stimulus / task sequence (high resolution)
matrix of time courses
for a slice
HRF (high resolution)
⊗
high resolution model
multiple high resolution models
shift + nT; n = 1, 2 …, L
total (2L+1)(# of slices) reference
functions at effective TR resolution
matrix of reference functions
for a slice
volume image
sequence
multiple detection process
cc or GLM
cc or t-score matrix
for a slice
optional max cc or t-score
optional assigning latencies optional filtering out low cc
activation map
time frames
thresholding
Figure 4.2. The block diagram of the dfMRI methodology.
33
Chapter 5 Simulation: Limit of the Temporal Resolution
To estimate the limit of fMRI temporal resolution, event-related (finger-
movement) fMRI data were simulated for the dfMRI detection. (Analysis using real
experimental data is in Ch. 6.) While the image sampling interval (T
R
) is equal to 1.2
seconds, the timing of the event was randomized to produce a mean of 15 seconds
with a standard deviation of 2 seconds, i.e., the inter-stimulus interval (ISI) was 15 +2
seconds. Because of the jittering in the ISI, the actual time at which the HRF (using
Eq. 2.5) was sampled (i.e., the slice acquisition time) was staggered from trial-to-trial
with respect to the HRF. Assuming a stationary HRF following each event, a plot of
the HRF and the different times at which it is sampled throughout the acquisition
sequence are presented in Fig. 5.1 [54]. It is apparent that almost the entire HRF is
sampled in less than 200 ms time-intervals with a significant sampling below 100 ms,
suggesting that a possible temporal resolution of 200 ms.
Figure 5.1. (left) Temporal sampling of the hemodynamic response function by combining
samples from multiple trials. As an example, the time-points when slice 1 is acquired are shown
by open circles. As shown at the right, the interval between successive temporal samples is under
200 ms for almost all samples and under 100 ms for a majority [54].
34
5.1 Monte Carlo Simulation
A simulation study was conducted to estimate the limit as a function of SNR in
the time-course of a pixel where the SNR was defined as the ratio of peak signal
height to the standard deviation in the baseline. The simulation relied on the protocol
used to acquire data for the finger flexion study described later in Ch 6. The
convolved HRF was shifted in steps of +100 ms and down-sampled at time-points
corresponding to the slice acquisition time to form a series of 61 modeled reference
functions per slice. The reference functions thus span a time-shift of +3 s with respect
to the zero-delay reference function, where the zero-delay reference function
represents no additional delay except that inherent to the peaking time of the HRF
with respect to the onset of the stimulus, to form a reference function matrix X (Eq.
4.1). These reference functions also represent noise-free time-courses at specific
delays.
It is well known that there are two main sources of noise in fMRI,
physiological noise and electronic noise, and the net noise is a quadrature summation
of the two. The physiological noise is a function of the precession frequency
0
ω ,
whereas the electronic noise is a function of
0
ω . The physiologic noise represents
fluctuations in the hemodynamic response of the brain from respiration and cardiac,
cerebrospinal fluid (CSF) and vascular pulsations. Any overlap in the physiologic
frequencies or their aliases with the stimulus presentation frequencies can lead to the
so-called “physiological artifacts” in fMRI [2, 27].
35
To achieve a SNR between 1 and 10, the effect of electronic noise was
simulated by adding random Gaussian white noise (GWN) to the noise-free convolved
HRF. Typically, the BOLD signal represents a 1-5% increase in signal intensity from
the baseline, and baseline fluctuations are approximately 0.5- 1%. Thus the SNR is in
the 1-10 range, where the higher SNR values are for signals originating from veins and
lower values from the microvasculature [18].
After adding electronic noise, the convolved HRF was down-sampled at a
time-point lying midway between two reference functions to obtain a time-course
simulating the worst-case scenario where the delay of a pixel might lie exactly
between two reference functions. Physiological noise was added to this time-course
by a sinusoidal function at the primary frequency of the stimulus presentation
sequence. The sinusoidal function was estimated from experimental data (in Ch. 6) by
comparing the power spectra of the time-courses of ~100 activated pixels lying within
the primary motor area, PMA, (referred to as set 1) to the power spectra of ~100
pixels lying within gray and white prefrontal brain regions where no activity was
expected during finger movements (set 2). It was found that for an average SNR =
2.0 in the time-courses of the PMA activated pixels, the ratio of the power between set
2 and set 1 pixels at the primary stimulus frequency ranged from 0.05 to 0.2.
Accordingly, the amplitude of the sinusoidal function was selected randomly such that
after repeated trials, the ratios of the sinusoidal power to that of the time-course (after
adding noise) at the primary stimulus frequency were distributed uniformly between
0.05 and 0.2. The phase of the sinusoidal function was also assigned randomly
36
between -π and π for each trial. For simulating time-courses at different SNRs, the
amplitude of the sinusoidal function was scaled appropriately to conform to the
selected SNR.
The process of adding random GWN and physiological noise was repeated
10,000 times to simulate 10,000 samples of time-courses per SNR. The temporal
delay for each sample was then computed from the series of reference functions and
the standard deviation in the 10,000 delay values was used as the estimate of the
expected uncertainty in the delay.
5.1.1 Simulation Diagram
The block diagram of the simulation procedure is shown in Fig. 5.2.
HRF*(task impulse train) = noise free activation model
(sampling interval, T
s
= 10 ms )
Delay for 50 ms
Shift the model from
–3000 to 3000 ms at 100 ms
step size Æ 61 delays
Down sampling the 61
delayed models at MRI
slice acquisition timing
61 reference functions
(T
s
= 1.2 s)
Noise free activated data
time course (T
s
= 1.2 s) GWN to simulate electronic
and measurement noise
at SNR = 1, 2, …, 10
Sinusoidal wave (stimulus
frequency) with random power
(0.05 Æ 0.2) and random phase
to simulate physiologic noise
Down sampling at MRI
slice acquisition timing
Synthetic time course
Compute cross correlation
with all 61 reference functions
61 correlation coefficients (cc)
Find maximum cc Assign delay Detected delay
x 10,000 times
+
+
Compute means & std’s
HRF*(task impulse train) = noise free activation model
(sampling interval, T
s
= 10 ms )
Delay for 50 ms
Shift the model from
–3000 to 3000 ms at 100 ms
step size Æ 61 delays
Down sampling the 61
delayed models at MRI
slice acquisition timing
61 reference functions
(T
s
= 1.2 s)
Noise free activated data
time course (T
s
= 1.2 s) GWN to simulate electronic
and measurement noise
at SNR = 1, 2, …, 10
Sinusoidal wave (stimulus
frequency) with random power
(0.05 Æ 0.2) and random phase
to simulate physiologic noise
Down sampling at MRI
slice acquisition timing
Synthetic time course
Compute cross correlation
with all 61 reference functions
61 correlation coefficients (cc)
Find maximum cc Assign delay Detected delay
x 10,000 times
+
+
Compute means & std’s
Figure 5.2. The block diagram of the simulation procedure.
37
5.2 Simulation Results
Examples of time-courses from three activated pixels in a typical human study,
where one pixel was within the PMA (SNR = 4.26), the second within the sensory
area, SA, (SNR = 1.84) and the third within the supplementary motor area, SMA,
(SNR = 1.25), are shown in the top three rows of Fig. 5.3 and a typical Monte-Carlo
simulated time-course at SNR = 4.0 is shown at the bottom. The simulated data
appears very similar to the experimental data.
The results of the Monte-Carlo study to estimate the uncertainty (i.e., the
standard deviation) in the temporal delay of a pixel are plotted in Fig. 5.4(Left) as a
function of the SNR in the pixel's time course. The results show that the uncertainty or
temporal resolution ranges from about 350 ms at SNR = 1 to about 50 ms at SNR =
10. The 50 ms limit is due to our assumption that the event occurs exactly mid-way
between two reference functions and represents a worst-case scenario.
Fig. 5.4(Right) shows the same std-SNR plot (without added simulated
physiological noise to the simulated signal) and the plots from the simulations with
added physiological noise (1, 3, 5, and 7 times of the estimated noise power from the
real data) to the simulated signal. At the physiological noise power level estimated
from the real data (+ PN x 1), the uncertainty of the temporal delay detection does not
degraded much. Moreover, the methodology is still robust at the 3 times physiological
noise power (+ PN x 3) , in which the detection uncertainty is increased by
approximately 50 ms.
38
Figure 5.4 (Left). A plot of the standard deviation (std) in ms from 10,000 trials of a Monte-
Carlo study, indicating the uncertainty in determining the temporal delay of a pixel as a function
of SNR in its time-course [54]. (Right). Plots from simulations with added physiological noise
(1, 3, 5, and 7 times of the estimated noise power from the real data) to the simulated signal.
Figure 5.3. (top 3 rows) Experimental time-courses from three pixels activated in the finger
movement protocol. The SNR in these time courses are (top-to-bottom) 4.26, 1.84 and 1.25
respectively. The bottom row shows a simulated time-course at SNR = 4.0 [54].
39
5.3 The Temporal Limit
In fMRI experiments, it is not difficult to achieve the SNR = 3-4. Thus the
temporal resolution in an event-related fMRI experiment can be as low as 100 ms by
the dfMRI detection method. In addition, the fMRI resolution can be improved by
designing the ISI and T
R
for a sub-T
R
resolution. For a standard fMRI experiment (a
box-car or epoch design), where there are not many epochs per experimental session,
we cannot expect such a fine temporal resolution. However, the dfMRI technique will
still be able to display the spatiotemporal information in the BOLD signal. (An
example of a picture naming study on a patient will be presented in Chapter 7.)
40
Chapter 6 Experimental Validation: Temporal Resolution
In Ch. 5, the simulation to estimate the temporal resolution of BOLD responses
was presented. To validate the results, human experiments were conducted. A
quantitative section is presented first, and dfMRI result is later presented with movie
frames.
Event-related fMRI studies were conducted during a finger flexion task for
three normal human volunteers. The subjects were instructed to flex either the right
index finger (for two right handed males) or the left index finger (for one left-handed
female) and place it on a response button fiber-optically coupled to a computer. Upon
hearing a verbal cue, the task was to immediately press the button (thereby
electronically marking the initiation time of the movement), extend the finger and
return to the original position with the flexed finger resting on the button. The
subjects were trained before entering the magnet to press the button, extend and return
to the original position in a smooth motion in approximately 700 ms. The timing of the
cue was randomized to produce a mean of 15 s with a standard deviation of 2 seconds,
i.e., the inter stimulus interval (ISI) was 15 +2 seconds. A GE 1.5 T LX 8.3 EPI
system, equipped with a quadrature head-coil and gradients of strength 25 mT/m and
125 T/m/s slew rate, was used to acquire time-series data with T
R
= 1.2 s (effective),
TE = 45 ms, 90 degree flip angle, 64x128 acquisition matrix, 20x40 cm
2
field-of-view,
and NEX (number of excitations) = 1. A total of 250 images were acquired per slice
from two 8 mm thick contiguous transaxial slices covering a substantial portion of the
41
motor regions. The total acquisition time was approximately 5 min and contained 19
cued flexions.
The fMRI data were processed by correlating the time-course of each pixel to a
series of reference functions. The stimulus presentation sequence was modeled in
steps of 700 ms wide rectangular pulses and convolved with a model of the HRF [10]
at a high temporal sampling rate of 10 ms and then down-sampled at time-points
corresponding to the actual image acquisition time per slice, yielding a different
reference function for each slice. We refer to this reference function as the zero-delay
reference function.
Tracking of the BOLD signal per pixel was achieved by correlating the time-
series of each pixel to a series of reference functions per slice. These reference
functions were generated by shifting the convolved response in steps of Δt before
down-sampling, where Δt was +100 ms for the simulation and validation study and
+200 ms for the pixel tracking study described below and in Ch.4, the dfMRI method.
(Note that this procedure is not equivalent to the commonly used procedure of shifting
a given reference function to model delay.) In this approach, not only the position of
peaks but the entire shape of the reference function is subject to change as a function
of Δt in accordance with the shape of the HRF. An example of a partial series of
reference functions composed of four reference functions with Δt = 400 ms is
presented in Fig. 6.1 [54].
42
Each reference function within a series of reference functions is a marker of a
specific time-delay, i.e., latency with respect to the stimulus onset time. The time-
series of every pixel in each slice was then correlated with the reference functions, and
through a logic test (Eq. 4.4), each pixel was assigned to that reference function which
produced the highest correlation coefficient. The latency of each pixel was thus
represented by the delay of the reference function producing the highest correlation
coefficient value.
6.1 A Real Human Experiment to Validate the Simulation
A human study was conducted to validate the Monte-Carlo simulations. Using
the MRI pulse sequence and protocol as described above, the subject performed a cued
flexion of the right index finger 19 times to complete one experiment in 5 min. Then
the experiment was repeated seven more times during the same MRI scanning session
Figure 6.1. An example of four reference functions with Δt = 400 ms. In the actual studies, Δt =
200 and 100 ms were used. Only the first four lobes of the 19-lobe reference functions are
shown here for clarity. Notice that in addition to the shift in the peak position, the shape of the
reference function also changes with Δt [54].
43
to obtain data from a total of eight trials. These fMRI data were processed as
mentioned earlier, using a series of reference functions spaced +100 ms apart to
determine the delay of each activated pixel with respect to the zero-delay reference
function in each experiment. Then the relative delay or temporal-shift between a pair
of pixels, where one pixel was selected from the left motor area and the other from the
right motor area, was determined separately from the 8 experiments, and the standard
deviation in the 8 temporal shift values was compared to the Monte-Carlo computed
standard deviations at the experimentally estimated SNR. This process was repeated
for several pairs of pixels. The temporal shift between two pixels, instead of the actual
delay of a pixel with respect to the zero-delay reference function, was used to avoid
any systematic errors due to our lack of knowledge about the precise MR slice
acquisition time with respect to the button press time.
6.2 Consistent Result (vs. Simulation)
The results of the experimental study, where the standard deviation from 8
experiments was compared to the Monte-Carlo simulation, are presented in Table 6.1
for eight representative pixel pairs. The experimental standard deviation σ
exp
arises
from the difference between two time courses whose SNR values are labeled SNR
1
and SNR
2
. To relate σ
exp
to the simulated value σ
sim
, the standard deviations σ
1
and σ
2
corresponding to SNR
1
and SNR
2
were estimated from Fig. 5.4 and added in
quadrature to obtain σ
sim
. Fig. 6.2 shows the block diagram of how to obtain σ
exp
44
from the experimental data and σ
sim
from the estimated SNR through the simulation
look up table (Fig. 5.4).
Dynamic fMRI detection 61 reference functions
Computing voxel pair delay
Select a pair of voxels
(one in ROI-1 and
another in ROI-2 )
σ(exp)
Estimating SNR in
time-course of 2 voxels
σ(sim)
LUT: SNR Æ σ
Dynamic fMRI detection 61 reference functions
Computing voxel pair delay
Select a pair of voxels
(one in ROI-1 and
another in ROI-2 )
σ(exp)
Estimating SNR in
time-course of 2 voxels
σ(sim)
LUT: SNR Æ σ
Figure 6.2. The block diagram of the approach to validate the temporal resolution.
45
The experimental standard deviations are in good agreement with Monte-Carlo
simulations. The assumption that events occur exactly midway between two reference
functions could increase the simulated standard deviation by as much as 50 2 , or
about 70 ms and it is seen that most of the simulated values are higher than the
experimental values consistent with this explanation.
6.3 The Case Studies of dfMRI
To demonstrate that the BOLD signal could be tracked on a pixel-to-pixel
basis, data from three human volunteers were processed by a series of 41 reference
functions per slice that were spaced by Δt of +200 ms and spanned a range of +4 s
with respect to the zero-delay reference function. The longer range was used to enable
tracking of pixels in the supplementary, primary and sensory areas, which may be
Table 6.1.
Comparison of experimentally determined standard deviation σ(exp) in ms,
to that obtained from the Monte-Carlo simulation σ(sim). The SNR1 and
SNR2 denote signal to noise ratios for the two pixels in each pair [54].
46
several seconds apart, and the 200 ms interval step was selected to reduce the
computation time. Shorter Δt values could be used in the future. Also, instead of
assigning a pixel to only that reference function which produced the highest
correlation coefficient, a pixel was assigned to all reference functions whose
correlation coefficient values with respect to the pixel's time course were within 1% of
the maximum value (Eq. 4.4). This 1% tolerance is consistent with the 200 ms
spacing because the correlation coefficient is reduced by approximately 1% between
two reference functions when one is displaced 200 ms with respect to the other on this
protocol. In addition, the tolerance allows a pixel to persist through different delays or
disappear and reappear through the activation time sequence, both of which may occur
in practice. Pixels whose correlation coefficient values exceeded a threshold of 0.3
(corresponding to p < 0.00001) were displayed in image frames as a function of Δt.
6.4 Movie Frames: The dfMRI Results
Examples of pixel tracking for two different subjects, where subject 1 moved
his right index finger and subject 2 moved her left index finger are presented in Figs.
4.1 and 6.3 respectively. Though this study was conducted with Δt = 200 ms
reference functions, selected images corresponding to reference functions or frames
delayed by Δt = 400 ms are shown for lack of space. The spatiotemporal relationships
among different regions of the brain are visible in these images. These relationships
are best depicted by a movie where one can actually visualize the path of the BOLD
signal as it propagates through the vasculature.
47
The results from all three subjects show that various SMAs are activated first,
followed by the contralateral and ipsilateral PMAs, and then the SAs. The start of
activation at certain locations (referred to as activation nodes) followed by propagation
of the activation through regions which geometrically and anatomically conform to
veins can be seen from individual frames. The frames also reveal that activation
within large veins is delayed by several seconds with respect to the activation nodes,
which is consistent with previous studies [34, 52].
Figure 6.3. Selected frames 400 ms apart showing the spatiotemporal relationships among pixels
activated within the SMA, PMA and SA during the finger movement experiment for subject2. The
subject moved her left finger [54].
48
6.5 Analyzing the dfMRI Results
Using the metric ‘A’ to quantify activation, where A = average z-score times
the number of pixels, the onset time and the peaking time of all activated clusters
(where a cluster was defined by a minimum of 4 contiguous pixels) were determined
to quantify their temporal relationships. An example of this analysis for subject # 2 is
presented in Fig. 6.4.
The activated SMA clusters lying within both hemispheres of the brain are
shown in the image displayed at the top in Fig. 6.4(a). This image is derived from a
composite of all frames. The various SMA regions have been labeled R1-R4 for the
right hemisphere and L1-L2 for the left hemisphere. The plots of the metric ‘A’
within each cluster as a function of latency are shown below the image. These plots
group the latencies of pixels within each cluster. Latencies corresponding to the
cluster onset time (which is presumed to correspond to activation within the
microvasculature) or the peaking time (which should correspond to activation within a
relatively large vein) can be obtained from these plots. Similar clusters and plots of
‘A’ in the PMA and SA are shown in Figs. 6.4(b) and (c) respectively. Since the
microvasculature is presumed to underlie regions of neuronal activity, the delays
between cluster onset times should reflect the temporal relationship between regional
neuronal activity.
As an example of determining temporal shifts between two activated regions,
the delay between the first activated cluster in the SMA and the first activated cluster
in the PMA was computed. The latency of the onset time of the first SMA cluster L1
49
in Fig. 6.4(a) was 2.9 s and the latency of the first PMA cluster R4 in Fig. 6.4(b) was
3.5 s, suggesting a temporal delay of 600 ms between these regions. Similar plots for
the other two subjects indicated delays of 1100 ms and 1000 ms respectively. These
delays are consistent with previous fMRI and evoked potential studies [57].
In addition to the first-activated SMA and PMA pixels, several other clusters
within the SMA, PMA and SA were activated in fMRI as seen in Figs. 6.2-6.3. This
work provides a tool to measure the temporal delays among regions based on the
assumption that the timing of neuronal activity is correlated to the onset timing of
pixels within a cluster.
Figure 6.4(a). Clusters in the SMA and plots of corresponding activation metric 'A' defined in
the text [54].
50
Figure 6.4(b). Clusters in the PMA and plots of corresponding activation metric ‘A’ defined in
the text [54].
Figure 6.4(c). Clusters in the SA and plots of corresponding activation metric 'A' defined in the
text [54].
51
6.6 Directly Estimating Temporal Resolution from 2-Event Studies
In previous sections, the fMRI temporal resolution is estimated from a
simulation and a multi-session motor study. It is possible to estimate the resolution
directly from an fMRI study with 2 types of event with varying timing as presented
here. The experiment was conducted on 2 subjects.
6.6.1 Visio-motor Tasks: The Method & Design
Under a 1.5T GE MRI scanner, two volunteers (5 experiments on subject # 1,
and 4 experiments on subject # 2) were scanned (TR = 800 ms, 512 images per
experiment) while performing a motor task of pressing a fiber-optically coupled button
as soon as they heard audio cues at 4000 + DT ms (where DT = +400, -300, -100,
+200, +0 ms) after a 1-second-interval visual stimulus (7.5 Hz pattern reversal) at
every 20.1 seconds. (By the design, events are staggered with 100 ms interval.) The
actual motor response differential timing (dt) was retrieved from the button-press
records. Approximately, DT = dt, but the accurate value, dt, was used in this analysis.
In using dt, the inter- and intra-subject variability of motor responses could be
avoided. The detected timing of the visual response (t
v
) of a voxel in the visual cortex
was considered as the reference timing (t = 0) in order to measure the detected timing
difference (Δt) of 2 motor response timings (t
m1
and t
m2
) of a voxel in the motor
cortex. Fig. 6.5 shows the design diagram of the experiment. There should be no
difference between Δt and dt., i.e., (Δt - dt) = 0. The standard deviation of (Δt - dt)
52
t
v
t
m2
dt
Visual
Stimulus
Motor Responses
Detected Timing
= Δt
t
m2
-t
m1
t
m1
Designed Timing
0
t
v
t
m2
dt
Visual
Stimulus
Motor Responses
Detected Timing
= Δt
t
m2
-t
m1 = Δt
t
m2
-t
m1
t
m1
Designed Timing
0
Figure 6.5. The design diagram of the experiment with visual and motor tasks to estimate fMRI
temporal resolution by using dfMRI.
from the detected activated voxels can represent the estimated fMRI temporal
resolution from these experiments.
6.6.2 Data Processing: dfMRI Processing
Separately, the visual and motor fMRI data were processed by the dfMRI
methods to obtain the timing of functional response for each voxel. In using the
detected visual response timing (t
v
) as the reference timing (t = 0) to compute
detected motor timing difference (Δt), inter-region hemodynamic response variability
can be ignored. The designed motor response differential timing (dt) and measured Δt
of pairs of early activated visual and motor voxels were analyzed.
53
Figure 6.6. Selected frames 400 ms apart (subject # 2) showing the spatiotemporal relationships
among voxels in the visio-motor experiment.
Fig. 6.6 shows the examples of dfMRI results of an experiment on subject # 2.
The differences of brain activations can be easily observed during a 400-ms time
interval. Only first 5% early activated voxels were used to calculate Δt’s because they
are spatially and temporally closest to the active neurons. As previously mentioned in
Chapter 3, the main and late activated voxels (as detected by a standard fMRI method)
are mainly associated with the big draining veins that are not only spatially and
temporally far away from the active neurons but also subject to anatomical variability.
6.6.3 Results and Discussion
In Table 6.2, we can observe that the means of (Δt - dt) are very close to zeros
(1.35 and -8.79 ms) in subject # 1 and subject # 2 as expected. The standard
deviations of (Δt - dt) = 126 and 93 ms, respectively. Therefore, the shift in timing
between 2 stimuli can be determined to the accuracy of approximate 100-150 ms.
54
Table 6.2.
The difference of designed differential timing and detected differential timing in the visio-motor
event-related experiment processed by the dfMRI methods.
55
Chapter 7 Applications of dfMRI
Because dfMRI can extend standard fMRI to one more dimension
(latency/time), it is useful in analyzing, understanding, and visualizing fMRI
experiments. Most importantly, spatiotemporal dynamics or networks of brain
activations according to a task can be easily conceptualized. A classic event-related
motor study was presented in Ch. 6. More of dfMRI benefits including some for
clinical investigations are presented in this chapter.
7.1 Motor Task in Patients
Even in a standard (box-car design) fMRI study, dfMRI can also improve the
detection of BOLD responses by not only allowing an investigator to observe early
and late BOLD responses but also accommodating the possible asynchrony between
the experimental design and the patient responses. Most importantly, motion artifacts
could be observed earlier than the BOLD responses.
Figure 7.1 is the resultant fMRI after performing the dfMRI detection, and then
avoiding motion artifact timings. A brain tumor patient moved his right hand for 40
seconds and rest for 40 seconds, and then repeated this sequence twice. For the left
hand movement study, the same paradigm was used. The scanning parameters were:
R
T =5.02s (effective), TE=45ms, FA=90°, # of slices=5, # of images per slice=100,
slice thickness=8mm, FOV=40x20mm
2
, acquisition matrix=128x64.
56
7.2 Picture Naming Task to Identify Language Areas
An fMRI picture-naming task can be utilized to identify language areas for
neuro-surgical planning. Because the brain dynamically uses many areas to perform
the task, dfMRI is a perfect tool to study brain responses. Typically, there are
activations in the visual cortex first, then in the Broca and the motor areas involving
speaking. There may be activations in the auditory cortex because the subject hears
own spoken words. During the experiment, every 18 seconds, the subject was shown
a picture to be named immediately.
In an example result (Fig 7.2), there were activations of speech areas, which
have been displaced by the tumor. Also activations in the visual cortex and the motor
cortex (due to speaking) could be observed.
Figure 7.1. A motor study in a brain tumor patient. Top row - activations associated with right
hand movement. Bottom row - activations associated with left hand movement.
L
57
7.3 Artifact Identification and Reduction
Even though a motion correction algorithm (e.g., the realignment process in
SPM) can reduce activation artifacts, it can hardly remove artifacts due to movements
associated with tasks, stimuli, or responses. Because dfMRI can locate when the
movements occur precisely, the time courses of the artificial voxels can be confidently
selected to fabricate a motion artifact model or a covariate to be used in a statistical
detection method. In addition, the artifact timing might help an investigator to
confidently identify the real fMRI activations or dynamics.
7.4 Working Memory
When an fMRI task becomes complicated or lengthy such as that in a working
memory study, a standard fMRI might not be adequate enough to detect or analyze the
brain activation dynamics. Sequentially, sub-events can be detected and displayed by
Figure 7.2. A study to identify speech areas by naming pictures.
R
L
58
dfMRI. The main effect of interest may be a differential comparison (> or <), while
another event may be stimuli or responses such as hearing or speaking.
Fig. 7.3 demonstrates the usefulness of dfMRI in detecting artifacts due to
speaking in a working memory study (sorting vs. repeating numbers). The artifacts
occurred 4 seconds before the sluggish BOLD response for the working memory
because the artifacts due to the verbal responses are prompt but the BOLD
hemodynamic response peaks approximately 4-5 seconds later.
Figure 7.3. In a working memory study, artifacts due to subject’s verbal responses were
detected at 7 seconds by a reference function designed to detect prompt effects, while the
working memory activations were detected 4 seconds later.
59
7.5 The Software at Work: Step by Step
Currently available software packages for processing fMRI data are relatively
slow and inadequate to address timing uncertainties in the stimulus on-set timing,
times at which specific task are performed, variability of individual hemodynamic
responses, and motion artifacts related to the task. As presented earlier, dfMRI can
alleviate the influences of time uncertainty and other factors, and then discover
normally undetected responses by adding a temporal dimension to the detection
process.
Thanks to the compactness of the formula (Eq. 4.1), the processing time of
dfMRI is not a problem even with the added dimension of time if the algebraic engines
of the software are optimized for vector-matrix operations. We have designed
software incorporating vector-matrix operations. An example of step by step
processing applied to a motor study of a tumor patient for neurosurgical planning is
summarized as below:
1. Generate multiple (in a range from negative to positive delays) reference
functions for multiple slices from the experiment design (Fig. 7.4(a).) while
waiting for image data being transferred from the MRI scanner
2. Import EPI data time-series into the program as a matrix (Fig. 7.4(b).)
3. Import T1 anatomical images on which the resultant activations will be
superimposed (Fig. 7.4(c).)
4. Register the T1 images to the EPI data time-series (Fig. 7.4(d).)
5. Make a mask in order to mask out regions outside the brain (Fig. 7.4(e).)
60
6. Start the dfMRI processing engine (Fig. 7.4(f).)
7. Visualize the timeframe results (Fig. 7.4(g). - Fig. 7.4(k).); normally, motion
artifacts appear before BOLD responses; spatiotemporal dynamics of the brain
responses can be observed; optimal results can be chosen (Fig. 7.4(l).).
61
Figure 7.4(a). Running the software: making multiple reference functions.
Figure 7.4(b). Running the software: importing EPI data time-series.
62
Figure 7.4(c). Running the software: importing T1 anatomical images.
Figure 7.4(d). Running the software: registering T1 images to EPI time-series.
63
Figure 7.4(e). Running the software: making a brain mask.
Figure 7.4(f). Running the software: running dfMRI engine.
64
Figure 7.4(g). Running the software: resultant timeframes 1-4.
Figure 7.4(h). Running the software: resultant timeframes 5-8.
65
Figure 7.4(i). Running the software: resultant timeframes 9-12.
Figure 7.4(j). Running the software: resultant timeframes 13-16.
66
Figure 7.4(k). Running the software: resultant timeframes 17-20.
Figure 7.4(l). Running the software: optimal views of sensorimotor activation for
(top 2 rows) right hand movement and (bottom 2 rows) left hand movement on 10 slices.
67
Chapter 8 Neuromagnetic Signal Detection with MRI
As mentioned earlier (Ch 3.6), a practical limitation of positive BOLD fMRI is
that the signal is not confined to the microvasculature but spreads into draining veins,
which may be several millimeters or up to few centimeters away from the actual site
of neuron activities. Thus, the identification of the actual site of activation and its
temporal behavior, i.e., the spatiotemporal resolution of BOLD fMRI is limited not
only by the hemodynamic response characteristics but also the vascular architecture.
The temporal response of the BOLD signal is modeled by a so-called hemodynamic
response function or HRF (Ch. 2.3.1). Though the HRF exhibits a relatively slow
peaking time of about 3 to 5 seconds, temporal shifts on the order of 200-300 ms have
been extracted from activated regions using the dfMRI, a specialized methodology
(Ch. 4, 5, and 6).
Partly based on the dfMRI, a new method to detect neuron activities with MRI
is presented in this chapter. The technique does not depend on the BOLD property as
fMRI and dfMRI, instead, it needs to avoid or to filter the BOLD effects out. Most
importantly, not only the direct detection of neuron activities is possible with MRI but
also the spatiotemporal dynamics with a very fine temporal resolution of the brain
activations according to a task can be studied and visualized. The crucial concept of
the technique is to control as well as to exploit the aliasing property of a signal.
68
8.1 Not a Blood Related Functional Study
In contrast to a blood-based detection method such as fMRI and PET,
Electroencephalography (EEG) and Megnetoencephalography (MEG) can also be
utilized to investigate brain functions non-invasively. However, both require solving
inverse problems to estimate the spatial information of the neuron activities. Accurate
spatial information can be obtained by MRI. Hence, if MRI could be utilized to detect
neuron activities (neuromagnetic signals) directly instead of the BOLD effect, the
accurate spatial information of neuron activities could be achieved. However, the
magnetic signals caused by electric signals in working neurons are very small for
signal detection with MRI [48, 49].
8.2 B
0
and Working Neurons
Mentioned in Ch. 2.1, the main magnetic field B
0
is the foundation of MRI. If
there is an interference with B
0
, there should be a change in the MR images. Working
neurons produce electrical currents generating a magnetic field change, interfering
with the B
0
.
An approach to image the electromagnetic activity of neurons directly by MRI
by detecting accumulated phase shifts in the B
0
was first reported by Singh, M. in
1991 [48] and a subsequent study published in 1994 [49] suggested that the expected
small signals were probably overwhelmed by physiological noise during visual and
auditory stimulation at 1.5T in humans. The basic idea draws on the well-known
relationship between the phase angle and the z-component of the external magnetic
69
field (which includes B
0
and the z-components from any other magnetic field) at a
given location [25, 29, 45]. The phase-shift at any location, i.e., the shift in phase after
removing the effect of B
0
, is directly proportional to the z-component of the evoked
neuromagnetic field, where the neuromagnetic field refers to the magnetic field
produced by neurons as mentioned in MEG [48, 49, 50]. If ‘ N
zi
’ is the z-component
of the neuromagnetic field at location ‘i’, then the phase-shift ‘φ
i
’ at location ‘i’ is
given by φ
i
= γN
zi
t where ‘γ’ is the gyromagnetic ratio and ‘t’ is the duration of the
neuromagnetic field. The phase shift may be positive or negative depending on the
orientation of the z-component of the neuromagnetic field, and if the neuromagnetic
fields within a voxel change with time, the phase shift represents an integration over
time ‘t’. Thus, the electrical activity of neurons could be imaged directly by mapping
phase-shifts with an appropriate MRI pulse sequence.
Test object and human studies reported by Singh [48, 49] using stimulus
locked spin-echo and gradient echo sequences suggested that successful neural activity
imaging would require the capability to detect a phase shift of 0.35 degree above the
noise. However, noise from cardiac, respiratory and movement-related effects
apparently overwhelmed the neural signals in the measurements made at that time, and
results of human studies including visual and auditory stimulation were equivocal.
The approach to image working neurons directly has been rediscovered
recently [6, 7, 30] but whether there would be adequate signal from phase-shifts in
human imaging to overcome physiological noise remains uncertain. An alternative
approach based on detecting phase dispersion or dephasing and the ensuing T
2
*
70
reduction by neuronal currents, rather than phase-shift, has shown promising though
unconfirmed results [60]. Dephasing within a voxel is likely due to different
orientations of neurons and variations in neural current among these neurons, leading
to intra-voxel fluctuations in the z-component of the neuromagnetic field. The
approach [60] relies on comparing images acquired in the first 1-second after stimulus
presentation to those acquired during the next 1-second. MRI acquisitions were
synchronized to a visio-motor stimulus so that these differences could be detected in
100 ms frames. Contamination from the BOLD signal, which would limit this
approach, was reduced by omitting the first 20-seconds of data under the assumption
that the BOLD signal would thereafter be in a steady-state. However, it is well known
that even in steady-state, there are significant fluctuations in the BOLD signal from
physiological effects such as cardiac pulsations and respiration, which could easily
mask the approximately 1% signal expected from neural currents. At
R
T = 1 second
used in their work, it is difficult to separate the cardiac and respiration related BOLD
signals from neural currents and thus the reported results [60] remain equivocal.
Another recent study where a fast stimulus presentation sequence was embedded
within a slow stimulus presentation to acquire both the BOLD signal and possible
neural current signals failed to find a significant neural current signal at 3T using a 16-
channel coil [11]. Correlation with the slow presentation was used to fit the relatively
slow BOLD signal, which was then separated from the fast signal by a second
correlation designed to fit the fast neural currents. A relatively short
R
T of 100ms was
used, which would allow separation of signals from cardiac and respiratory
71
fluctuations in the data, but at the expense of signal to noise ratio (SNR). Though the
reduction in SNR was partially compensated by the increased sensitivity of the 16-
channel coil, the SNR reduction may be a possible reason for not being able to detect
neural currents. Also, the exact shape of the BOLD HRF necessary to fit and remove
the slow BOLD signal in this approach is not known and any discrepancy between the
modeled and the actual BOLD signal could influence the detection of the neural
current signal.
8.3 Fine Resolution Modeling and Aliasing
Similar to BOLD fMRI in detection processes, however, a model for detecting
neuronal activity was based on an average event-related potential (ERP) of EEG data
of flash stimuli at 4 Hz (Fig. 8.1) resulting in an approximately inverted skewed
narrow Gaussian curve with a peak at around 150 ms as the model to be detected.
Figure 8.1. Representative EEG using the same 4 Hz visual stimulus. Two main peaks can be
observed, a sharp peak at 50 ms and a broad peak at 125-175 ms [53].
72
Because MRI sampling rate cannot be as high as that of EEG or MEG (ms level), it
seems that it is nearly impossible to detect any fast signal as those in EEG or MEG
studies. However, if the fast signal to be detected is periodic, with an appropriate
sampling rate (even lower than the Nyquist frequency), a detected signal can be the
same shape with the fast signal but at different frequency, i.e., an aliasing signal. The
two crucial points of the experimental design are: (1) the fast functional task or
stimulus must be periodic and (2) as described later, the sampling rate must be suitable
to exploit and to control aliasing signals.
To illustrate the design concept graphically, we assume that there is only one
inverted peak curve signal per period of the 4 Hz flash stimuli (250 ms), and we set an
MRI scanner to scan at TR = 275 ms (i.e., the sampling interval at each anatomical
slice). Then, in every sampled data point interval (275 ms), the fast signal (250 ms
period) has traveled for 275-250 = 25 ms plus 1 period (250 ms). Thus, it will take
250/(275-250) = 10 sampling data points (10x275 = 2750 ms) to record the shape of
the fast signal for 1 period. It is as if the fast signal is expanded from 250 ms to 2750
ms (i.e., a factor of 2750/250 = 11). Fig. 8.2 illustrates the design concept with 2
anatomical sampling slices. The second slice sampling time is at 125 ms after the first
slice. Therefore, the sampled signal of the second slice is delayed by 11x125 = 1375
ms comparing to the first slice sampled signal.
73
In signal processing, because the sampling rate (1/.275 = 3.636 Hz) is less than
the twice of signal frequency (2x4 = 8 Hz), an aliasing effect is produced. Due to the
aliasing property, after sampling, the signal is periodic with the frequency of 4-1/.275
= .3636 Hz (i.e., the period = 2.75 s, or 10 sampled points).
Figure 8.2. (a). A conceptual depiction of how a gradual shift is obtained between the stimulus
presentation and slice acquisition (or slice sampling) time; (b). Formation of reference
functions by a model of the expected neural current or neuromagnetic signal (black line) and its
sampling by the two slices (o and *) [53].
74
8.4 Experiment Design and Formulae: Exploiting Aliasing
In the previous section, a design example to detect fast neuromagnetic signals
by MRI is presented. Ironically, an alias signal is needed instead of being avoided.
However, the alias signal in the sampled data must be controlled by the design to have
the same shape as the original signal. Thus, the number of sampled points to contain
the shape of the fast signal for one period,
T
s
T
n
TT
=
−
(8.1)
must be an integer number. T = the fast signal period (= ISI),
s
T = sampling interval
(= effective
R
T in MRI). If
T
n is positive, the sampled signal has the same shape as
the original fast signal but is expanded by a factor,
Ts
a
nT
s
T
= (8.2)
If
T
n is negative, Eqs. 8.1 and 8.2 are still valid. However, the sampled signal is the
expanded mirror image of the original fast signal.
To recover the detail of the signal shape,
T
n > 5. This is equivalent (in signal
processing) to the requirement that to recover the detail of the signal shape, the
sampling frequency,
s
f , should be much higher than the Nyquist frequency,
N
f ,
where
N
f = twice the signal frequency f .
If
T
n is not an integer number, the sampled signal will be complicated (e.g.,
there are low frequency signals modulating a high frequency signal). It is crucial to
choose a
s
T for a T to get an adequate
T
n . From Eq. 8.1, to get a high
T
n in order to
75
reconstruct the original signal shape,
s
T must be very close to T . However, T is
constrained by a functional task or stimulus in the experiment. Only a function
experiment with a periodic design can utilize this methodology. Moreover, the
task/stimulus interval, T , should be short enough so that the sampling period,
s
T or
R
T , can be short, then a band-pass filter can be applied to the sampled MRI data to
remove some periodic physiological noises, e.g., breathing, pulses from heart beats.
The example in Ch. 8.3 was one of our early experimental designs. In that
experiment,
s
T = 275 ms, T = 250 ms,
T
n = 10 and the scaling factor,
a
s = 11.
Because the sampled signal frequency = 1/(11x.250) = .36 Hz, while, approximately,
the heart beat frequency = 1 Hz and the breathing frequency = .3 Hz, a band-pass filter
can be designed for the sampled MRI time-series data.
8.5 Slice Timing Correction: A Special Case
In general, an fMRI experiment is conducted with a sampling frequency,
s
f or
1/
s
T , which is much higher than a task or stimulus frequency, f or 1/T . Even in an
event-related experiment or a jittering design, the
s
T or
R
T is still relatively short.
Thus, in a slice timing correction (Ch. 2.3.6), one can simply shift the time-series data
from a different slice or a model by a fraction of
R
T . However, in the new design to
detect a fast neuromagnetic signal, the sampled signal is as if expanded with a factor
of
a
s . Thus, a slice timing correction must include the scaling factor, not just a
fraction of
R
T anymore. The effect can be clearly observed in Fig. 8.2.
76
8.6 The Experiment
Human data were acquired on a 1.5T GE EXCITE scanner from two 8 mm
thick contiguous slices bracketing the calcarine fissure using 512-image time-series
EPI per slice with
R
T = 275 ms, TE = 45 ms, flip angle = 30 deg. A robust visual
stimulation paradigm where a checkerboard was illuminated by a flash at 4 Hz was
used [51]. A series of reference functions to detect the neural activity as a function of
latency were designed [54] in accordance with the EEG measured latencies from the
same stimulus. A representative EEG recording is presented in Fig. 8.1 showing a
narrow peak with a latency of about 50 ms and a prominent wide peak at a latency of
125-175 ms for this stimulus paradigm. Though peaks may be positive or negative, the
phase dispersion will reduce the T
2
* weighted signal in either situation. Thus the
series of reference function was designed to detect a negative peak as a function of
latency. An example of how a particular reference function within the series was
formed is presented in Fig. 8.2(b). The black line represents the effect of neural
currents on the T
2
* weighted signal, showing a dip corresponding to latency. Here we
have selected the latency of the prominent EEG peak (125-175 ms) to create the dip
(inverted skewed Gaussian curve). The blue open circles represent the time at which
slice 1 is acquired and the red stars the time at which slice 2 is acquired. Thus the blue
and red lines trace the modeled MRI signal time-course for slices 1 and 2 respectively,
and become representative reference functions for these slices. Since the latencies of
the peaks are in general unknown, these representative time courses were shifted in
time to form a series of reference functions. Consistent with the choice of
R
T and TI, a
77
series of 10 reference functions, where each successive reference function was shifted
by 25 ms with respect to the previous function, was used in this study. The resultant
shifting will be expanded as mentioned in Ch. 8.5.
SPM (http://www.fil.ion.ucl.ac.uk/spm) was used to detect the activity
correlated to each reference function at a p < .05. Thus each detected pixel also had a
temporal marker to identify its latency depending on the particular reference function
that detected it [51]. As indicated earlier, any temporal sampling interval for the
latency could be achieved in principle by adjusting the difference between TR and TI.
However, shorter the temporal sampling, more images would be required to complete
one cycle or epoch and thus fewer cycles would be contained within a given
acquisition time. Thus the latency sampling is limited by the statistical power to detect
signals in this approach.
The stimulus was turned on for 30 s before the MRI measurements began, and
left on during the entire time-series acquisition. This was to ensure that the
hemodynamic response (i.e., the BOLD signal) was driven into a steady-state and
could be considered as a baseline. Likewise, the initial "negative dip" and any flow-
related effects would also be driven into the baseline due to the steady state 4Hz
stimulation. However, fluctuations in the baseline, including those in the resting
BOLD signal, and physiological effects arising from cardiac and respiratory pulsations
represent potential sources of error to be avoided.
78
8.6.1 Results
Three subjects were scanned. The first and second subjects (subject 1 and 2)
were each scanned four times to obtain four independent measurements using the 4 Hz
flashing checkerboard stimulus. The third subject was scanned three times. An
example of images showing pixels detected at p < .05 (using SPM) by the 10
successively shifted reference functions is presented in Fig. 8.2. These images
represent activity detected as a function of latency in steps of 25 ms. The top two
rows show images obtained for the two slices for subject 1, study 1, and the next two
rows show the results for the same subject for a repeat scan. Results for similar two
studies for subject 2 are shown below subject 1 in the same order. Activity is clearly
seen in the visual cortex at different latencies. For comparison, an example of
conventional event-related BOLD imaging for the same stimulus but TI = 20 s. at p <
.01 is shown in Fig. 8.3. Ample activity is produced in the visual cortex by this
stimulation, but as expected, most of the activated pixels lie in the macrovasculaure.
The images in Fig. 8.4 are confined to gray matter and likely correspond to neural
activity co-localized with the microvasculature.
The t-scores were computed and summed within the visual cortex at p<0.05
from images shown in Fig. 8.4 for subjects 1 and 2 including the additional two
studies per subject not shown here. The t-scores were similarly summed for the three
studies for subject 3. All t-scores are listed in Table 8.1. A plot of the t-score for the 3
subjects, averaged over all scans per subject, as a function of latency is presented in
79
Figure 8.3. An example of conventional event-related BOLD imaging for the flash checker
board stimulus with TI = 20 s. at p < .001 for the same two subject shown in Fig. 8.4. Strong
activation is in the macrovasculature as expected [53].
Fig. 8.5. Signals were typically ~ 1%. These plots are consistent with the EEG where
peaks are seen at latencies of approx 50 ms with a prominent peak at 125-175 ms.
An example of the time-course from one prominently activated pixel is
presented in Fig. 8.6. The selected pixel from subject 2 at a latency of 125 ms is
indicated in Fig. 8.6(a), the average time-course with averaging computed over all
epochs (similar to what is usually performed in EEG) is shown in Fig. 8.6(b), the full
time-course is shown in Fig. 8.6(c) and the power spectrum is shown in Fig. 8.6(d).
80
Figure 8.4. Images representing activity detected as a function of latency in equivalent steps of
25 ms (mapped to 275 ms) [53].
81
Table 8.1.
Sum of t-score for activated voxels lying within the visual cortex at p < .05 [53].
Figure 8.5. Plot of t-score averaged over all scans. They are consistent with EEG latencies [53].
82
8.6.2 Discussion
A novel strategy to present stimuli and acquire MRI data has been presented
which allows us to obtain images as a function of latency following the stimulus on-
set. By adjusting the stimulus presentation interval TI and the MRI repetition rate
R
T ,
it becomes possible to acquire images as a function of latency in the milliseconds
range. The sampling time, i.e., the latency resolution, is limited by the number of
repeated cycles or epochs that would be contained within a given acquisition time.
Thus the statistical power to detect the signal would limit the latency resolution. In this
work we selected TI = 250 ms (4Hz) and
R
T = 275 ms as representative values to
evaluate the feasibility of this approach.
The value of
R
T and the number of images contained within a cycle or epoch
also determines the frequency of the signal of interest. Since cardiac, respiration and
very slow (<0.1 Hz) fluctuations are potential confounding sources, the frequency of
interest must be selected to allow differentiation of these confounds from the effect of
interest. Cardiac fluctuations range from about 0.8 to 1.3 Hz and respiration
fluctuations are usually centered around 0.3 Hz. The
R
T , TI combination selected in
this work yield a signal of interest at 0.36 Hz, which is separable from these
confounds. Filtering techniques specifically tailored to extract the frequency of interest
could be used in the future to avoid these confounds.
The T
2
* signal from neural currents is expected to be ~ 1%, which is very
close to baseline fluctuations attributed to MR instrumentation noise and noise due to
the fluctuations mentioned above. In this work we used the latency of the detected
83
signals and its consistency with measured EEG signals from the same stimulus to
suggest the feasibility of this approach. There remains a finite chance that despite our
methodology to avoid confounds, the detected signals are not due to neural currents
but occur through some combination of cardiac and respiratory fluctuations. Future
studies using different combinations of
R
T and TI and filtering techniques to reduce
the contributions from other frequencies are needed to convincingly demonstrate the
potential of this approach.
8.7 Spatiotemporal Dynamics Once Again
Similar to the spatiotemporal dynamic fMRI in previous chapters, the proposed
methodology can provide the spatiotemporal dynamics of brain activities according to
a task or stimulus. However, it does not require the indirect effects in the blood
supply systems in the brain (BOLD property) as fMRI or dfMRI, in contrast, it detects
magnetic field interferences due to electrical currents in working neurons directly.
Moreover, a fine temporal resolution in fMRI (tens ms level) can be achieved with a
scaling factor (1/
a
s ).
84
Chapter 9 Neuromagnetic Signal Detection Difficulties
Direct detection of neuromagnetic signals with MRI is difficult because the
signals are very small and fast, relative to the capacities of MRI scanners. However,
the proposed method in the previous chapter shows that it is possible by a special
design to exploit and to control the aliasing signals. Much more than those in a
standard fMRI, this method poses several problems and difficulties, addressed below.
9.1 Noise I: Low Frequency
Because the fast signal to be detected is an aliasing signal, the sampled signal
frequency becomes much lower (Eqs. 8.1 and 8.2). As in the example experiment
(4Hz flash stimulus, effective
R
T = 275 ms) in the previous chapter, the sampled signal
frequency = .36 Hz, which is very close to the breathing frequency. Therefore, the
breathing can easily cause fault positives in the detection. Moreover, a noise such as a
subject movement can cause artifacts in the detection of low frequency signal. Thus it
is very important to choose the signal interval, T , and the sampling interval,
s
T , to
achieve an optimal aliasing signal frequency. Then a filter can be designed to remove
the low frequency noises.
85
9.2 Noise II: High Frequency
A high frequency noise such as that from heart beats (approximately, 1+ Hz)
can reduce the detectivity of the small signal of interest. However, it can be filtered if
the
s
T or
R
T is short enough (e.g., 200-400 ms) with a low-pass filter.
9.3 Sampling and Task/Stimulus Rates
The sampling rate (
s
f , 1/
s
T , or 1/
R
T ) must be very close to the task/stimulus
frequency, 1/T , to produce the desired number of sampling points per signal period
T
n , and the scaling factor
a
s . A jittering or random task/stimulus design cannot be
applied in this method, because a periodic signal is needed to produce a controlled
aliasing signal.
9.4 The Possible Solutions
As mentioned, it is crucial to choose the stimulus interval T, and the sampling
interval
s
T or
R
T . Moreover, a suitable filter should be applied in preprocessing the
sampled time-series image data.
9.4.1 Design
To achieve a wider gap between the breathing frequency and the aliasing
frequency, an experimental design with T = 200 ms and
s
T = 220 ms would be more
suitable as the aliasing frequency = .45 Hz, comparing to .36 Hz in the previous
experiment.
86
9.4.2 Filters
A good band-pass digital filter can be designed and applied to the data. The
wider gap between the aliasing frequency and the breathing should lead to better
results. The lower cutoff of the band-pass filter can be set at the mid point of the
breathing frequency and the aliasing frequency, while the higher cutoff would be just
below 1 Hz.
87
Chapter 10 Neuromagnetic fMRI: An Optimal Design
Because the presented neuromagnetic functional MRI (nfMRI) is based on
directly detecting the aliasing signals of working neurons by design, the signals to be
detected in image time-series are in fairly low frequency band. The signal in the pilot
experiment (Ch. 8) was at 0.36 Hz., which was close to the breathing rate
(approximately 0.3 Hz.). The sampling rate was at 4 Hz., which was good enough in
order not to have aliasing signals from the heart rate. However, if the subject’s heart
rate is faster than 109 per minute, there will be aliasing signals. If a subject is a little
excited during the experiment in an MRI scanner, her or his heart rate can possibly
reach 109 per minute.
Fig. 10.1 shows the power spectrum of a voxel time-course from a subject in
Ch. 8. The aliasing signal of interest peak of 0.36 Hz. can be observed, so can the
breathing rate (0.3 Hz.) peak and its harmonic peak (0.6 Hz). The heart rate frequency
peak is distinct at 0.9 Hz.
Figure 10.1. Power spectrum of a voxel time-course showing the aliasing signal of interest peak
and physiological signal peaks.
88
10.1 Choosing Aliasing Frequency of the Signal: The Design
A desired aliasing frequency should be far from breathing frequency and its
harmonic frequency. The sampling rate, 1/
R
T , should be high so that there will be no
aliasing signals from the heart rate. However, if
R
T is too short such as less than 200
ms, the SNR in the data image time-series will be quite poor.
Table 10.1 shows examples of design parameters where the stimulus
frequency, 1/T and sampling interval,
R
T are in the first two columns. The third
column is the aliasing frequency of the signal to be detected due to the sampling
rate,1/
R
T . The fourth column is the number of samples per one period of the aliasing
signal, n
T
, computed by Eq. 8.1. The last column, based on Nyquist’s theory, shows
the maximum frequency of a signal that will not cause aliasing frequency. The
parameters of the experiment in Ch. 8 are in the 5
th
row.
Table 10.1.
Examples of design parameters for neuromagnetic signal detection with MRI.
Stimulus TR (eff) Aliasing signal # of samples Max. unaliased
(Hz) (ms) frequency (Hz) per signal period frequency (per min)
1 1200 0.17 5 25
2 600 0.33 5 50
3 400 0.50 5 75
4 300 0.67 5 100
4 275 0.36 10 109
5 240 0.83 5 125
5 220 0.45 10 136
6 200 1.00 5 150
8 150 1.33 5 200
10 120 1.67 5 250
89
The higher the integer number of samples per signal period, n
T
, the better
shape of the signal can be reconstructed. From the table, the 7
th
row parameters seem
to give an optimal design. The aliasing signal to be detected is at 0.45 Hz., the mid
point between the breathing frequency and its harmonic frequency. The
R
T is 220 ms
which is not very short so that the image data will not have a poor SNR. Moreover,
the subject’s heart rate can be as high as 136 per minute without having aliasing
signals. Finally, there are 10 samples per aliasing signal period, which should give a
fairly detailed reconstructed signal.
10.2 Designing A Digital Zero-Phase Bandpass Filter
Because the frequency of the signal of interest is known, a digital bandpass
filter can be designed to remove or reduce signals of non-interest such as breathing, its
harmonic and heart beating. The filter must be a zero-phase filter so that it will not
interfere with the detection process and the interpretation of results.
The digital filter was designed by using MATLAB [39]. In filter designing,
the Least Square method for the filter order of 21 were used to get a linear phase filter
with passing band between .405 and .505 Hz. The resultant filter gain at 0.45 Hz is
0.5, thus the final gain for zero-phase filtering is .25 for the aliasing signal of interest.
10.3 The Experiment
The experiments were performed on 5 subjects. The settings and equipments
were very similar to those in Ch.8, but the
R
T = 220 ms, and the flash frequency = 5
Hz (T = 200 ms).
90
10.4 Data Processing Methods
Every voxel data time-series was filtered by the designed linear phase filter,
and then SPM2 was used to detect negative signals (compared to the data points
synchronized with the flash). The first reference function is an impulse train with an
impulse every 10 data points. The impulse was at the time point that the flash light
was synchronized with the scanner scanning a particular slice. The next reference
function is simply the first reference function but the impulses were shifted by one
data point. Because the signals to be detected were periodic with the period of 10 data
points, there are 10 reference functions to be used in data processing. Nine contrasts
were set to detect the negative signal compared to the signal at the synchronized flash
timing. Even though the data point interval = 220 ms, it represents 20 ms interval in
the original fast neuromagnetic signal (s
a
= 11, explained in Ch. 8.4). Finally, the data
need to be processed slice-by-slice because the standard slice timing correction cannot
be performed. If the fast signal is shifted by 100 ms (approximately,
R
T /2), the data
have to be shifted by 100x11 = 1100 ms, which is equivalent to 1100/220 = 5
R
T .
10.5 Results and Discussion
At corrected p < .001, the detected neuromagnetic brain activities (responding
to flash light at 5 Hz stimulus) of 5 subjects are shown in Fig. 10.2. Activities in the
visual cortex can be observed. Most activities are at 80-120 ms after the flash
stimulus. Sum of t-score in the visual cortex confirms the observations (Fig. 10.3).
The preprocessing with the bandpass filter may explain why the detection threshold
91
Figure 10.2. Detected neuromagnetic signals at phases: 20, 40, …, 180 ms on 2 selected subjects
(S1-S2), corrected p < .001. S1 and S2 show activity superimposed on echo planar images.
Figure 10.3. Plot of t-score in the visual cortex ROI of 5 subjects.
92
can be set at a very high value (corrected p < .001) in these experiments compared to
those in Ch.8 (uncorrected, p < .05). Unfortunately, we did not have EEG data
available to us for a 5-Hz stimulus. However, it can be expected that the EEG signal
would be an inverted Gaussian curve similar to that observed in the EEG experiment
with 8-Hz flash stimulus, which would be consistent with Fig. 10.3.
93
Chapter 11 Future Work and Contributions
Even though dfMRI enables us to detect and visualize brain activation in
greater details than the standard fMRI by adding the temporal dimension, still further
research should be conducted on improving spatiotemporal relation in the information-
rich MRI data. High field MRI enables high SNR, therefore the time-course of
activation detected by dfMRI should lead to even better spatiotemporal localization at
higher field strengths. Recently, physiologic brain networks have been uncovered by
fMRI. We hope that the temporal dimension of dfMRI will expand the possibility to
better understand brain spatiotemporal networks.
Using fMRI as a spatial constraint, EEG and MEG inverse solutions have been
improved lately. By their nature, EEG and MEG have enormous temporal
information. Thus, with the spatiotemporal information from dfMRI to further
constrain the solvers, accurate spatiotemporal solutions are expected, especially in
situations where EEG and dfMRI are conducted simultaneously.
The need of a periodic stimulus or task by nfMRI is its limitation to be used in
neuro-scientific and physiologic research. However, its ability to directly detect a
brain signal at a certain voxel may be of significant use in many situations where
periodic stimuli are used. Moreover, its equivalent temporal resolution can be as fine
as 20 ms as presented in the last few chapters.
The nfMRI results should be improved with high field MRI and the ability to
control the MRI scanner to precisely synchronize with the stimulus.
The challenge will never end.
94
11.1 Original Contributions
During the past years toward my Ph.D., I have come with ideas and
implemented them to improve fMRI. In the previous chapters, I have categorized
them into dfMRI (dynamic functional MRI) and nfMRI (neuromagnetic functional
MRI). At the end of my thesis, I would like to itemize them and give brief
explanations.
11.1.1 Accurate Slice Timing Correction
In general, a series of volume images is acquired in an fMRI experiment. In an
image volume, there are many slices that are acquired at different time. Most fMRI
processing programs shift the (noisy) time-series data to a reference slice timing that
will be used in the HRF modeling. For better accuracy, we use very fine resolution
(noise free) HRF model, and then perform sampling at the actual acquisition timing of
each slice.
11.1.2 Multiple Reference Functions to Address Variability
Instead of using only one HRF model to be the reference function for the fMRI
detection, we shift the fine resolution model to address the brain response variability,
and then, again, perform sampling at the actual slice timings. Therefore, there are
multiple reference functions for the fMRI detection. Finally, the function response
timing of a voxel can be determined.
95
11.1.3 Making Complicated Computation Simple
Because multiple reference functions are used in the dfMRI detection, complex
and massive computation are required. However, with the compact formula (Eq. 4.1),
the dfMRI data processing becomes fast and simple if using a vector-matrix oriented
program.
11.1.4 Spatio-temporal fMRI Visualization
The results of dfMRI can be displayed as frames and time delays. Thus they
can be visualized as a movie to display the space and time relations of the brain
responses.
11.1.5 Addressing Difficulties in Certain fMRI Experiments
In some fMRI such as pre-sugical planning fMRI studies, the subject has to
respond or move her/his extremities, causing motion artifacts. With dfMRI timing
results, the motion artifacts can be identified and discarded. Hence the final fMRI
results are reliable and accurate.
11.1.6 Detecting Neuro-magnetic Signals Directly with MRI
Studies were done to optimize parameters for direct detection of the weak and
fast neuro-magnetic signals. Optimization is based on choosing the sampling
frequency and the signal of interest (stimulus) frequency to control the frequency
range of interest.
96
References
[1] Aguirre, C. K., Zarahn, E., and D’Esposito, M. (1998). The variability of
human, BOLD hemodynamic responses. NeuroImage, 8: 360-369.
[2] Al-Dayeh, L., Kim, T., Sungkarat, W., and Singh, M. (1999). Multi-frequency
reference function to reduce noise in functional MRI. IEEE Trans. Nucl. Sci.,
46(3): 513-519.
[3] Ances, B. M. (2003). Coupling of changes in cerebral blood flow with neural
activity: What must initially dip must come back up. Journal of Blood Flow &
Metabolism, 24: 1-6.
[4] Arthurs, O. J. and Boniface S. (2002). How well de we understand the neural
origins of the fMRI BOLD signal? Trends in Neuroscience, 25(1): 27-31.
[5] Bandettini, P. A., Jesmanowics A., Wong, E. C., and Hyde, J. S. (1993).
Processing strategies for time-course data sets in functional MRI of human
brain. Magnetic Resonance in Medicine, 30: 161-173.
[6] Bodurka, J., and Bandettini, P. A. (2002). Toward direct mapping of neuronal
activity: MRI detection of ultraweak, transient magnetic field changes.
Magnetic Resonance in Medicine, 47: 1052-1058.
[7] Bodurka, J., Jesmanowicz, A., Hyde, J. S., Xu, H., Estkowski, L., and Li, S.J.
(1999). Current-induced magnetic resonance phase imaging. J. Magn. Res.,
137:265-271.
[8] Brady, T. J., Goldman, M. R., Pykett, I. L., Buonanno, F. S., Kistler, J. P.,
Newhouse, J. H., Burt, C. T., Hisnshaw, W. S., and Pohost, G. M. (1982).
Proton nuclear magnetic resonance imaging of regionally ischemic canine
hearts: Effect of paramagnetic proton signal enhancement. Radiology, 144:
343-347.
[9] Buxton, R. B., Wong, E. C., and Frank, L. R. (1998). Dynamics of blood flow
and oxygenation changes during brain activation: The Balloon model.
Magnetic Resonance in Medicine, 39: 855-864.
[10] Cohen, M. S. (1997). Parametric analysis of fMRI data using linear systems
methods. NeuroImage, 6: 93-103.
97
[11] Chu, R., de-Zwart, J. A., Gelderen, P. V., Fukunaga, M., Kellman, P., Holroyd,
T., and Duyn, J. H. (2004). Hunting for neuronal currents: absence of rapid
MRI signal changes during visual-evoked response. NeuroImage, 23: 1059-
1067.
[12] Detre, J. A., Leigh, J. S., Williams, D. S., and Koretsky, A. P. (1992).
Perfusion imaging. Magnetic Resonance in Medicine, 23:37-45.
[13] Duong, T. Q., Yacoub, E., Adriany, G., Hu, X., Ugurbil, K., and Kim, S. -G.
(2003). Microvascular BOLD contribution at 4 and 7 T in the human brain:
Gradient-Echo and Spin-Echo fMRI with Suppression of Blood Effect.
Magnetic Resonance in Medicine, 49: 1019-1027.
[14] Friston, K., Fletcher, P., Josephs, O., Holmes, A., Rugg, M., and Turner, R.
(1998). Event-related fMRI: Characterizing differential response.
NeuroImage, 7(1): 30-40.
[15] Friston, K. J., Holmes, A. P., Worsley, K. J., Poline, J. P., Frith, C. D., and
Frackowiak, R. S. J. (1995). Statistical parametric maps in functional imaging:
A general linear model approach. Hum. Brain Mapp., 2: 173-181.
[16] Friston, K., Jezzard, P., and Turner, R. (1994). Analysis of functional MRI
time-series. Human Brain Mapping, 1: 153-171.
[17] Friston, K., Mechelli A., Turner, R., and Price C. J. (2000). Nonlinear
response in fMRI: The Balloon model, Volterra kernels, and other
hemodynamics. NeuroImage, 12: 466-477.
[18] Gati, J. S., Menon, R. S., Ugurbil, K., and Rutt, B. K. (1997). Experimental
determination of the BOLD field strength dependence in vessels and tissue.
Magnetic Resonance in Medicine, 38(2):296-302.
[19] Giltelman, D. R., Penny, W. D., Ashburner, J., and Friston, K.J. (2003).
Modeling regional and psychophysiologic interactions in fMRI: the important
of hemodynamic deconvolution. NeuroImage, 19: 200-2007.
[20] Glover, G. (1999). Deconvolution of impulse response in event-related BOLD
fMRI. NeuroImage, 9(4): 416-429.
[21] Goutte C., Nielsen F. A., and Hansen L. K. (2000). Modeling the
haemodynamic response in fMRI using smooth FIR filters. IEEE Trans. Med.
Imaging, Vol. 19, No. 12, 1188-1201.
98
[22] Hall, D. A., Goncalves, M. S., Smith, S., Jezzard, P., Haggard, M. P., and
Kornak, J. (2002). A method for determining venous contribution to BOLD
contrast sensory activation. Magnetic Resonance Imaging, 20: 695-706.
[23] Heeger, D. J. and Ress, D. (2002). What does fMRI tell us about neuronal
activity. Nature Reviews (Neuroscience), 3: 142-151.
[24] Hernandez, L., Badre, D., Noll, D., and Jonides, J. (2002). Temporal
Sensitivity of Even-Related fMRI. NeuroImage, 17: 1018-1026.
[25] Holz, M. and Muller, C. (1980). NMR measurements of internal magnetic
field gradients caused by the presence of an electric current in an electrolyte
solution. J. Magn. Res., 40: 595-599.
[26] Hu, X., Le, T. H., and Ugurbil, K. (1997). Evaluation of early response in
fMRI in individual subjects using short stimulus duration. Magnetic
Resonance in Medicine, 37: 877-884.
[27] Hu, X., Le, T., Parrish, T., and Erhard, P. (1995). Retrospective estimation and
correction of physiological artifacts in fMRI. Magnetic Resonance in
Medicine, 34:201-212.
[28] Josephs, O., Turner, R., and Friston, K. (1997). Event-related fMRI. Human
Brain Mapping, 5: 243-248.
[29] Joy, M., Scott, G. C., and Henkelman, R. M. (1989). In-vivo detection of
applied currents by magnetic resonance imaging. Mag. Res. Imag., 7: 89-94.
[30] Kamei, H. I. K., Yoshikawa, K., and Ueno, S. (1999). Neuronal current
distribution imaging using magnetic resonance. IEEE Trans. Magn., 35: 4109-
4111.
[31] Kwong, K. K., Belliveau, J. W., Chesler, D. A., Goldberg, I.E., Weisskoff, R.
M., Poncelet, B. P., Kennedy, D. N., Hoppel, B. E., Cohen, M. S., Turner, R.,
Cheng, H. -M., Brady, T. J., and Rosen, B. R. (1992). Dynamic magnetic
resonance imaging of human brain activity during primary sensory stimulation.
Proc. Natl. Acad. Sci. USA, 89: 5675-5679.
[32] Lange, N. and Zeger, S. L. (1997). Non-linear fourier time series analysis for
human brain mapping by functional magnetic resonance imaging. Applied
Statistics, Journal of the Royal Statistical Society, Series C, 46(1): 1-29.
[33] Lauterbur, P. (1973). Image information by induced local interactions:
Examples employing nuclear magnetic resonance. Nature, 242: 190-191.
99
[34] Lee, A. T., Glover, G. H., and Meyer, C. H. (1995). Distribution of large
veinous vessels in time-courses spiral blood-oxygen-level-dependent magnetic
resonance functional neuroimaging. Magnetic Resonance in Medicine,
33(7):745-754.
[35] Liu, H. L., Pu, Y., Nickerson, L. D., Liu, Y., Fox, P. T., and Gao, J. H. (2000).
Comparision of the temporal response in perfusion and BOLD-based event-
related functional MRI. Magnetic Resonance in Medicine, 43: 768-772.
[36] Mandeville, J. B., Marota, J. J. A., Ayata, C., Zaharchuk, G., Moskowitz, M.
A., Rosen, B. R., and Weisskoff, R. M. (1999). Evidence of a cerebrovascular
post-arteriole windkessel with delayed compliance. J. Cereb. Blood Flow
Metab., 19: 679-689.
[37] Mandeville, J. B. and Rosen, B. R. (2002). Functional MRI. In Toga, A. W.
and Mazziotta, J. C. (Eds.), Brain Mapping The Methods, 2/e, pp. 315-349.
Academic Press.
[38] Mansfield, P. and Grannell, P. K. (1973). NMR `diffraction' in solids?
J. Phys., C6, L422-L426.
[39] MATLAB is a software by MathWorks, Cambridge, Massachusetts, USA
[40] Menon, R.S., Luknowski, D.C., and Gati, J.S. (1998). Mental Chronometry
using latency-resolved functional MRI. Proc. Natl. Acad. Sci. USA, 95:
10902-10907.
[41] Nakai, T., Matsuo, K., Kato, C., Takehara, Y., Isoda, H., Moriya, T., Okada,
T., and Sakahara, H. (2000). Post-stimulus response in hemodynamics
observed by functional magnetic resonance imaging—difference between the
primary sensorimotor area and the supplementary motor area. Magnetic
Resonance Imaging, 18: 1215-1219.
[42] Ogawa, S., Lee, T. M., Kay, A. R., and Tank, D. W. (1990). Brain magnetic
resonance imaging with contrast dependent on blood oxygenation. Proc. Natl.
Acad. Sci. USA, 87: 9868-9872.
[43] Price, C. J., Veltman, D. J., Ashburner, J., Josephs, O., and Friston, K. J.
(1999). The critical relationship between the timing of the stimulus
presentation and data acquisition in blocked designs with fMRI. NeuroImage,
10: 36-44.
[44] Saad, Z. S., Ropella, K. M., DeYoe, E. A., and Bandettini, P.A. (2003). The
spatial extent of BOLD response. NeuroImage, 19: 132-144.
100
[45] Scott, G. C., Roy, M., Armstrong, R., and Henkelman, R. M. (1991).
Measurement of nonuniform current density by magnetic resonance. IEEE
Trans. Med. Imag., 10(3): 362-373.
[46] Silva, A. C., and Kim, S. G. (1999). Pseudo-continuous arterial spin labling
technique for measuring CBF dynamics with high temporal resolution.
Magnetic Resonance in Medicine, 42: 425-429.
[47] Silva, A., Lee, S-P., Yang, C., Iadecola, C., and Kim, S-G. (1999).
Simultaneous BOLD and perfusion functional MRI during forepaw stimulation
in rats. J. Cereb. Blood Flow Metab., 19: 871-879.
[48] Singh, M. (1991). Feasibility of in-vivo neuromagnetic field detection using
MR phase shifts. Proc. Soc. Magnetic Resonance in Medicine, Vol. 2: 837.
[49] Singh, M. (1994). Sensitivity of MR phase shift to detect evoked
neuromagnetic fields inside the head. IEEE Trans. Nucl. Sci., 41(1): 349-351.
[50] Singh, M., Doria, D., Henderson, V. W., Huth, G. C., and Beatty, J. (1984).
Reconstruction of images from neuromagnetic fields. IEEE Trans. Nucl. Sci.,
NS-31(1).
[51] Singh, M., Kim, S., and Kim, T-S. (2003). Correlation between BOLD-fMRI
and EEG signal changes in response to visual stimulus frequency in humans.
Magnetic Resonance in Medicine, 49: 108-114.
[52] Singh, M., Patel, P., Khosla, D., and Kim, T-S. (1996). Segmentation of
functional MRI by k-means clustering. IEEE Trans. Nucl. Sci., 43(3): 2030-
2036.
[53] Singh, M., and Sungkarat, W. (2005). A novel approach to image neural
activity directly by MRI. Proc. SPIE, 5746: 109-118.
[54] Singh, M., Sungkarat, W., Jeong, J., and Zhou, Y. (2002). Extraction of
temporal information in functional MRI. IEEE Trans. Nucl. Sci., 49(5): 2284-
2290.
[55] Smith, S. M. (2001). Preparing fMRI data for statistical analysis. In Jezzard,
P., Matthews, P. M., and Smith, S. M. (Eds.), Functional MRI: An
Introduction to Methods, pp. 229-241. New York: Oxford Univ. Press Inc.
[56] Veltman, D. and Hutton, C. (2001). SPM99 manual.
http://www.fil.ion.ucl.ac.uk/spm/course/manual/manual.pdf
101
[57] Wildgruber, D., Erb, M., Klose, U., and Grodd, W. (1997). Sequential
activation of supplementary motor area and primary motor cortex during self-
paced finger movement in human evaluated by functional MRI. Neuroscience
Letters, 227: 161-164.
[58] Woolsey, T. A., Rovainen, C. M., Cox, S.B., Henegar, M. H., Liang, G.E.,
Liu, D., Moskalenko, Y.E., Sui, J., and Wei, L. (1996). Neuronal unit linked
to microvascular modules in cerebral cortex: Response elements for imaging
the brain. Cereb. Cortex, 6: 647-660.
[59] Worsley, K. J., Liao, C. H., Aston, J., Petre, V., Duncan, G. H., Morales, F.,
and Evans, A. C. (2002). A general statistical analysis for fMRI data.
NeuroImage, 15: 1-15.
[60] Xiong, J., Fox, P. T., and Gao, J. H. (2003). Directly mapping magnetic field
effects of neuronal activity by magnetic resonance imaging. Human Brain
Mapping, 20: 41-49.
[61] Yacoub, E. and Hu, X. (2001). Detection of the early decrease in fMRI signal
in the motor area. Magnetic Resonance in Medicine, 45: 184-190.
[62] Yacoub, E. and Hu, X. (1999). Detection of the early negative response in
fMRI at 1.5 Tesla. Magnetic Resonance in Medicine, 41: 1088-1092.
[63] Yacoub, E., Le, T. H., Ugurbil, K., and Hu, X. (1999). Further evaluation of
the initial negative response in functional magnetic resonance imaging.
Magnetic Resonance in Medicine, 41: 436-441.
[64] Yacoub, E., Shmuel, A., Pfeuffer, J., Moortele, P-F. V., Adriany, G., Ugurbil,
K., and Hu, X. (2001). Investigation of the initial dip in fMRI at 7 Tesla.
NMR Biomed., 14: 408-412.
Abstract (if available)
Abstract
An important challenge in functional magnetic resonance imaging (fMRI) is to achieve the most spatially accurate results, i.e., to localize activation as close as possible to the actual site of working neurons. Because fMRI detection methods are based on the blood oxygen level dependent (BOLD) property of brain microvascular and venous systems, any model-based fMRI detection method will provide inaccurate results if the model is biased to the BOLD response of a draining venous system, e.g., a vein. Here, a novel detection method is proposed to achieve accurate localization by measuring the spatiotemporal dynamics of brain activations in fMRI to separate microvasculature from venous activations. The method utilizes precisely down-sampled multiple models (with positive and negative delays). Thus conventional 3D fMRI results are expanded into 4D (with an additional temporal dimension). Example fMRI studies underline the usefulness of the approach.
Linked assets
University of Southern California Dissertations and Theses
Conceptually similar
PDF
Theory and simulation of diffusion magnetic resonance imaging on brain's white matter
PDF
Functional magnetic resonance imaging characterization of peripheral form vision
PDF
Functional models of fMRI BOLD signal in the visual cortex
PDF
Unveiling the white matter microstructure in 22q11.2 deletion syndrome with diffusion magnetic resonance imaging
PDF
Statistical signal processing of magnetoencephalography data
PDF
Fast flexible dynamic three-dimensional magnetic resonance imaging
PDF
Applications of graph theory to brain connectivity analysis
PDF
Correction, coregistration and connectivity analysis of multi-contrast brain MRI
PDF
Optimization methods and algorithms for constrained magnetic resonance imaging
PDF
High-dimensional magnetic resonance imaging of microstructure
PDF
New theory and methods for accelerated MRI reconstruction
PDF
Estimation of cognitive brain activity in sickle cell disease using functional near-infrared spectroscopy and dynamic systems modeling
PDF
Localization of multiple targets in multi-path environnents
PDF
Signal processing methods for interaction analysis of functional brain imaging data
PDF
Functional connectivity analysis and network identification in the human brain
PDF
Zero-power sensing and processing with piezoelectric resonators
PDF
Quantitative MRI for oxygenation imaging in cerebrovascular diseases
PDF
Switching dynamical systems with Poisson and multiscale observations with applications to neural population activity
PDF
Statistical inference for dynamical, interacting multi-object systems with emphasis on human small group interactions
PDF
3D vessel mapping techniques for retina and brain as an early imaging biomarker for small vessel diseases
Asset Metadata
Creator
Sungkarat, Witaya
(author)
Core Title
Dynamic functional magnetic resonance imaging to localize activated neurons
School
Viterbi School of Engineering
Degree
Doctor of Philosophy
Degree Program
Biomedical Engineering
Defense Date
12/06/2006
Publisher
University of Southern California
(original),
University of Southern California. Libraries
(digital)
Tag
aliasing,dfMRI,drect neuromagnetic signal detection,fMRI,magnetic source,msMRI,neuroelectric signal,nfMRI,OAI-PMH Harvest,spatiotemporal dynamics
Language
English
Advisor
Singh, Manbir (
committee chair
), Wolf, Walter (
committee member
), Yen, Jesse T. (
committee member
)
Creator Email
sungkara@usc.edu
Permanent Link (DOI)
https://doi.org/10.25549/usctheses-m274
Unique identifier
UC1315936
Identifier
etd-Sungkarat-20070215 (filename),usctheses-m40 (legacy collection record id),usctheses-c127-161474 (legacy record id),usctheses-m274 (legacy record id)
Legacy Identifier
etd-Sungkarat-20070215.pdf
Dmrecord
161474
Document Type
Dissertation
Rights
Sungkarat, Witaya
Type
texts
Source
University of Southern California
(contributing entity),
University of Southern California Dissertations and Theses
(collection)
Repository Name
Libraries, University of Southern California
Repository Location
Los Angeles, California
Repository Email
cisadmin@lib.usc.edu
Tags
aliasing
dfMRI
drect neuromagnetic signal detection
fMRI
magnetic source
msMRI
neuroelectric signal
nfMRI
spatiotemporal dynamics