practiCal fMRI: the nuts & bolts

New user training guide/FAQ

2012-03-08T18:00:00.000-08:00

I've just uploaded a new user training guide/FAQ that we use at Berkeley to initiate newbies into the ways of the dark side. It is Siemens-specific, for a Trio/TIM.

As last time, the guide is a bit rough. Sorry for English-isms and typos. It's worth exactly what you pay for it. It's free. Use and abuse it however you like. It's a Word document so that you can reorder things, add your own notes, etc. I would appreciate constructive feedback, especially if you find mistakes or have suggestions to improve it, but there's no need to ask permission to use it, change it, replicate it, sell it...

The most recent version of the training guide/FAQ is available from this web page:

http://bic.berkeley.edu/scanning

Locate the file attachment towards the bottom of the page, it's called 3T_user_training_FAQ_08Mar2012.doc. The most recent contents and a list of changes since the last version (April, 2011) appear below.

Caveat emptor.

The document is only a component of user training, don't expect to learn how to scan by reading it! Rather, use the tips to extend your understanding, refine your experimental technique and so on. Note also that this document is for a Siemens TIM/Trio (with 32 receive channels) and running software VB17. There may be subtle or not-so-subtle differences for the Verio and Skyra platforms, for software VB15, VD11, etc. so keep your wits about you if you're not on a Trio with VB17!

You may have local differences, e.g. custom pulse sequences, that allow you to do things that contradict what you find in this user guide. Talk to your physicist and your local user group before taking anything you find in this guide/FAQ too literally.

Finally, you wont find many (any?) references in this guide/FAQ. It's for the training of newbies, not a comprehensive literature review! If you are seeking further information on something I mention in the guide and you can't find a suitable reference yourself, shoot me an email and I'll do my best to point you in a useful direction.

--------------------------------------

Update Notes (8^th March, 2012):

Updated with new operating modes available under software syngo MR version B17.

General tweaks to improve readability.

Further recommendations on using the 32-channel coil for fMRI.

Added a description of the new AutoAlign procedure, AAHScout.

Added a new section: “I have an existing protocol that uses the old AutoAlign (AAScout). How do I get and use the new AutoAlign (AAHScout)?”

Added a new section: “I want to add a new acquisition and acquire exactly the same slices as this other EPI acquisition I just acquired. How do I tell the scanner to do that?”

Extended the discussion on the relative merits of PACE versus using an offline realignment alone, in the section on the ep2d_pace sequence.

Fixed a typo concerning the slice ordering for descending slices.

Added a new section: “What is a field map and how does it fix EPI distortion?”

Added a new section: “I want to try to fix my distortion with a field map. What do I need to acquire?”

Updated the sections on partial Fourier for EPI, noting that Siemens simply zero fills the omitted portion of k-space rather than doing a conjugate synthesis.

Extended checklists.

-------------------------------------

CONTENTS (as of 8th March, 2012):

Sections added or modified since last version are highlighted in yellow.

SETTING UP AND ACQUIRING SCANS:

What is the practical difference between the 12-channel and 32-channel head coils? Which one is best for fMRI?

I have a subject who has a lot of dental work. Is this person okay to scan?

Why does the scanner instruct me that the patient bed might move when I start the first scan in my session (usually a localizer)?

I can’t hear anything happening? How can I tell what the scanner is doing right now?

Why do I sometimes get a message that the subject might experience peripheral nerve stimulation? Should I tell the subject?

How does the AutoAlign feature work? Should I use it?

NEW: I have an existing protocol that uses the old AutoAlign (AAScout). How do I get and use the new AutoAlign (AAHScout) instead?

I don’t want to trust AutoAlign. How should I define my slice positions manually?

NEW: I want to add a new acquisition and acquire exactly the same slices as this other EPI acquisition I just acquired. How do I tell the scanner to do that?

When does shimming happen and what is actually done?

I want to re-shim my subject’s brain midway through my session. How do I do it?

How do I know whether I should re-shim or not?

I want to know how long my scan will take. Where is the scan time shown?

What is the difference between the Scan and Apply buttons for starting a scan?

Help! What pulse sequence am I using?

EPI: BASIC PARAMETER AND SEQUENCE ISSUES

I’ve been told not to use echo spacing between 0.6 and 0.8 ms for EPI. How come?

How many dummy scans happen before the first real (saved) volume of EPI in my time series?

I want 200 volumes in my EPI time series. How do I do that?

On the BOLD card, what is Motion Correction? How do I turn it on or off?

My protocol has TE set at 28 ms for EPI. But I saw somebody else’s protocol that uses a TE of 22 ms. How come?

I am using ep2d_bold. What are the specifics of using this sequence?

I am using ep2d_pace. What are the specifics of using this sequence?

I am using ep2d_neuro. What are the specifics of using this sequence?

What flip angle should I use for fMRI?

What TR should I use for fMRI?

Should I use interleaved or sequential slices for fMRI?

In what order does the scanner acquire EPI slices?

EPI: ARTIFACTS

I hear a lot about ghosting when people talk about EPI. What is a ghost and what causes them? How do I get rid of them?

On the Contrast tab I notice that fat suppression is enabled for EPI. What does it do?

What is the origin of signal dropout in EPI? Can it be fixed?

What is the origin of distortion in EPI? Can it be fixed?

NEW: What is a field map and how does it fix EPI distortion?

NEW: I want to try to fix my distortion with a field map. What do I need to acquire?

Whoa! I’m watching my EPIs on the Inline Display window and I’m seeing all sorts of weirdness. What’s going wrong?

How much subject movement is too much?

EPI: ADVANCED PARAMETER AND SEQUENCE ISSUES

What the hell is iPAT? Last time I checked, grappa was a strong Italian drink! It makes no sense!

Is GRAPPA a good technique to use? What are the caveats?

What is “partial Fourier” and why might I want to consider it for EPI?

Is partial Fourier a good technique to use? What are the caveats?

It looks like I will need to use either partial Fourier or GRAPPA to get the spatial resolution and coverage that I want. Which method should I use?

FINAL ISSUES:

I want to scan overnight. Is there anything I need to watch out for?

I hear we have a research agreement with Siemens. Why should I care?

APPENDIX 1: CHECKLISTS

Normal operation checklists:

Experimenter prep
Lab prep
Subject prep
Subject setup
Start of scan
Experimental protocol
End of scan

Emergency checklists:

Unexpected image feature
Panicked subject
Magnetic object accident
Fire
Earthquake

Common persistent EPI artifacts: Distortion and dropout

2012-02-28T09:24:00.000-08:00

The origins of distortion and dropout in EPI were covered in PFUFA Part Twelve, and both of these artifacts have been mentioned in passing in the previous articles concerning abnormally high ghosting. In some instances these artifacts are "co-morbid" because certain issues that cause abnormally high ghosting - such as a poor shim because of asymmetric placement of the subject's head in the magnet - are likely to increase distortion and dropout effects at the same time. Except that it can be very difficult to evaluate distortion and dropout by inspection, during an experiment. The ghosts can be used as a fairly independent "barometer" of the experiment's quality if, as is often the case, some of them fall into an image region that is otherwise noise. Not so with distortion and dropout. By definition these artifacts plague signal regions in the brain, and even an experienced operator can have a tough time determining when either issue is worse than it might otherwise be.

So I'm afraid I don't have a whole lot of new information to offer on either distortion or dropout, from the perspective of diagnosing and potentially changing (improving) your experiment on the day. Other than very obvious deficiencies, as might happen if the subject has a highly conductive hair product, for example, I don't spend much time evaluating distortion or dropout by inspection. Ghosts can be a good surrogate for all that ails distortion and dropout, so I focus on those.

Where you can potentially improve the situation for distortion and dropout is with parameter selection when you are establishing your experimental protocol. Distortion and dropout will generally change with slice prescription, as we already saw in the "good data" posts. And it may be that reduction of dropout leads you to use a particular slice direction, e.g. coronal slices for improved frontal lobe signal. After that, the other common tactics to minimize dropout are to use the thinnest possible slice thickness, possibly using higher in-plane spatial resolution, and perhaps decrease TE. These are protocol/parameter questions that are covered somewhat in my user training guide/FAQ, and I will expand on those sections below. Be warned, however, that it is very difficult to provide general guidelines for all fMRI experiments. Instead, the parameter choices tend to be dictated by your primary requirements. You might select very different parameters for a study that is primarily interested in orbitofrontal cortex than you would use for a sensorimotor task. It's horses for courses.

Approaches to tackling distortion

The level of distortion in the phase encoding dimension is a function of the echo spacing, as explained in PFUFA Part Twelve. Tactics to reduce the distortion level involve making fundamental changes to the phase encoding k-space scheme, e.g. multi-shot segmented k-space, or parallel imaging methods. In each approach the essential idea is to increase the k-space step size, thereby increasing the bandwidth of the phase encoding dimension.

Distortion reduction techniques such as segmented EPI and GRAPPA are beyond the scope of this artifact series. For now I am focusing exclusively on single-shot EPI because that's what most people use. But I would note in passing that each of these distortion reduction schemes comes at a cost - usually increased motion sensitivity. I already discussed the motion sensitivity of GRAPPA in a previous post. Segmented EPI gets very little use for fMRI at 3 T these days, but at some point I may do a PFUFA post on it. Until then, Stuart Clare's PhD thesis is online and has a good description of the method as well as some of the issues.

That leaves another approach entirely: trying to fix the distortion. At this juncture the most common way to attempt a fix is with a magnetic field map. I'm not going to get into the specifics of field maps today because I see this as a peripheral issue to the theme of this post series, which is on the acquisition of the EPI data. You can find a basic introduction to the acquisition and use of a field map in my most recent user training guide/FAQ (see Note 1). Your choice to use a field map to try to correct distortion has no direct influence on what you do in the EPI; the parameters you select will produce a certain level of distortion, and your decision on whether to try to remedy the distortion is one that can be made independently. That said, I may try to do a dedicated post on distortion correction at some point, because there are some experimental considerations, such as when to acquire a field map for distortion correction, and whether more than one field map acquisition is required per scan session.

Until that dedicated post I will direct you to a short review article, and references therein, courtesy of Peter Jezzard:

"Correction of geometric distortion in fMRI data." P. Jezzard, NeuroImage Epub (2011).

I would point out the youthful Peter Jezzard who appears in Figure 1, except that the bugger looks exactly the same today!

The paper doesn't address worse than normal/typical distortion, it's about the problem and how to potentially fix it. The distortion you get in your experiment should only be a function of the parameters you've selected, and the shim (i.e. the residual magnetic susceptibility gradients). And you have already seen how to diagnose a poor shim using the N/2 ghosts.

Advanced methods to try to undo distortion are (at this juncture) beyond the scope of this series because it's not something that is likely to be a function of your skill as an experimenter. Furthermore, few of the methods are (yet?) available on commercial scanners. Do a Pubmed search for "epi distortion correction" to pull up a few dozen methods. It's interesting to note that most of the methods are aimed at fixing distortion for EPI-based diffusion imaging. This may well be because people doing tractography are more aware of the need to match to undistorted anatomical space, but it could also be because changing the pulse sequence for diffusion imaging doesn't have the same statistical implications as it does for fMRI. In fMRI we seek to make each sample - essentially, each TR period - as independent as possible. Many schemes aimed at reducing distortion cause a temporal (and possibly spatial) smoothing to be imposed on the time series. But, as I say, this is far beyond the scope of this post or this blog right now.

Approaches to tackling dropout

There are modified pulse sequences, such as Z shim methods (do a Pubmed search for "EPI z-shim") that seek to use regionally optimal refocusing gradients and then pool the results into a final, improved image. None of these sequences has gained widespread acceptance yet, I suspect because there tends to be a temporal penalty in the acquisition of different pieces of the final image, and because there is additional complexity in the data processing (even if the pulse sequences were available on commercial scanners, which most aren't). Instead, there are techniques to tease out improved performance with a series of tweaks to TE and voxel size, such as this one from the FIL:

"Optimized EPI for fMRI studies of the orbitofrontal cortex: compensation of susceptibility-induced gradients in the readout direction." Weiskopf et al., MAGMA 20(1), 39-49 (2007).

The general principles are easy to comprehend: lower TE or higher spatial resolution (smaller voxels) tend to reduce the susceptibility-induced dephasing, leading to signal recovery in some brain regions. The costs of these approaches should be relatively obvious by now, too. A shorter TE will sacrifice some BOLD sensitivity in well-shimmed brain regions - but probably not enough to be of concern - while attaining higher spatial resolution requires higher gradient performance and acquisition time, potentially limiting coverage to a portion of the brain. However, in decreasing the signal loss through higher resolution - a process that increases the effective T2* for the problem voxels - it can then be useful to increase the TE to regain optimum BOLD sensitivity. This was the conclusion of a study that investigated optimum parameters for fMRI of amygdala:

"The impact of EPI voxel size on SNR and BOLD sensitivity in the anterior medio-temporal lobe: a comparative group study of deactivation of the default mode." Robinson et al., MAGMA 21(4), 279-90 (2008).

There are several other papers concerning sets of tweaks for specific applications; the above paper covers many of them. But I'm going to quit here because this post is already rambling. Feel free to post a question on particular brain regions or applications and I'll do my best to track down prior work for you to follow.

_________________

Notes:

1. I'll be releasing an updated user training guide/FAQ in the next week. The acquisition of a field map for distortion correction is one of the new sections.

As previously, the most recent version of the training guide/FAQ is available from this web page:

http://bic.berkeley.edu/scanning

Locate the file attachment towards the bottom of the page. The last version was called 3T_user_training_FAQ_19April2011.doc.

Terminology change: characterizing EPI artifacts

2012-02-19T10:43:00.000-08:00

After considering the rest of the topics that I want to cover in the series of posts on EPI artifacts, I've decided to change the terms "static" and "dynamic" to "persistent" and "intermittent," respectively. I think the new terms better reflect the dominant temporal character of each artifact. The idea is to sort them based on whether they are likely to plague every frame of an EPI time series, or come and go.

Take external RF interference, for example. Say you fail to close properly the RF-sealed door to the magnet room and the door opens slightly during your scan, leading to RF contamination from "environmental" sources. In this instance the RF interference itself isn't likely to be static - it will vary with whatever sources of RF happen to be in your scanner environment - but it will persist (at some level) until the door is closed. With a good eye or some appropriate diagnostics it would be possible to show the persistence of the problem throughout an acquisition. Contrast this situation with a static electrical discharge somewhere within the scanner room; tiny sparks that cause a broad range of electromagnetic frequencies, including radiofrequencies. These can arise if the humidity of the magnet room air becomes too low. Depending on the source of the sparks, the humidity, etc. you might find that only one or two TRs of a time series are contaminated, or you could find the entire time series is affected. I will therefore characterize static electrical discharges as an intermittent artifact.

Pedants will spot that the first example can be modulated by the position of the magnet room door, rendering the artifact intermittent, while in the second example the propensity for static electrical discharges will persist as long as the source exists, while the humidity remains low, etc. So, yeah, in some ways the distinctions I'm making are subjective. FMRI, like life, is complicated! Still, I'm hopeful that a more practical characterization of artifacts will assist you in differentiating and diagnosing them when it matters: during your experiment. And once you're an expert you will find it easy to comprehend the nuances of temporal behavior, when my artificial distinctions will be all but irrelevant to you.

So, there you have it. I'll be going back to edit the existing posts in this series over the next couple of days. Apologies for any confusion the switch creates.

Physics for understanding fMRI artifacts: CONTENTS

2012-02-16T18:21:00.001-08:00

Figured it might be useful to have some summary/contents pages. I'll do similar admin posts as the other series mature, too. And I will label these pages with "Contents" to make them easier to find via the sidebar.

Part One

An introduction to the series, followed by an introductory video courtesy of Sir Paul Callaghan: What is NMR and how does it work?

Part Two

Further videos explaining the principles of nuclear magnetic resonance - how the intrinsic spin of certain atomic nuclei interacts with applied magnetic fields to yield useful information.

Part Three

Videos showing the anatomy of a miniature scanner, a basic NMR experiment, why shimming is important for NMR (and MRI), how and why a spin echo works, and the relaxation of spins back towards their ground state.

Part Four

Mathematics of oscillations: an introduction to imaginary and complex numbers, and frequency and phase.

Part Five

An introduction to the Fourier transform - what it does and how it works. Includes a description of Fourier analysis, and explains conjugate variables and (Fourier) domains.

Part Six

Practical issues arising from the use of the Fourier transform in MRI:

Fourier pairs
Convolution in pictures
Time-decaying signals
Finite acquisition periods and signal clipping
Digitization
The Nyquist frequency
Aliasing

Part Seven

Magnetic field gradients and one-dimensional MRI:

A magnetic field alters the local resonance frequency
Acquiring an MR signal in the presence of a gradient
One-dimensional imaging in pictures

Part Eight

Gradient-recalled echoes:

An explanation of how GRE works
The benefits of acquiring a gradient echo

Slice selection:

How slice selection works
Using a GRE with slice selection

Part Nine

K-space - conjugate variables redefined:

Conjugate variables revisited
Representing pictures in reciprocal space
One-dimensional imaging as seen in k-space
Tracing kx through time
Getting off axis (into 2D)

Part Ten

K-space in two dimensions:

A useful pictorial representation of imaging pulse sequences
The goal revisited
Gradients along the x direction (again)
Gradients along the y axis
The equivalence of frequency and phase encoding
Gaining an intuitive understanding of phase encoding

Part Eleven

Resolution and field-of-view as seen in k-space:

Spatial frequencies in k-space: what lives where?
Why does the signal level change across k-space?
Defining parameters in k-space to yield the image you want
Image field-of-view
Image resolution

Part Twelve

The echo planar imaging (EPI) pulse sequence:

Gradient echo EPI pulse sequence
Processing 2D k-space for EPI

EPI artifacts:

Ghosting
Distortion
Signal dropout

(Part Thirteen)

To come: ramp sampling.

(Part Fourteen)

To come: partial Fourier.

Common persistent EPI artifacts: Abnormally high N/2 ghosts (2/2)

2012-02-15T08:43:00.005-08:00

In the previous post I covered sources of persistent ghosts that arise as a result of some property of the subject, such as the orientation of the subject's head in the magnet. These are what I'm categorizing as subject-dependent effects. In this post I will review the most common sources of persistent ghosts attributable to the scanner, either from an intrinsic property that you might encounter inadvertently, or from mis-setting a parameter in your protocol. As I mentioned last time, I am restricting the discussion to factors that you have some control over as the scanner operator. Ghosts that arise because of a scanner installation error, such as poor gradient eddy current compensation or inaccurate gradient calibration, are issues for your facility physicist and/or your service engineer.

Scanner-dependent conditions:

Rotated read/phase encode axes

GLOBAL - affects all slices to some extent.

This is an insidious problem that we could categorize as pilot error, except that it's very easily encountered without realizing it. When you set up your slice prescription you are primarily concerned with capturing all those brain regions you need for your experiment. Or you might be concerned with setting a particular slice angle relative to the brain anatomy, e.g. parallel to AC-PC. Now, if the subject's head is precisely aligned such that the read and phase encode axes of your imaging plane are matched perfectly with the gradient set axes (i.e. with the magnet's frame of reference), then for axial slices the readout dimension will be attained using pure X gradient (subject's left-right) while the phase encode dimension uses pure Y gradient (subject's anterior-posterior). (See Note 4 in the post on "Good" coronal and sagittal data for an explanation of why the gradients are established this way, for subject safety/comfort reasons.) But, if the head is twisted slightly, or you're a little sloppy with your slice positioning, then it is quite easy to have a readout gradient that is mostly X with a little bit of Y, and a phase encoding gradient that is mostly Y with a little bit of X. This in-plane rotation ought not be a problem if the X and Y gradients performed equivalently, but they're only similar and not identical. There tend to be small differences in the response time of the gradients, which means that when the scanner tries to drive the read gradient to its desired k-space trajectory, one component (say the X component) can respond faster than the other. This produces a slight mismatch between the target (ideal) k-space trajectory and the trajectory that's actually achieved by the gradients, thereby leading to a source of zigzags that will produce N/2 ghosting.

Now the good news. You've got to rotate the image plane by quite a lot before the ghosting starts to become apparent. It's common to have rotations of 1-2 degrees and these will generate almost no additional ghosting. Once the rotation gets much larger than 5 degrees (depending on the specifics of your scanner) then you might start to see additional ghosting. Below on the left is an ideal prescription, while on the right I've intentionally rotated the image plane by 8 degrees, leading to a small but noticeable increase in ghost level:

(Click to enlarge.)

Note that the subject's head is properly aligned in the magnet in both conditions; I didn't change the head position at all. The ghosting on the right arises purely as a result of rotating the image plane. Thus, were the subject's head rotated in the magnet, you might decide to "straighten things out" by rotating the image plane in the opposite direction to yield images that appear like those on the left, but that are actually a result of offsetting rotations - a rotation of the subject's head one way and an equal but opposite rotation of the image plane the other. Both the rotation of the subject's head (see last post) and the rotation of the image plane are likely to lead (through separate mechanisms) to increased ghosting.

As a general rule I would advocate setting the in-plane rotation to zero, even if the subject's head is significantly rotated (about the subject's longitudinal axis). This will at least avoid one source of ghosts although, since the shim will already be degraded by such a rotation (relative to the magnet/gradient axes), you're usually better off repositioning the subject and negating any need to rotate the image plane whatsoever. (Siemens users: see Note 1 on how to set the in-plane rotation to zero by hand.) But if you must be sloppy and you decide you're not going to reposition the subject's head, don't compound one problem with another.

No fat suppression

GLOBAL - all slices containing scalp fat will be affected.

I covered the need for fat suppression in my user training guide/FAQ, which is available via this post. See the section entitled "On the Contrast tab I notice that fat suppression is enabled for EPI. What does it do?" (Note: I will be uploading an updated version of the guide/FAQ soon. The current one dates from July, 2010.)

What I didn't show in that document was what happens if you decide to disable the fat suppression. Below, on the right, is an example of what happens if you disable the fat suppression compared to leaving it enabled (left):

(Click to enlarge.)

Disabling fat suppression (aka "fat sat") leads to intense ghosts from the subcutaneous lipid. These ghosts overlap the brain in just about every slice. The variance in the overlapped regions can be expected to be significantly higher than other regions, leading to serious complications of interpretation in fMRI.

What, you don't think these ghosts look so bad, and you're tempted to run without fat sat? Hmmm, maybe I should have mentioned that this subject was a fairly serious cyclist with body fat in the single digits. He also happens to have just about the most spherical brain I've ever scanned. It's darn hard to get ghosts of any type on this guy! Your subjects will show much higher ghosting, I guarantee it.

On some scanners there is an alternative to using fat suppression. Instead, the RF excitation pulse is designed to be what's called a spatial-spectral pulse. It aims to excite water in a slice while simultaneously avoiding excitation of fat resonances. GE scanners use the this approach by default, I'm told. (Siemens users, see Note 2.)

There are naturally pros and cons to the avoidance (spatial-spectral excitation RF pulses) or presaturation (fat sat) strategies. Presaturation works fairly well and permits the slice excitation to be implemented separately, and possibly optimally, but at the temporal cost of some 15 ms per slice. Let's do the arithmetic. That's an overhead of nearly half a second for thirty slices. Ouch! One might be tempted to disable the fat suppression in order to acquire more slices per TR, but if you do then ghosts of the type shown here will be the result. The principal disadvantage of the spatial-spectral approach is that the slice profile tends to be broader - less rectangular, more trapezoidal - and the minimum attainable slice thickness may be higher than for the corresponding fat suppression option. My advice is to use what's optimized on your scanner platform, which appears to be fat suppression for Siemens and spatial-spectral excitation for GE. I don't know what Philips offers.

Mechanical resonances

GLOBAL - affects all slices equally.

I covered the basics of this problem in the user training guide/FAQ, too. See the section, "I’ve been told not to use echo spacing between 0.6 and 0.8 ms for EPI. How come?" (Note: I will be uploading an updated version of the guide/FAQ soon. The current one dates from July, 2010.)

In the last version of the guide/FAQ document I didn't given examples of mechanical resonances, so here is a video showing the temporal stability of ghosts arising from mechanical resonance in the X gradient (the default readout gradient for axial slices):

It's quite remarkable how temporally stable the ghosts are. Vibration of the gradient coil is in a steady state - mechanical resonance - and there is a constant mismatch between the desired and actual k-space trajectories, leading to sufficiently large offsets between odd and even lines of phase-encoded k-space that the ghost correction scheme, using three navigator echoes, isn't well matched to it.

On Siemens Trio scanners the mechanical resonances are really only a problem for axial or axial oblique slices - prescriptions that use the X gradient for readout - and they primarily affect echo spacings in a range of about 0.6-0.8 ms. Most facilities are aware of the problem and have characterized specific "no go" zones for that particular scanner, setting up "default" protocols that avoid the issue entirely. If you do decide to go "off piste," e.g. to try to push a few extra slices in your TR, you would want to talk to your facility techs/physicists to run some quick empirical tests and check that you haven't stumbled inadvertently into a higher ghost level than is possible by using a different echo spacing. I'll leave the issue here for today because it's bordering on an issue that's one for your physicist.

Excessive ramp sampling

GLOBAL - affects all slices equally.

Disclaimer: Let me begin by stating that I haven't fully investigated this phenomenon, I'm simply going to report what I see on my scanner and do my best to explain how it looks in data. Therefore, please don't take this at anything but face value. My diagnosis may be incorrect in any number of ways. For example, it is entirely possible (even probable) that the effect I've termed "excessive ramp sampling" is actually another manifestation of the aforementioned mechanical resonances, i.e. a k-space instability generated by a mechanical resonance that is modulated in some fashion by ramp sampling, when the extent of ramp sampling fraction becomes large. (See Note 3.) I've not had the time (or, to be honest, the imperative) to investigate further than what you'll see here. I just want to be complete and make sure you can avoid any similar problems that might arise on your scanner.

This issue is another one that will likely require your physicist's input to avoid properly on your scanner; it's in the gray area between being under your control and being a system/scanner feature. But for the sake of your education I'm going to give you an example and describe one way the problem can be encountered somewhat inadvertently, if you are setting up a new EPI protocol.

Let's look at an example of the ghosting that arises when I set the echo spacing time to be the shortest possible (along with the highest readout gradient bandwidth), commensurate with a fixed TE of 28 ms and a 64x64 matrix. On the left is the echo spacing set at 0.5 ms using a readout gradient bandwidth of 2298 Hz/pixel, on the right is the echo spacing set at the minimum attainable echo spacing 0.43 ms, with readout gradient bandwidth at 2790 Hz/pixel:

(Click to enlarge.)

I don't think I need to point out the ghosts on the right, eh? And at first glance - assessing a single TR - the ghosts appear consistent with the mechanical resonance ghosts in the previous section. Except that, unlike the ghosts arising purely from a mechanical resonance, these ghosts are temporally unstable:

The ghosts appear to swirl in a periodic fashion. Any small perturbation in the power supply - or any other slight variation in the scanner's stability - causes a sizable effect in the ghost; the correction algorithm (using the three navigator echoes) produces a slightly different result for each TR.

Of course, in this instance the temporal properties are arguably moot. You don't want to be running with such large ghosts at all, whether the bulk of the ghost is temporally stable or not! Still, let's consider the temporal stability independent of the large magnitude of the ghosts, because there might be instances when the ghosts are much less evident - lower magnitude - but still have a temporal instability arising from "excessive ramp sampling."

Why use quotation marks? Because of my initial disclaimer. I think the effect arises because the degree of ramp sampling is high; in this case roughly half of the readout data points are being acquired on the readout gradient ramps. Thus, small instabilities in either the readout gradient, or mismatches in timing between the digitizer and the readout gradient waveform, can lead to temporal changes in the ghosting level from TR to TR. On my scanner I think that the source is the relatively poor electrical power - I don't have a power conditioner for the scanner, so my electrical stability is whatever I get fed from campus - leading to a small but noticeable instability in the readout gradient amplitude whenever the readout gradient is driven extremely hard. The use of a high degree of ramp sampling seems to exacerbate the gradient instability in the ghost level. But, as I said at the outset, this is just a tentative, working hypothesis. It may be something else entirely.

On my scanner I can quite easily avoid the worst effects of what I'm calling "excessive ramp sampling" by applying a simple rule of thumb: I use an echo spacing that is not less than 0.02 ms longer than the minimum attainable echo spacing. This means I fix all other acquisition parameters except echo spacing and bandwidth, then search for the shortest echo spacing that can be achieved at all the allowable bandwidths. Once I determine the minimum echo spacing I then avoid using any echo spacing//bandwidth combination that sets the echo spacing at less than the minimum echo spacing plus 0.02 ms. (This sounds complicated, it's actually very simple sitting on the scanner console!) Thus, in the example shown above, I found that the minimum echo spacing was 0.43 ms, attainable at a bandwidth of 2790 Hz/pixel. By using an echo spacing of 0.45 ms or longer - at whatever bandwidths are now attainable - I reduce the percentage of data points acquired on the ramps and impart temporal stability to the ghosts. I have simultaneously to make sure that I also stay away from the mechanical resonances, of course! That said, I rarely drive the system even this hard, preferring instead to stay 0.04 ms or more above the minimum echo spacing (and away from the mechanical resonance "zone" of 0.6-0.8 ms).

The next video shows how using moderate ramp sampling, with an echo spacing of 0.5 ms and readout gradient bandwidth of 2298 Hz/pixel leads to small and temporally stable ghosts. All other parameters, including the TE and matrix size, were held constant from the previous example:

So, in conclusion, you should beware the effects of overly aggressive reductions in echo spacing, e.g. as might happen if you attempt to reduce the degree of distortion in the phase encoding dimension. (The shorter the echo train length the lower the inherent distortion in the phase encoding dimension, as discussed in PFUFA Part Twelve.) Be careful that you don't introduce a temporal instability and/or increased ghosting. Check not only for the magnitude of the ghosts at the echo spacing you select, but also their temporal stability. A quick phantom experiment will show you what you need to know.

I will try to update this post in a couple of weeks' time, after I've had a chance to test on my scanner the potential for similar instabilities with coronal and sagittal slice prescriptions, i.e. those using the Z gradient (head-to-foot) for readout rather than X (left-to-right), as for the axial and axial oblique slices considered here. (See Note 4 of the "Good" coronal and sagittal data post for an explanation of why sagittal and coronal prescriptions preferentially use the Z gradient for the readout axis.)

Next post: not sure yet, but possibly an addendum to the PFUFA series, which would be Part Thirteen, on ramp sampling. Either that or the next post in this series, which will be on either distortion or dropout and tactics to reduce them.

____________________

Notes:

1. On the Routine tab, click the three dots to the right of the Phase enc. dir. field. This opens the Inplane Rotation window, shown below. Assure the Rotation angle is zero. Set it to zero if needed.

Pay particular attention to the in-plane rotation when using AutoAlign. It has a habit of introducing non-zero rotations.

2. I'm not going to get into it here but Siemens users should note that the "spatial-spectral" option in the Fat Suppr. field of the Contrast tab is implemented incorrectly. The actual slice thickness is some three times larger than the number reported in the Slice Thickness parameter field. I wouldn't recommend using the "Spatial-spectral" option unless you're a highly experienced MR person and you have some programming experience. You're gonna need to modify some source code for things to work properly... Note also that the fat suppression module resides outside of the particular EPI pulse sequence you're using, so there is no easy workaround by simply changing from, say, ep2d_bold to ep2d_pace. The problem resides in the fat suppression module itself.

3. I will do a separate post on ramp sampling soon. In brief, ramp sampling takes advantage of the fact that the k-space trajectory through the 2D plane can follow any direction and speed we like, provided we end up with a fully sampled 2D plane of (in this case) 64x64 equally spaced k-space points. Thus, our gradient shape isn't restricted. We can sample on the trapezoidal readout gradient, using the ramps up and down in addition to the flat top of the gradient. To do this we need to digitize much more quickly than the Nyquist frequency established by the field-of-view we've selected in the final image, then we simply discard the extra data points on the flat portion, and ensure that delta-k is the same area under the ramps as for the flat portion. And if this is all sorts of confusing, don't worry about it now, wait for the new post! I'll write it next, before I move on to the next post on common persistent EPI artifacts.

Common persistent EPI artifacts: Abnormally high N/2 ghosts (1/2)

2012-01-29T09:04:00.000-08:00

In this and a subsequent post I am going to cover some common situations when the N/2 ghosts can become abnormally high, i.e. higher than it is possible to achieve with comparatively small tweaks to the setup. For now I am going to restrict the discussion to temporally static, or persistent, ghosts. Furthermore, I will restrict the discussion to situations over which you can exert some control, usually through the subject setup and via EPI parameter selection. I'll cover the origins of dynamic ghosts later on in this series, once you've got a better grasp of the common persistent ghosting sources and are in a position to differentiate between a source that is intermittent and a (persistent) ghost that is being modulated by subject motion.

Before we get into the different experimental conditions that can lead to abnormally high ghosting, it is important that you are familiar with the reason why N/2 ghosts arise in EPI in the first place. So, if the following section sounds like Swahili (and you don't ordinarily speak Swahili) then I would encourage you to spend twenty minutes reviewing the section on N/2 ghosts in PFUFA Part Twelve before continuing here.

Origins of N/2 ghosts: review

Anything that causes a zigzag offset between the odd and even lines of the phase encoding dimension of k-space will lead to an N/2 ghost. Typically, one strives to minimize the zigzag at source, e.g. via shimming, and then a software tweak is applied (using navigator echoes) to clean up the residue as far as possible. The navigator echoes (Siemens uses three) are acquired immediately prior to the EPI readout (the 2D k-space). A phase correction term is derived from the navigator echoes and this term is applied to the k-space prior to 2D FT to reduce the zigzag phase differences. A successful correction is thus dependent on a good match between the phase offset obtained from the navigator echoes and the zigzag through the echo train. Discrepancies will lead to uncorrected ghosting.

There can be many reasons why the navigators don't capture accurately the zigzag across k-space. For example, the navigator echoes measure phase discontinuities arising over just 1-2 milliseconds whereas the echo train takes 20-30 milliseconds to play out. Physical effects that lead to phase evolution over the longer timescale of the echo train won't be corrected - because they won't have been measured - by the navigator echoes. However, it's not essential to understand exactly how the navigator echoes might not be able to capture the phase evolution across k-space. Instead, we can simply take an empirical approach and assess situations when the ghosting is higher than compared to some normal, or best case, level. It's then clear that the navigator echo correction was insufficient, and we should seek to eliminate or reduce the source of the phase difference.

Categorizing the origins of persistent N/2 ghosts

There are several sources of persistent ghosts that can arise from poor scanner installation or tuning, e.g. from gradient calibration errors, or from residual eddy currents. As a general rule there is little that you, as an experimenter, can do to affect these sources. These are issues for your facility physicist and your service engineer. So, from this point forward I'm going to focus attention on the factors that you have at least partial control over. Thus, for simplicity I will assume that your scanner has been installed correctly and is working nominally (as they say at NASA). The question, then, is whether you're doing everything properly. (Siemens users, see Note 1.)

As this is a blog dedicated to improving experimental techniques I am going to organize the different ghost origins along a couple of practical lines. I'll separate the effects that can be attributed to the specifics of the subject and/or the way the subject is oriented in the magnet, from those that are essentially "scanner misuse" effects. The scanner-specific effects would be seen in EPIs of a phantom whereas the subject-specific effects might require a human head to manifest themselves. This distinction is purely to assist you in making appropriate diagnostics should you (I mean when you) encounter an unexpectedly high ghost level.

Finally, while this post will focus on temporally stable ghosts, it is useful to consider the degree of temporal stability for some of the effects we'll see below. Note, however, that this can be a rather artificial distinction when considering brain data because any movement whatsoever of the subject will couple with the root cause of a major ghost and render it temporally unstable to some extent! It is thus important to recognize that in real situations the distinction between "persistent" and "intermittent" may itself be in flux. If you want me to be pedantic then I'm saying the origins of the different effects we'll consider here are temporally constant, or persistent. Later on in this series I will consider sources of ghosting that will vary temporally whether your subject is moving or not. At that point you'll have a complete overview of common ghost sources and you should be able to formulate a systematic approach to diagnosis, and remediation.

Global or local effects?

As was emphasized in Note 5 of PFUFA Part Twelve, the spatial distribution of ghosts may be global (e.g. a badly calibrated on-resonance adjustment due to scanner heating) or local (e.g. regions of high magnetic susceptibility in the frontal lobes) relative to your subject's head and to the EPI slices across it. Put another way, every single slice may exhibit essentially the same bad ghosting or you may just see some regions affected in just one or a few slices. In what follows I will do my best to suggest when an effect is more likely to be global than local, but bear in mind that there may not be a clear distinction in practice. It all depends what's screwed up and how.

Subject-dependent conditions:

Asymmetric orientation of the subject's head leading to a poor shim

GLOBAL - likely to affect all slices to some extent.

Talk about a bad way to start your session. You're in a hurry, already watching the clock, so you rush getting the subject onto the patient bed, ram a couple of pieces of foam down the sides of the subject's head, pop on the top of the head coil and it's go time!

Slow down! Take a moment to check that the bridge of the subject's nose is pointed towards top dead center of the magnet and try to ensure, as best you can by eye, that the head isn't yawed or rotated relative to his looking directly ahead. (Another clue might be that the subject reports not being able to see all of the display, e.g. a coil-mounted mirror.) If you don't you may find that you have to waste more time by repacking the subject's head once you discover that the scanner can't shim the head properly. (See Note 2.) The lower the axes and planes of symmetry through your subject's brain, the harder you have just made it for the shimming algorithm to find a good solution that will maximize brain signal and minimize N/2 ghosts.

In the left image below the subject's head was positioned symmetrically in the magnet, resulting in low ghosts for all EPIs because the shim was good. On the right, the same subject was placed in the magnet with a 5-10 degree rotation. (The subject rotated his head to look to his left very slightly.) A new shim was performed, then the EPI slices were prescribed across the brain to match those in the first case. (This has the effect of nullifying the rotation of the head in the magnet reference frame, so you no longer see the rotation in the images.) At first glance these new EPIs look acceptable, until you compare the ghost level to the left-hand image. As always, I've brought the background intensity up to highlight the ghosts:

Left: good ghost level resulting from the head placed symmetrically in the magnet, i.e. following a good shim. Right: poor ghost level arising from head rotation and a poor shim.
(Do not adjust your set! Yeah, I managed to change the volume on my Mac while I was doing a timed screen shot. Oh well.)

Distortion and dropout would also likely be higher for the rotated head, too. But the more useful diagnostic - the thing that tells you to STOP and check the subject's head orientation, is the high degree of ghosting.

Why is the ghosting higher when the head is rotated? There are usually only eight shim terms on your scanner, three with linear spatial dependence (X, Y, Z) and five with second-order spatial dependence (e.g. Z2, XZ, YZ, XY, X2-Y2) in the magnet frame of reference. The magnetic field across a brain is being affected by venous sinuses, skull, ear cavities, etc. and you are trying to minimize the variation with a low number of spatial terms. Those shims that act with an X and/or Y dependence are operating across the magnet bore, with X being left-right and Y being up-down on most scanners. (See Note 6 in PFUFA Part Seven for a typical magnet frame of reference.) The best combination of shim terms is computed from a field map. It's a tough problem to begin with, having so few shims and such a complex magnetic field. Reducing the already low intrinsic symmetry makes it harder still to compute a good solution to the field map, and reduced magnetic field homogeneity is the result. This leads to bigger k-space offsets and higher ghosting (amongst other things) in EPI.

Poor shim as a result of subject motion during/immediately after shimming

GLOBAL - likely to affect all slices to some extent.

Even when you're careful setting up your subject and you get a good initial shim, the subject may well be working against you! A sneeze or perhaps adjusting body position for comfort can easily displace the subject's head from the initial, well-shimmed position.

In the following figure the images on the left were acquired immediately following a good shim with the head placed symmetrically, and on the right is the result of a head rotation following the shim (with the subject now stationary in the new position) :

(Click to enlarge.)

The slice prescription wasn't changed so you can see clearly the rotation of the brain in addition to the higher ghost level on the right. These twin observations are your cue to have the subject return his head to the starting position (or you'd go do it for him) and re-shim. (Siemens users, see Note 3 in the post, Tactical approaches to (re)shimming for the procedure to trigger a re-shim at any point during a scan session.) A less sophisticated user might just re-acquire a localizer scan and redefine the slice prescription. That would ensure the brain coverage is as desired, but note that the shim is still sub-optimal and therefore the ghost level will remain high. You'll have to take my word for this, I'm afraid I didn't think to acquire a second set of EPIs with a straight prescription. I can guarantee that the ghost level would not change much!

Any time during an fMRI experiment, if you know or suspect that the subject has moved his head since the (last) shim was performed, do another shim! Whether or not you interrupt the session to realign the subject's head is up to you. I would reposition the subject's head only for large deviations that lead to ghosting that isn't fixed by a new shim. Re-shimming is fast - 30-40 seconds - and does no harm even if it doesn't generate a massive benefit, e.g. because the rotation was small enough not to perturb the shim significantly.

Presence of FOD or an implant causing a poor shim

LOCAL - will usually affect some slices more than others.

This one's on you again: pilot error. You can't blame the subject, not even if she forgot to take that last hair clip out. (Well, you could try, but in my lab you're in charge of screening, not the subject!)

In this example I used a spherical gel-filled phantom against which I placed a small push-pin; low-grade steel with a plastic top. The steel pin is about 15 mm long and about 1 mm in diameter. It is considerably smaller than most hair clips but might contain a similar amount of metal to some types of "scrunchie" or metallic hair band. Here's what it does to what should be rather circular slices through a spherical phantom:

(Click to enlarge.)

It should be obvious that the presence of the FOD (foreign object/debris) is massively disrupting the local magnetic field homogeneity here, leading to extreme ghosting amongst other things. Indeed, you could argue that the ghosting is the least of your worries! You've got a rather large hole (dropout) accompanied by some bizarre distortions.

The particular combination of ghosting, dropout and distortion will obviously depend on the size and composition of the conductive foreign object, as well as its location in/near the subject's head. (In an earlier post on FOD I described a couple of incidents when subjects had "missing brain.") You might also see a diffuse effect on ghosting, e.g. because of the presence of a conductive hair product or hair coloring. Or it could be multiple focal effects, as I have seen before with wigs. (You should screen for wigs, by the way!)

If the artifact is the result of an implant that somehow got through screening, e.g. a metal plate in a cheek bone, then it's nearly impossible to tell what the effect(s) will be on the resulting EPI ahead of time. Now, on occasion you might have a valuable subject who has, say, a titanium screw in his upper jaw. Assuming you can't easily replace him, and you've suitably assessed the safety risks, my usual approach is to run a quick pilot session and assess the EPI quality (including the TSNR) before wasting time with subject training, etc. I take a similar approach to subjects with a large number of metal amalgam fillings, especially in the upper jaw.

Next post: Scanner-dependent conditions leading to abnormally high (persistent) N/2 ghosts.

____________________

Notes:

1. The magnitude and temporal (in)stability of the ghosts may have an RF coil dependence. In this post all data were acquired with the standard 12-channel head coil on a Trio/TIM system. I will deal with the 32-channel head coil in a dedicated post. It turns out that the 32-channel coil isn't a TIM coil and this can lead to different ghost behavior than for the 12-channel coil.

2. Another likely consequence of hasty subject setup is subject discomfort. The subject's very first position on the patient bed is unlikely to be his most comfortable one. Give your subjects time to make small adjustments before inserting them into the magnet, otherwise you're very likely to find that the subject will make these adjustments anyway - during the scan. Best case, there might not be any significant effect if the subject moves between EPI runs and manages to keep his head in the same position. Or, worst case, you might end up with a restless subject who seeks comfort throughout the runs, trashing your experiment with excessive motion. All because you didn't want to "waste" a couple of minutes at the start of the session getting the subject positioned comfortably and correctly.

New blog! MathematiCal Neuroimaging

2012-01-27T17:54:00.001-08:00

My colleague, DS has begun a new blog, mathematiCal Neuroimaging, dedicated to exploration of the mathematical principles underlying neuroimaging methods. Here are some excerpts from the section labeled About:

As the name of the blog partially implies, the topic here will be the mathematics and physics of neuroimaging. In particular the focus will be upon functional imaging of the brain.

The style of this blog will range from tutorial-like expositions of present functional neuroimaging technology to whimsical explorations of how we might create better functional neuroimaging technology.

An entry on parallel imaging has just been posted. See you over there!

New stats/analysis blog, with MATLAB examples

2012-01-05T10:08:00.000-08:00

New blog from Kendrick Kay on stats and analysis together with some supporting MATLAB
scripts. Looks VERY useful!

http://randomanalyses.blogspot.com

Common persistent EPI artifacts: Gibbs artifact, or ringing

2011-12-30T05:48:00.000-08:00

Don't ask me why there's no apostrophe, it looks possessive to me. Perhaps it's (the) Gibbs artifact rather than Gibbs (his) artifact. Most people simply refer to the effect as ringing anyway, so let's move on. This post concerns a phenomenon that, like aliasing last time, isn't unique to EPI but is a feature of all MRIs that are obtained via Fourier transformation.

In short, ringing is a consequence of using a period of analog-to-digital conversion in order to apply a (discrete) FT to the signals and produce a digital image. Or, to put it another way, we are using a digital approximation to an analog process and thus we can never properly attain the infinite resolution that's required to fully represent every single feature of a real (analog) object. Ringing is an artifact that results from this imperfect approximation.

We had already encountered one consequence of digitization in the Nyquist criterion in PFUFA Part Six. However, for our practical purposes, ringing isn't a direct consequence of digitization like the Nyquist criterion, but instead results from the duration of the digitization (or ADC) period relative to the persistence of the signals being measured. In principle, a signal decaying exponentially decays forever, which is rather a long time to wait for the next acquisition in a time series, so we instead enable the ADC for a window of time that coincides with the bulk - say 99% - of the signal, then we turn it off. This square window imposed over the exponentially decaying signal causes some degree of truncation, and it's this truncation that leads to ringing. (See Note 1.)

An example of ringing in EPI of a phantom

Let's start with an unambiguous example of ringing by looking at the artifact in a homogeneous, regular phantom. Below is a 64x64 matrix EPI acquired from a spherical gel-filled phantom. You're looking for the wave-like patterns set up inside and outside the edges of the main signal region:

In the left image, which is contrasted to highlight ringing artifacts within the signal region itself, the primary ringing artifact appears as a series of concentric circles, each with progressively smaller diameter and lower intensity as you move in from the edge of the phantom. One section of the bright bands is indicated with a red arrow, but you should be able to trace these circles all the way around the image. Also visible is a strong interference pattern (blue arrows) that arises between the aforementioned ringing artifact and the overlapping N/2 ghosts. This is because the ghosts maintain the contrast properties of the main image; they are, after all, simply weak (misplaced) clones of the main image.

The ringing pattern in the ghosts is very easily discerned once the background contrast is brought up, in the right image above. The orange arrow indicates a clear ringing artifact in the ghost, and a little bit of mental interpolation should be sufficient to comprehend how that artifact will persist in those sections of the ghosts that overlap with the primary phantom image, leading to the interference pattern just discussed.

With the background contrasted high we can now also see that the ringing artifact extends outside of the main phantom image; it's not just circles and patterns within the signal region. The top yellow arrow on the right image highlights one clear region of ringing extending into what should be noise. There are other similar regions easily discerned on all sides of the main signal region. What's more, the ringing artifact that extends away from the sample also generates its own N/2 ghosts, as indicated by the lower yellow arrow. (The N/2 ghosts will replicate all the features of the parent image.)

So, what's happening here? As a general rule you can expect to see ringing artifacts whenever spatial resolution is low and/or there is a strong contrast boundary, as is the case for the bright circular signal region in the image above. However, in the absence of any sort of filtering (see below), it is actually the case that all images, whether 2D or 3D, will exhibit ringing. The question isn't whether it's there or not, it's whether you can see it with the naked eye or not. The reason is that although we think of the individual pixels within an image as being square (or rectangular), this isn't actually the spatial response that comes out of the Fourier transform.

What causes ringing?

As we saw in PFUFA Part Six with an FID, if we're going to use a digitizer (an analog-to-digital converter, or ADC) to acquire k-space samples then at some point it's going to have to be turned on, and at some later point it's going to have to be turned off. It's a window in time. However, in principle, once a coherent magnetization (an observable MR signal) has been generated in your subject's brain, this signal will persist (theoretically) for infinity as it decays with time constant T2*. I'm sure you see the potential conflict. As we saw in the fifth cartoon in PFUFA Part Six, turning off the ADC "too soon," while the signal is still appreciably above the noise level, results in wiggles that spectroscopists call feet; hence the term apodization for getting rid of them.

Exactly the same phenomenon occurs in imaging. Now, to comprehend the situation for imaging we first need to recall that each signal in the EPI echo train is a gradient echo, i.e. a signal peak. Either side of the signal peak the signal changes (exponentially) as some function of the dephasing imparted by the imaging gradient plus the inherent T2* relaxation. (See PFUFA Part Eight, first figure, for a review of a single gradient echo.) Still, the signal hasn't necessarily been driven down into the noise level at either side of the echo. Any signal that persists at either side of the ADC period will be "clipped."

Here's an illustration of ringing for a one-dimensional image (a readout profile) through a homogeneous cubic phantom, resulting in what should be a square 1D profile. For the purposes of the illustration I've also included cartoon noise so as to differentiate the ringing effect from random wiggles. In the left column are the k-space (or time) domain representations while in the right column are the image (or frequency) domain responses produced from 1D FT of the left partners. The top row shows the ideal case, the bottom row shows the response from the ADC period by itself, and the middle row shows the experimental result: a time domain signal that has been "clipped" at either end by the ADC period, resulting in a 1D profile that has sinc-like wiggles at edges:

Please note that the slight discontinuity in the center of the top-right profile was a slip of my pen! It's supposed to be a flattish profile with noise. Consider the discontinuity as a large noise spike!

You will remember from PFUFA Part Six that multiplication in one domain - here the k-space domain - is equivalent to convolution in the conjugate domain (image space). The sinc response from the ADC period convolves with the ideal profile to generate pronounced waves at the edges of the signal, both extending into the signal itself and into the adjacent "noise" (or no signal) regions around the phantom.

Now, there's something subtle going on here that can help you understand why you see ringing in certain regions of an image and not others. In PFUFA Part Eleven you saw how different spatial frequencies within an object, such as a brain (or a Hawker Hurricane aeroplane), "live" in different parts of reciprocal, or k-space. High spatial frequencies - edges, fine details - live in the periphery of k-space. Now think about the k-space trajectory for the 1D gradient echo that you saw in PFUFA Part Nine. The center of k-space - where all the low spatial frequencies live - is by definition in the center of the k-space trajectory, and in the center of the ADC period (if we're being smart and digitizing the entire k-space line!). What lives at the either end of the k-space trajectory, the parts that are going to get chopped off (clipped) first, should there be appreciable signal as the ADC switches on or off? That's right, the high spatial frequencies! Thus, the act of clipping the signal affects the high spatial frequencies - edges - more than it does the low spatial frequencies. And so the strongest effects of ringing tend to be detectable where there are fine details or boundaries between high and low intensity signals.

An analogous situation arises for the phase encoding dimension. We define an effective digitization window and, as in the frequency encoding axis, increasingly large values of k produce signal attenuation as well as spatial encoding. If infinite k values in the phase encoding dimension were obtainable then, as for the readout direction, there would be no "clipping" and hence no ringing. But that isn't feasible in practice and so we have ringing in the phase encoding dimension in just the same way as in the readout dimension. That's why the real example of ringing at the top of this post was so symmetric.

Can ringing be reduced or fixed?

In PFUFA Part Eleven we saw how extending k-space increased image resolution, and that restricting the extent of k-space resulted in low resolution images. It's important to remember that the k-space matrix is, as its name implies, a discrete space but we are using it to map something that is continuous in space. Real brains aren't made up of pixels!

Increasing image resolution makes the absolute distance affected by ringing lower, because the pixels are getting smaller. But this doesn't actually fix the problem; it tends to restrict the spatial extent (which may be all you need to do). Let's take a look at the effects of resolution on the spatial extent of ringing by comparing three different resolutions of EPI in a spherical phantom:

Left: 64x64 matrix. Middle: 96x96 matrix. Right: 128x128 matrix. (Click to enlarge.)

Parameters other than the resolution were held constant as far as possible, but to get to the higher resolution it was necessary to increase the echo spacing (thereby increasing the distortion) and also increase the TE (thereby decreasing the signal-to-noise ratio) slightly.

Still, increasing the resolution to 128x128 hasn't eliminated the ringing, just reduced the spatial extent and made it harder to detect by eye. If you wanted to reduce the ringing further then even higher resolution would help, or you could apply a smoothing function to the data. I don't want to get side-tracked today with smoothing because this post is on recognizing the artifact, not fixing it. I may do an exhaustive post on smoothing another day but for now see Note 2 for some thoughts on smoothing.

What does ringing look like in real data?

Compared to the slice through a spherical phantom above, detecting ringing in even low-resolution (64x64 matrix) EPI of brains can be quite a challenge, perhaps leading you to wonder why I'm bothering to do an entire post on the subject. Well, it's because as I've already said: ringing is ubiquitous and I want you to be able to differentiate the relatively benign ringing artifact from other, more important artifacts for fMRI, those where timely intervention might make the difference between good and not so good data quality.

To ease the transition from EPI of phantom to EPI of brain, I'm going to take a short detour via regular "spin warp" imaging of brain. In these anatomical scans the SNR is generally a lot higher than for EPI with the same matrix (because TE is much lower and each phase encoding line is acquired after a new RF excitation, generating a signal averaging effect for the final image) and this tends to make the ringing artifact more easily detected. (There is also less inherent smoothing arising from T2* than for EPI, but I'll deal with this fact below in more detail.)

Here is an anatomical (gradient echo spin warp) 128x128 matrix 2D image through a brain with contrast to highlight anatomy (left) and background noise (right):

(Click to enlarge.)

The ringing both inside and outside of the brain is pretty easily recognized without the yellow arrows. The most prominent ringing artifacts arise from the black-white contrast border at the sides of the brain and from the scalp fat around the head.

Now you're tuned into seeing the wave-like patterns, see if you can make out similar patterns within the brain in these 64x64 matrix axial slices of EPI:

(Click to enlarge.)

It's considerably harder to see in the EPIs than in the (higher resolution) anatomicals above! Why should this be, given that earlier on in this post I showed that the low resolution should make the effects of ringing more visible, not less? There are two major reasons. Firstly, the SNR in the EPIs is five to ten times lower than in the anatomical scan. Thus, all other factors being equal, our ability to detect such a small feature is greatly reduced. We often don't have the SNR in EPI to tease out fine details, not even those that are artifacts! Secondly, in the phase encoding dimension of the EPI scans there is an effective smoothing function being "applied" by virtue of the T2* decay that is occurring simultaneously with the readout echo train. Recall from PFUFA Part Twelve that the practical limit for the number of echoes in our train is set by T2*; we're in a race to recycle and reuse the magnetization as many times as we can until T2* renders the signal level too low. The phase encoding dimension might well experience reduced ringing as a result. But not always.

Inferior regions of the brain suffer extensively from magnetic susceptibility gradients, leading to signal dropout as well as higher ghosting and distortion of any remaining signal. For this remaining signal, then, the presence of the susceptibility gradients will tend to interact with the desired, or ideal, k-space matrix and may have concomitant effects on ringing. We know that the signal maximum should occur at the center of k-space as defined by the readout and phase encoding imaging gradients. But, if a region experiences high magnetic susceptibility gradients these additional dephasing processes can cause the signal for that local region to attain its maximum early (if the susceptibility gradient adds to the phase encoding gradient) or late (if the susceptibility gradient tends to offset some of the phase encoding gradient). Furthermore, the displacement causes one side of the 2D signals to shift towards the edge of the assumed (ideal) k-space matrix until they start to "fall off" the plane, like this:

Illustration of the effects of background magnetic susceptibility gradients on k-space. In this example the susceptibility gradients have caused the early refocusing of signals along the phase encoding dimension, ky.

It should be fairly obvious that signals located towards the bottom of this k-space matrix have been "pushed off" the sampled (ideal) k-space plane by magnetic susceptibility gradients and would therefore be truncated (clipped) more severely than they should have been, had the signals been properly centered at the k-space origin. This will cause parts of the resulting image - those having appreciable k-space representation in the clipped signal - to exhibit increased ringing.

Now, in my cartoon example above I shifted the signals for the entire slice earlier (down) in ky. That may happen in an experimental situation, but it's not common. What is more likely to occur is the shift of signals from just a localized region, with the remainder of the signals staying closer to the intended k-space origin. This results in some parts of the slice exhibiting pronounced ringing but not others, as in this example:

An inferior axial oblique EPI slice showing localized ringing.

Understanding why this ringing should be localized should be relatively straightforward. After all, dropout is usually localized, too. The magnetization residing in these poorly shimmed regions has managed to survive the aggressive dephasing of the susceptibility gradients to the point where it will yield appreciable (non-zero) signal in the final image, but the environment is so contaminated by poor magnetic field homogeneity that the signals may be very disfigured indeed, often being characterized by excessive distortion, high N/2 ghosting and localized ringing.

Is ringing a problem for fMRI?

As with many artifacts in fMRI, the answer to this question depends on where you are in the brain, but it also depends on what you subsequently do with your EPI time series in the fMRI analysis.

We have already noted that the ringing effects are quite subtle, both because of the inherent smoothing and because of the relatively low SNR (compared to anatomical scans). But we have been considering (raw) EPI images, not the spatial response that you use in your fMRI analysis necessarily. For example, there are several reasons to smooth EPI data that have nothing to do with ringing. It is common to want to coalesce "active" pixels in order to get reasonably sized activation "blobs" that pass a statistical threshold. (We usually assume that truly active pixels don't appear in isolation in image space, a reasonable assumption given the broad spatial response inherent to BOLD imaging.) Another common reason to smooth is to reduce the total number of pixels in the data set, thereby allowing a valid reduction in the aggressiveness of a Bonferroni correction for multiple comparisons. Whatever the motivation, even moderate smoothing will likely eliminate the mild ringing that might be barely visible in the raw EPIs.

For the most part, then, we tend not to worry specifically about Gibbs artifact in fMRI. If you're smoothing your data then you almost certainly have no need to be concerned with ringing further. But if you find that your EPIs are exhibiting substantially worse ringing than the examples presented here then it's time to go talk to your friendly facility physicist. There could be a problem with your acquisition parameters or even a problem with the imaging gradients (the hardware) themselves. (See Note 3.)

Ringing and motion

In this post I've dealt with the static effects of Gibbs ringing. In later posts we will be looking at the way motion tends to interact with artifacts such as ringing. Until those in-depth posts I will leave you with a brief observation. Contrast of all sorts, whether it's "real" anatomical contrast arising from T1/T2*/spin density differences in the brain tissues, or artifact banding from Gibbs ringing, will tend to exacerbate the effects of motion to some extent. So, if you wanted to be pedantic (and why not?), you could argue that the banding exhibited along the inside edges of the brain in the axial EPIs above would lead to enhanced degradation of temporal SNR compared to a situation where those artifacts weren't present. True enough. Except, of course, there are lots of other valid edges - most obviously the margins of the brain - that must be present and that will exhibit similar degradation of TSNR in the presence of motion. All you can do is apply motion correction, perhaps, and hope that the TSNR is returned.

Thus, it seems to me that to worry about ringing as a source of motion sensitivity is to overdo it. Motion is a concern for other, more critical reasons than the way it interacts with ringing. At least, that is my opinion today. If I find any information to the contrary I shall be sure to post it immediately!

_______________________

Notes:

1. There is another way to comprehend the ringing phenomenon: as a series of sinc functions building up ever sharper features of the real object. This is a perfectly valid, and technically correct, way to approach the phenomenon. However, I decided against using it as the description here because I've always thought that consideration of the square ADC window and the relative duration of the MR signals is a more useful approach for experimental technique. Once you have a sense of the relative durations of the sampling and signal periods then you can use that intuition to predict qualitatively the effects of increasing the magnitude of a spatial encoding gradient, say (it causes more rapid MR signal decay), or the effects of magnetic field heterogeneity (which also causes accelerated signal decay). There is a good explanation of the so-called "point spread function" and the true spatial resolution of an image in Rick Buxton's book, Introduction to Functional Magnetic Resonance Imaging, 2nd edition pages 220-231.

2. To fully understand the effects of smoothing requires consideration of the "point spread function," the mathematical relationship that defines the actual shape of pixels in an MRI. We'd like pixels to have ideal, rectangular shapes but they aren't. In the absence of smoothing the pixels are actually sinc-shaped, i.e. the majority of signal is located at the pixel's nominal position, but some signal from that (real) position extends away from the pixel into neighboring pixels. There is a little bit of spatial overlap.

As mentioned in Note 1, there is a good explanation of the PSF and the true nature of spatial resolution of an image in Rick Buxton's book, Introduction to Functional Magnetic Resonance Imaging, 2nd edition pages 220-231. Other textbooks deal with the PSF in some fashion, too. I encourage you to find one and read the relevant sections because I'm afraid I can't find a nice overview (for MRI) online. Perhaps I'll try to write one and append it to the PFUFA series as Part Thirteen.

Until a more comprehensive treatment, then, let's quickly go through some of the ways that ringing can be manipulated beyond simply changing the acquired matrix (which isn't usually an option in EPI for fMRI). Another alternative is to apply the imaging equivalent of "apodization" that was mentioned in PFUFA Part Six. In that post the ringing was shown for a one-dimensional example. But the same principles can be applied to 2D images. When it's applied to images we tend to refer to it as smoothing because in addition to reducing the ringing the effect is to blur, or smooth, the image, i.e. to reduce the absolute resolution in the image plane. That's the price of the fix.

Now, it may well be that you want to apply smoothing to your fMRI time series data anyway, e.g. to increase functional SNR, or to reduce the multiple comparisons problem. Thus, if you know ahead of time what level of smoothing you're going to use then you could, if you insisted, include a smoothing filter at the data acquisition stage. My personal opinion is that this isn't a good idea; just apply the smoothing offline. (Multiplication of a smoothing function in the k-space domain is equivalent to convolution of the FT of the smoothing function in the image domain. You can do either.) Then you have the luxury of testing different smoothing kernels to find an optimum, or deciding to run some data-driven analyses that could benefit from the highest possible spatial resolution, etc.

"Zero filling" is another commonly used trick to reduce the apparent effects of ringing. It's a form of interpolation in the image domain. Here, instead of using a 64x64 matrix of k-space values, the k-space matrix is extended, or padded, up to some higher matrix, say 128x128, by appending matrix values having zero signal level (and zero noise) to the acquired data. This extension of zeros adds no new information to the matrix, and it doesn't change the nominal resolution of the resulting image (i.e. the maximum acquired k-space values). But it does produce (in this example) twice as many pixels with which to define the image plane, thereby lessening the ringing effects detectable to the naked eye. Here's an illustration of the effects of zero filling on an image:

Both images have the same nominal resolution, but the right image has the illusion of better resolution because it doesn't appear pixellated.

A final note about zero filling: it's applied in the k-space domain, prior to Fourier transformation. Thus, it's unlikely that you will have much facility to use zero filling on the scanner, unless your center has the capacity to pipe the raw (k-space) data off the scanner and you can do all the signal processing steps, including 2D FT, offline. That's not the case for the vast majority of routine fMRI labs, and I wouldn't think that most sites would want to install such a pipeline for the purposes of zero filling. Something a lot more useful would need to provide the motivation! But given that any sort of interpolation in the spatial domain isn't common (as far as I know) in fMRI, it's more likely that you will use a simple smoothing function, such as a Gaussian smoothing kernel, offline and be done with the issue.

Incidentally, the difference between smoothing and interpolation? With interpolation you're adding data points, with smoothing you're throwing them away.

3. This discussion concerns fully sampled EPI, not partial Fourier or parallel accelerated (e.g. GRAPPA) variants of EPI. I will deal with these options separately, noting in passing that increased Gibbs ringing is one possible negative consequence of these data reduction methods. Thus, if your EPI does exhibit severe Gibbs ringing the very first thing to check is whether you're using a partial k-space technique and it's actually working as it should.

Another brief explanation of decoding

2011-12-27T10:49:00.000-08:00

Here's another short video produced by UC's media people in which Jack Gallant explains in broad terms how his group's recent decoding experiment was conducted:

A good place to go next for more details is the Gallant Lab website. Read the FAQ on that page to gain a basic understanding of what the experiment was, and what it wasn't. Then go read the paper, it's written very accessibly!

Common persistent EPI artifacts: Aliasing, or wraparound

2011-12-14T22:18:00.000-08:00

In Part Eleven of the series Physics for understanding fMRI artifacts (hereafter referred to as PFUFA) you saw how setting parameters in k-space determined the image field-of-view (FOV) and resolution. In that introduction I kept everything simple, and the Fourier transform from the k-space domain to the image domain worked perfectly. For instance, in one of the examples the k-space step size was doubled in one dimension, thereby neatly chopping the corresponding image domain in half with no apparent problems. At the time, perhaps you wondered where the cropped portions of sky and grass had gone from around the remaining, untouched Hawker Hurricane aeroplane. Or perhaps you didn't.

In any event, you can assume from the fact that this is a post dedicated to something called 'aliasing' that in real world MRI things aren't quite as neat and tidy. Changing the k-space step size - thereby changing the FOV - has consequences depending on the extent of the object being imaged relative to the extent of the image FOV. It's possible to set the FOV too small for the object. Alternatively, it's possible to have the FOV set to an appropriate span but position it incorrectly. (The position of the FOV relative to signal-generating regions of the sample is a settable parameter on the scanner.) Overall, what matters is where signals reside relative to the edges of the FOV.

Now, on a modern MRI scanner with fancy electronics, aliasing is a problem in one dimension only: the phase encoding dimension. (Yeah, the one with all the distortion and the N/2 ghosts. Sucks to be that dimension!) The frequency encoding dimension manages to escape the aliasing phenomenon by virtue of inline analog and digital filtering, processes that don't have a direct counterpart in the phase encoding dimension. Instead, signal that falls outside the readout dimension FOV, either because the FOV is too small or because the FOV is displaced relative to the object, is eliminated. It's therefore important to know what happens where and when as far as both image dimensions are concerned. One dimension gets chopped, the other gets aliased.

I will first cover the signal filtering in the frequency encoding dimension and then deal with aliasing in the phase encoding dimension. Finally, I'll give one example of what can happen when the FOV is set inappropriately for both dimensions simultaneously. At the end of the process you should be able to differentiate the effects with ease. (See Note 1.)

Effects in the frequency encoding dimension

Below are two sets of EPIs of the same object - a spherical phantom - that differ only in the position of the readout FOV relative to the phantom. In the top image the readout FOV is centered on the phantom, whereas in the bottom image the FOV is displaced to the left, causing the left portions of the phantom signal in each slice to be neatly, almost surgically, removed:

Readout FOV centered relative to the phantom.

Readout FOV displaced to the left of the phantom, resulting in attenuation of the signal from the left edge of each slice.

I set the contrast so that you could see the background noise and verify that there's no sneaky trace of the removed signal somewhere else in the image. It's really gone. But where? And why?

The readout dimension is the one for which the k-space sampling is actually performed, using an analog-to-digital converter (ADC). As you will recall from the various posts on k-space, it's only necessary to sample the rows of k-space to fully sample the 2D plane; the 2D k-space matrix is simply a stack of these rows. What this means practically is that, as the name suggests, the readout dimension (a.k.a. frequency encoding dimension) is the one for which data points are obtained from the ADC, the rows of analog signals having been passed to the ADC from the receiver RF coil electronics.

Now, whenever you have a stream of data points you can filter them, e.g. with a passband that rejects frequencies above and below some bandwidth. This is what happens to each row of k-space data points as they are processed by the ADC. By itself this would tend to attenuate signal (and noise) from frequencies (which we know correspond to spatial positions) that are outside the passband. Except that analog filters are imperfect. So, to clean things up your fancy scanner is also equipped with digital filters that trim the passband until it performs as a near perfect square set at the FOV. This step is achieved by a high degree of oversampling (multiples of the specified Nyquist frequency), and it's all done invisibly to you, the operator. What we end up with is a neatly trimmed image. No signal (or noise) from outside of the defined readout FOV survives the filtering process. Pretty nifty, huh?

Effects in the phase encoding dimension

The columns of phase-encoded data points aren't detected in a continuous stream like the frequency-encoded data points, but are instead built up from a succession of readout data points, each row possessing a slightly different phase. (See PFUFA Part Twelve for a review.) The analog-to-digital conversion occurs for just one phase encoding value at a time, for each entire readout period, then the digital results are stacked up to produce the final 2D k-space matrix ready for Fourier transformation. Thus, the properties of the ADC, including analog filtering and in-line digital signal processing, only get an opportunity to operate on one dimension - the frequency encoding dimension - of k-space. The digital properties of the phase encoding dimension can be thought of as a "synthetic" construction of a stack of digital readouts. And this leads to different properties for the two dimensions when it comes to dealing with extraneous signal; signal that lies outside the defined FOV. (See Note 2.)

Below is a typical example of aliasing in the phase encoding dimension. Displacement of the sample relative to one edge of the phase encoding FOV (and vice versa) leads to aliasing, or wraparound, of part of the signal to the opposite side of the image:

Left: phase encoding FOV centered relative to the sample. Right: the phase encoding FOV has been displaced such that part of the sample falls outside, causing that part of the signal to alias to the other side of the image. Note that the N/2 ghosts are aliased as well.

It's almost as if the entire image is on a conveyor-belt. What falls off one side seems to be cycled around and deposited on the other side of the image. (See Note 3.)

Why does aliasing occur in this dimension when, in the readout dimension, the signal outside of the FOV would have been removed by filtering? The algebra for the relationship between the k-space step size (in the phase encoding dimension) and the FOV was given in Note 4 of PFUFA Part Eleven. In brief, for each k-space increment there is a total of 360 degrees of phase, +/- 180 degrees about the center on-resonance (or "carrier") frequency, imparted across the sample in the direction of the applied gradient. (This is true for readout and phase encoding gradients, but let's focus on the phase encoding gradient from now on.) Put another way, at the extreme spatial positions of the FOV - the very edges of the image in the phase encoding dimension - the phase evolution is exactly +180 degrees on one side and exactly -180 degrees on the other. At all spatial positions in between the phase is somewhere in between this range.

But what if there is some sample residing outside of these extreme positions? Magnetization there still "feels" the effect of the applied gradient, and it still gets a phase imparted to it. The problem is that if, say, the signal is left of the left-hand edge of the FOV then the imparted phase will be greater than +180 degrees (if we assume left edge is +180, right edge is -180). But, because phase is modulo(360) it means that a phase value of +200 degrees (say) is indistinguishable from a phase of -160 degrees. And you should immediately recognize where that signal past the extreme left position will end up: at the -160 degree position, which is actually on the right-hand side of the image FOV!

Note that it's not just the signal regions that get cycled (aliased) around the FOV; the N/2 ghosts do, too. Indeed, the relative position of the ghosts to the sample hasn't changed because the zigzag phase modulation across the phase encoding k-space (see PFUFA Part Twelve) is simply phase-shifted by a constant amount; the zigzag itself is unchanged.

Effects in both dimensions when the image FOV is too small

Finally for this post let's look at what happens when the FOV is made smaller than the signal extent in both the readout and phase encoding dimensions. Here I've simply made the FOV smaller than the phantom:

By now I hope you can immediately determine which dimension is frequency-encoded and which is phase-encoded. You can look for the aliasing in the phase encoding dimension or, more generally, just look for the N/2 ghosts. (Distortion can be harder to see, and it's impossible to see in this image!) Okay, so phase encoding is the vertical dimension in this image, and it's easy to identify the two overlapped signal regions top and bottom. Note how the aliased signal simply adds to the correctly located signal, increasing the SNR in the overlapped regions. You've already identified the N/2 ghosts (right?), but in case you haven't, the easiest parts to see are the two crescents (a smiley and a frowny) right through the middle of the image. The FOV is so small that the ghosts themselves overlap! (You might just be able to make out an interference pattern, established by the aliased ghosts, in a horizontal line across the center of the image.) The readout dimension is neatly cropped at the edges and no aliasing occurs.

Pretty easy stuff. Don't worry if you didn't understand why the frequency encoding dimension gets cropped while the phase encoding dimension aliases. It's fine to commit these features to memory and be ready to interpret the differences whenever signals get close to the edges of your image FOV. Be especially aware if you are using small FOVs, such as 192 mm FOV for axial EPI. It's critical to know if scalp fat, say, will alias or get cropped based on where it falls relative to the image dimensions. You don't want to alias fat signals onto brain!

Next up: Gibbs ringing. No, nothing to do with hand bells or Christmas caroling.

_____________________

Notes:

1. What you're reading here applies equally well to other forms of MRI; these aren't just phenomena affecting EPI. However, as this is a blog dedicated to experimental fMRI I am going to use EPIs as the example images and discuss only related issues (such as the aliasing of ghosts) that pertain to EPI.

2. Siemens has a parameter called Phase Oversampling (on the Routine tab) which might look like a way to circumvent aliasing, but it's not really. It simply acquires additional k-space steps in the phase encoding dimension - by the percentage as set in the parameter - and in concert increases the FOV invisibly to the user. Then, once the acquisition is complete and the data has been transformed into image space the image is simply cropped to leave the FOV and resolution you specified in the Resolution tab. In other words, the extraneous data from the enlarged FOV is simply discarded.

Now, given that in EPI the echo train must be extended to increase the number of phase encode steps, this "trick" has a cost. It increases the echo train length, thereby increasing the minimum TE and altering (generally increasing) the blurring due to T2* during the echo train. It also decreases the number of slices permissible in the TR period. That said, it is sometimes possible to set the Phase Oversampling parameter greater than zero without having to change any other parameters, but your images will still have the increased T2* blurring.

This oversampling feature is equivalent to explicitly increasing the desired FOV in the phase encoding dimension, while (explicitly) increasing the number of points in that dimension to maintain a constant nominal pixel resolution and then, once the data is acquired, simply zooming (or cropping) the image back to the originally desired FOV. Thought about this way it becomes obvious that there's no free lunch here, and all that the Phase Oversampling feature really achieves is an indirect way to increase the FOV and phase encode points in concert, then save a tiny amount of space in the database by not storing a larger image. Overall, I don't see much utility in Phase Oversampling for circumventing aliasing in EPI.

3. It is possible to reconstruct an unaliased image from an aliased image, but it's not good practice. You may find that heating effects or motion in the phase encoding direction leads to "interesting" statistics for the lines of pixels that fall along the edges of the aliased image. But, in a pinch, I'd probably try it to salvage an important data set. Of course, the point of reading this blog is that you become sufficiently skilled that you don't accidentally alias your images!

Twitter. Damn.

2011-11-28T10:08:00.000-08:00

@practiCalfMRI

I thought I could resist, I really did. (I've been off Facebook* for more than two years!) But when Neuroskeptic took the plunge in June I started thinking that maybe I should suck it up, too. I mean, Neuroskeptic blogs ten times more frequently than I do and he still has time to tweet.

Not sure exactly how it will go. I'm going to treat it as an experiment. I can guarantee that there won't be daily tweets let alone hourly ones. I'm not going to bring anyone's cellular network to its knees. But I do come across little things related to fMRI that aren't worth a full blog post. A micro-blog ought to fit the bill, eh? We'll see...

* Okay, so technically I am still on fb. I maintain an account so that I can post comments to websites, merge with other online media, etc. Turns out it's really, really hard not to have fb unless you don't mind registering separately for every online newspaper and music service yet invented. But I never actually look at my fb page, so I apologize if you have a "friend" request suspended somewhere in cyberspace.

Understanding fMRI artifacts: "Good" coronal and sagittal data

2011-11-27T11:09:00.000-08:00

Front, back, side to side

Now that you have an appreciation of "good" axial EPI time series data we should be able to zip through a review of "good" coronal and sagittal EPIs. This isn't the post to get deep into the reasons why you might want to acquire these prescriptions instead of axial or axial-oblique slices, but here's a short list (and some music) for you to be going on with:

Pros

coronal slices tend to exhibit less dropout of frontal and temporal lobes compared to axial slices.
coronal slices might permit a smaller field-of-view and higher spatial resolution without signal aliasing than achievable with other prescriptions, assuming your gradient performance and other pulse sequence parameters can be driven sufficiently hard.
sagittal slices may also show some improved signal in frontal and temporal lobes compared to axial slices, but the real benefit is the unique coverage afforded. You could acquire a single hemisphere, for instance; could be useful in a handful of situations. Alternatively, if you are interested in the whole brain, including cerebellum and perhaps even brain stem, these structures are naturally included in sagittal slices.
sagittal slices tend to make the most common type of head motion - chin to chest rotations - an in-plane phenomenon which might lead to improved motion correction in post-processing.

There are, naturally, drawbacks to coronal and sagittal slices, just as there are for axial slices. I'll mention some of these in more detail below, as we consider the individual artifacts, but here's another brief list:

Cons

safety limits on gradient switching (to avoid peripheral nerve stimulation) tend to force the phase encoding direction to be left-right for coronal slices, rendering the EPIs strongly asymmetric. While the absolute level of distortion may actually be very similar to that present in axial slices, the disruption of left-right symmetry can be a shock to your aesthetic sensibility.
bizarre distortion is also a "feature" of sagittal slices where, as you'll soon see, the distortion can make the frontal lobes look like a duck's bill! But, as before, the absolute level of distortion may not be significantly different to that in axial slices; it's really the unnatural appearance that shocks us. (We ought to be just as outraged at the symmetric distortions in axial slices!)
perhaps the biggest limitation to both coronal and sagittal prescriptions is the number of slices required to cover the entire brain in the given TR. Slicing along the longest axis of the brain, as done for coronal slices, is clearly the least efficient way to do it. The efficiency of sagittal slices falls somewhere between coronal and axial. And, of course, anything that leads to more (fixed width) slices means that TR might have to get longer. It all depends on your application.

Okay then, that's the introduction over with. Let's now put aside the justification for using one prescription over another and look at what constitutes "good" data in the case of coronal and sagittal slices. The features should be immediately recognizable from what you saw in the axial data of the last post.

Parameters

The data we will consider in this post used parameters comparable to those for the axial images in the last post. I used the same single shot, gradient echo EPI sequence on a Siemens Trio/TIM scanner, using the 12-channel head RF coil and a pulse sequence functionally equivalent to the product sequence, ep2d_bold. (See Note 1.) Parameters were: 33 slices, 3 mm slice thickness, 10% slice gap, TR=2000 ms, TE=28 ms, flip angle = 78 deg, 64x64 matrix over a 22.4 cm field-of-view yielding 3.5 mm resolution in-plane, full k-space with phase encoding oriented either left-right (coronal slices) or anterior-posterior (sagittal slices).

To accommodate the scanner's safety software the gradient switching rate had to be reduced slightly compared to the axial prescription, making the echo spacing slightly longer (0.5 ms instead of 0.47 ms) and permitting a maximum of 33 (rather than 34 axial) slices in TR=2000 ms for the coronal and sagittal acquisitions. I acquired just 20 volumes in each time series; enough to demonstrate the temporal features. See Note 2 if you'd like to download the raw data (or the movies and jpeg images) you see below.

Coronal EPI

Start out by reviewing the cine loop through the twenty frames a few times. Click the 'YouTube' icon on the embedded video to launch an expanded version in a separate tab/window if you prefer (or to download all the videos and jpegs from this post, see Note 2):

The most obvious feature in these coronal images is the left-right distortion which gives the brain a rather perverse, windswept appearance. As I've already mentioned, most of what you're seeing is, in quantitative terms, not appreciably worse than what you deal with in axial slices. Axial slices (with phase encoding set anterior-posterior) do at least maintain a semblance of left-right symmetry, even though we know that many brain regions (particularly frontal lobe) are considerably displaced from their true locations and, indeed, the entire brain appears more oval (stretched) than it should. What disturbs our aesthetic sense in these coronal slices is the disruption by a left-right shearing of what ought to be moderately symmetric anatomy. The problem is especially noticeable for the temporal lobes. (See Note 3.)

If you're wondering why I didn't simply make the phase encoding direction head-foot, the answer is scanner performance and subject comfort. Swapping the read and phase encoding axes means that the rapidly switching read gradient would become left-to-right, from head-to-foot. This might cause more peripheral nerve stimulation, so the scanner's safety software imposes stricter limits on the gradient timing and magnitude parameters in that configuration, and that would translate into longer echo spacing, higher distortion (albeit symmetrically in the brain!) and fewer slices per TR. (See Note 4 for more explanation on stimulus limits.)

Note that throughout the twenty frames of the movie all signal regions, whether strongly distorted (such as temporal lobes) or apparently undistorted (such as the left and right sides of the occipital lobe), remain stable; it is difficult to discern much movement of brain anatomy. The only obvious movement is happening in the neck, because of the pulsation of the carotid arteries (in which blood can be seen flowing as bright vertical lines in a slice in the second row, right-hand side).

What about the ghosts? The N/2 ghosts can be seen in the bottom few rows, but to really see what's going on we want to crank up the background intensity. Here is the same loop of twenty frames but slowed down to 5 frames/sec:

You'll have to take it from me that the ghosts have an intensity about 5% of the brain signal. (You don't really want ghosting more intense than this.) Importantly, the ghosts from the brain signal remain stable over the series.

Let's try to confirm the stability we think we're seeing in the movies by looking at the temporal SNR (TSNR) and standard deviation (SDEV) images:

TSNR image. (Click to enlarge.)

Standard deviation image. (Click to enlarge.)

The SDEV image reveals the carotid arteries like an angiogram! The superior sagittal sinus is also visible, particularly in the bottom row, as a bright spot atop the midline. In the brain we should expect to see CSF with high SDEV because of cardiac-driven pulsations, just as we saw for the axial slices. And, as before, because gray matter is more metabolically active (and more vascularized) than white matter, GM has higher SDEV than WM. All perfectly normal.

There's an interesting edge that highlights the brain in every slice of the SDEV image. The edge is right along the brain's surface and connects to sulci frequently, strongly suggesting that it arises from pulsating CSF in the meninges, not (whole) head motion or gradient heating. I spent a lot of time inspecting the movies at a host of contrast and zoom combinations, and I couldn't detect by eye a substantial difference for the edges versus the bulk of brain signal. The edges look physiological to me.

Assessing the TSNR image doesn't add much new insight to the picture, except to confirm that there are no brain regions that suffer from a low TSNR when they exhibited high signal in the raw data series. This might happen if, say, there were intense ghosts overlapping portions of the brain that would serve to reduce the TSNR of the overlapped region. (Of course, we would also see these regions having a correspondingly high SDEV, and we don't see any high SDEV ghosts.)

Sagittal EPI

The distortion in the sagittal EPIs makes the brain almost as peculiar as in the coronals. Magnetic susceptibility gradients mutate the frontal lobe into what looks like a peak on a cap. As before, the phase encoding direction is preferentially set anterior-posterior (rather than head-to-foot) to minimize the potential for peripheral nerve stimulation arising from the read gradient train. (See Note 4 again.) Anyway, loop through the movie a few times and see what changes:

Other than the direction of distortion, the features in the movie are quite similar to the previous coronal data and our previous "good" axial data. You now immediately recognize the carotid arteries pulsing away in the neck. If you're wondering why the blood in the carotids appears bright it's because blood is flowing up from the heart where it hasn't experienced the slice-selective RF pulses of prior TR periods. Thus, the "apparent T1" for blood in the carotids is much shorter than the actual T1 for arterial blood. Once the blood has been present in the head for a few seconds, however, it achieves a more conventional T1 steady state, and its signal level is reduced. (The T1 for blood is about 2 seconds at 3 T, so there's only partial relaxation at TR=2 sec.)

Interestingly, in this particular example of sagittal data the subject (me!) seems to be maintaining very good fixation because the eye signal is stable. It is only 40 seconds of data, however, so in a typical fMRI run of a few minutes you should expect to see more eye movements, as we saw previously in the axial data.

There's not much else to talk about, which is a good thing. All signal boundaries appear stationary by inspection. It is difficult to discern with certainty regions of N/2 ghosts so let's follow standard procedure and replay the movie with the background contrasted high:

Now we can clearly identify ghosts - they're most easily identified in the first and last few slices, from the left and right sides of the head, because the volume of brain signal is lower and the background area is larger. These ghosts do not appear by inspection to fluctuate greatly. Most of what's changing is signal in the neck, which one expects from the pulsatile effects of the carotids, perhaps swallowing, small jaw movements and so on. The brain signal appears to be quite stationary.

Okay, time to confirm our cine loop diagnoses with the TSNR and SDEV images:

TSNR image. (Click to enlarge.)

Standard deviation image. (Click to enlarge.)

Unsurprisingly, the carotid arteries are a main feature in the SDEV image. Another recognizable feature is the edge enhancement of the brain, much as we saw above in the coronals. Sulci are discernible along the edges again, too, so it's reasonable to assume that this feature arises from pulsatile CSF flow (and for the central few slices perhaps, from the superior sagittal sinus as well). The SDEV image confirms that fluctuations in the N/2 ghosts are low; their standard deviation is barely above that of the background.

The TSNR image shows that in spite of the high degree of distortion the signal in the frontal lobes is viable for fMRI when (as here) the head motion is low. Ugly doesn't imply useless.

Some final thoughts about good data

Now that you have an appreciation of what good time series EPIs look like, and noting that we have only discussed the principal slice prescriptions and a conservative set of parameters in each case, you are ready to start looking at the effects on individual images and time series data when things go wrong; what your physicist will tell you is ‘bad’ or sub-optimal data.

Be warned, though: many of the artifacts in EPI time series have remarkably similar appearance in spite of radically different origins. It can take some time to be able to discriminate between ‘good’ and ‘bad’ data, let alone to become accomplished at discriminating between the different problems when they do arise. That, however, is your task because the better you can diagnose an artifact during your session the faster you can correct it! And with that in mind it's time to shift the focus onto bad data, starting with the next post.

__________________________

http://dl.dropbox.com/u/26987499/Good_EPI_Coronal_64x64.zip
http://dl.dropbox.com/u/26987499/Good_EPI_Sagittal_64x64.zip

If you don’t already have a DICOM viewer, check out Osirix for Mac OSX (available via a link in the sidebar). ImageJ from NIH - also in the sidebar - has some nice features for ROI analysis, too.

If you want the movies and jpeg images that appear above, these are available here:

http://dl.dropbox.com/u/26987499/Good_EPI_Cor_Sag_jpg_movies.zip

You can use the data, movies and jpegs for any educational purpose you like. No need to acknowledge the source unless you really want to, in which case please cite practiCalfMRI.blogspot.com and the Henry H. Wheeler, Jr. Brain Imaging Center at UC Berkeley.

3. Is the asymmetric L-R distortion a problem for coronal slices? Strictly speaking, it's no worse a problem than it was with axial slices! If you weren't overly concerned with distortion (and its correction) before, with axial slices, then it's a bit feeble to start worrying only when the problem becomes more conspicuous via asymmetry! That said, however, if you apply conventional coregistration algorithms in your analysis pipeline then it's possible that the disruption of left-right symmetry (and in particular the shearing) might cause worse results for coronal than axial slices. I couldn't tell you offhand. Distortion (and its correction) is a complex subject for another post. All I would say at this point is that I wouldn't be dissuaded by distortion alone from exploring a coronal prescription if there were measurable benefits in other domains, such as reduced dropout. I used to use coronal slices for some studies on our old Varian 4 T scanner where it was near impossible to get proper frontal coverage with axial slices. For most of you doing fMRI at 1.5-3 T, though, chances are the reduced coverage along the A-P axis would be the principal reason not to pursue coronal fMRI.

4. The amount of current induced in the subject's (electrically conductive) body is proportional to the cross-sectional area in the plane perpendicular to the switched gradient direction. If the current induced in the body becomes too large it is possible that peripheral nerves will be activated, causing twitching in the subject's muscles.

Given that the dominant gradient is the read (or frequency encode) gradient - the read gradient is generally larger and switched faster than the slice-select and phase encode gradients, as shown in the EPI pulse sequence diagram - then it is this gradient that concerns us the most from a safety perspective. Once the slice selection direction is established it leaves us with just two options for the read gradient direction, the other axis becoming the phase-encoded axis by default.

In understanding the safety limits for switched gradients it is useful to consider the body's three planes as if they act like pick-up coils; loops of wire that can sense changing magnetic fields by having an electric current induced in them. (It's a situation not unlike the RF coils we use in practice to detect the oscillating magnetization! But I digress.) Consider this cartoon showing the effective current loops formed in a subject's body when a gradient is switched in one of three cardinal axes:

The relative areas of effective current loops (in black) produced by gradient switching. An effective current loop is induced in the plane perpendicular to the switched gradient axis. The three principal switched gradient axes are anterior-posterior (A-P), left-right (L-R) and head-foot (H-F), corresponding to effective current loops in the subject's coronal, sagittal and axial planes, respectively.

So let's consider our options for a coronal EPI slice prescription. Now, the body's cross-sectional area in the plane perpendicular to the H-F axis, i.e. the subject's axial plane, is less than the cross-sectional area in the plane perpendicular to the L-R axis, i.e. the subject's sagittal plane. It's as if the axial plane contains a pick-up loop that's smaller than one in the sagittal plane, as illustrated by the black loops in the cartoon. The induced currents in the subject will be lower if we choose H-F for the read gradient direction assuming, as is virtually always the case for conventional fMRI, that the read gradient is larger than the phase encoding gradient. Thus, the phase encoding axis becomes L-R.

Suppose that you don't particularly like the asymmetric distortion of your coronal EPIs, and you decide to try and swap the read and phase encode gradient directions. What are the consequences? The good news is that you do indeed return the L-R axis of the EPIs to being distortion-free, and any stretches/compressions would now be H-F. But there's a price: we run a higher risk of peripheral nerve stimulation. Accordingly, the scanner's software reduces the gradient timing and amplitude parameters (i.e. the hardware performance is reduced) to maintain a safe situation for the subject, thereby reducing our ability to achieve a given spatial resolution and brain coverage (via increased minimum TE, for example) than were achievable in the previous situation, with the read axis H-F. What's more, by slowing down the gradients the minimum echo spacing actually increases, making the absolute level of distortion higher with the phase encoding H-F than it was when it was L-R. So much for symmetry!

That's coronal slices dealt with. A quick glance at the cartoon reveals the preferred read gradient axes for the other two cardinal slice prescriptions. For axial slices, the preferred read gradient direction is L-R, making the phase encode axis A-P. For sagittal slices the preferred read gradient direction is H-F, making the phase encode axis A-P as well.

Another note on the subject of peripheral nerve stimulation: note that the effective current loops we've considered assume that the subject has his hands by his side and his legs uncrossed. Linking hands/feet will increase the size of the effective current loop, especially in the coronal plane where the current loop was already largest! Crossed arms and/or legs may mean that the effects of switching the slice select or phase encode gradients - gradients that we ignored in determining which axis to use for the read gradient - may start to induce peripheral nerve stimulation in the subject. The scanner software assumes that you're not creating artificially large current paths! (I've done these tests myself, as it happens. We have a perfusion sequence with an aggressive EPI readout that is just about able to twitch the trapezius muscles in my back. But if I clasp my hands together I can get reliable twitching, no problem at all! The twitching stops the moment I release my hands.)

What's the bottom line, here? Well, other than the implications for (a)symmetric distortion that depends on the slice prescription, you should now be able to understand how a set of parameters that can be used to acquire EPIs in one orientation may be precluded (by the scanner) for another. What might have seemed a trivial step of attempting to swap the the read and phase encode directions actually has safety implications. Your subject sure appreciates it, even if you don't!

Understanding fMRI artifacts: "Good" axial data

2011-11-16T17:29:00.000-08:00

Good EPI data has a number of dynamic features that are perfectly normal once a few basic properties of the sample - a person's head - are considered. The task is to differentiate these normal features from abnormal (or abnormally high) artifacts and signal changes. We'll look at axial slices first because these are the most common slice prescription for fMRI. (Axial oblique slices will exhibit much the same features as the axial data considered here.)

The data we will consider in this post were acquired with a single shot, gradient echo EPI sequence on a Siemens Trio/TIM scanner, using the 12-channel head RF coil and a pulse sequence functionally equivalent to the product sequence, ep2d_bold. (See Note 1.) Parameters were typical for whole cortex coverage (the lower portion of the cerebellum tends to get cut off): 34 slices, 3 mm slice thickness, 10% slice gap, TR=2000 ms, TE=28 ms, flip angle = 90 deg, 64x64 matrix over a 22.4 cm field-of-view yielding 3.5 mm resolution in-plane, full k-space with phase encoding oriented anterior-posterior. (See Note 2 for advanced parameters.) The entire time series was 150 volumes in duration but in the movies and statistical images that follow I've considered only the first fifty volumes. (See Note 3 if you want to download the entire raw data and/or the movies and jpeg images.)

Let's start by simply looping through the volumes with the contrast set to reveal anatomy. Play this through a couple of times to familiarize yourself with it, then read on (click the 'YouTube' icon on the video to launch an expanded version in a separate tab/window):

Other than movement of the eyes and some large blood vessels in the inferior slices, at this resolution it's difficult to determine with certainty which regions are fluctuating and which are stationary. So let's zoom in on some of the central slices and replay the cine loop:

Now we can see that there's quite a bit of brain pulsation going on. Indeed, nothing appears stationary now! However, the edges of the brain don't appear to be moving very much so we can be reasonably confident that the pulsation is due to normal physiology and not a fidgety subject.

Now let's replay the same zoom but with the images re-contrasted to make the ghosts and the background noise regions visible:

The pulsations in the brain are still just about visible, and the edges of the brain still appear to be stationary. Good. Now, though, we can clearly see the N/2 ghosts and these are stationary also. Even better! This strongly suggests that the subject movement is very low. If the subject were moving even a little bit it would tend to perturb the magnetic field homogeneity over the head (the "shim") and cause the ghost intensity to fluctuate. Furthermore, if the subject had moved a lot - say, he sneezed - in the time between the shim and the acquisition of these EPIs, we might expect the ghost level to be high but temporally stable. The ghost level is acceptably low, also indicating an absence of significant movement since the start of the session. (See Note 3 in the post, Tactical approaches to (re)shimming for instructions on triggering a shim at any point during a scan on a Siemens.)

Also in this view we see that the background regions that aren't contaminated with ghosts tend to fluctuate in a random manner, as they should. Nothing coherent pops up, it's just speckly, grainy noise. Also good.

But what about the inferior slices? The quality of the shim will be considerably worse for the frontal and temporal lobes, and the inferior surface of the brain. Then there are the eyes to consider. So let's take a look at a zoom into some of these slices:

Yup, the eyes are dancing around alright! And the muscles around the eyes generate some dynamic signal, too. Normal physiological fluctuations, especially large blood vessels (white spots) are present, but otherwise the edges of the brain appear relatively stationary.

The ghosts were already visible in the previous contrast, but let's re-contrast and highlight the background, as we did previously:

The brain's edges remain stationary by inspection - good - and the N/2 ghosts also appear to be quite consistent - also good. And finally, the background noise is relatively uniform - no obvious structure. If you're wondering what the faint horizontal streaks are that are especially pronounced in the bottom row of slices, it's a consequence of the interaction of a sharp contrast boundary and ramp-sampled k-space. I'll deal with this artifact in detail in a separate post. At this point you should content yourself with it being a normal feature of EPI.

Simple statistical images

The sort of zooming and contrasting that I've done with these movies is the sort of thing that you can (and should) do with your scanner's inline (real time) image display tool. Rather than catching up on your fiction reading while your fMRI experiment proceeds, spend the time watching the real time display to ensure that the images' appearance remains similar to what you've just seen.

To reinforce the features in the movies I want to change gears and show you how this fifty-volume data set appears in simple statistical images. Unless your scanner has real time analysis you're not going to be able to look at summary statistical images like these on-the-fly, but it is instructive to see how the fluctuations you've seen in the movies translate into statistical features. We will also be in a position to make an assessment on the likely impact on fMRI statistics because we're quantifying the background fluctuations that your experiment has to overcome to reach significance. (See the post, Comparing fMRI Protocols for detailed explanation on the usefulness of the following images a proxy for fMRI performance.)

Consider the temporal SNR (TSNR) image, which is the pixelwise mean divided by the pixelwise standard deviation of the fifty image series:

TSNR image. (Click to enlarge.)

Let's also consider the standard deviation (SDEV) image. It helps confirm what we think we're seeing in the TSNR images. (I am in the habit of using these as a pair for diagnostic purposes. It means I get to see the entire dynamic range of fluctuations.)

Standard deviation image. (Click to enlarge.)

The eyes are easily the brightest feature in the SDEV image - all those saccades we saw in the movie! - and this translates into relatively low TSNR. Favorably, however, there are no intense ghosts arising from the eyes apparent in the SDEV image. If there were, we should expect to see twin stripes (or pairs of oval-shaped ghosts) in anterior-posterior planes that encompass the eyes. (Recall from physics Part Twelve that the ghosts will be located in the A-P direction with the phase encoding axis oriented A-P.) We will look at ghosts from eye movement later on in this series.

Another bright feature in the SDEV image is large blood vessels, particularly in the lower slices (towards top left in the matrix view). This is good news: one hopes the subject has a cardiac cycle! These bright spots are perfectly normal, and acceptable.

The next most intense features in the SDEV image are the edges of the brain and regions that contain CSF, the ventricles and sulci. This is consistent with cardiac-driven pulsatile flow of CSF, as well as pulsations of the brain tissue in general. We might have predicted high SDEV (and low TSNR) for CSF-filled spaces because of the contrast between gray/white matter and CSF in the original EPIs; CSF appears quite a lot brighter than tissue. Coupled with the expected pulsations, the features make sense in the statistical images.

Looking more closely at the outside edges of the brain in the SDEV image, there is also a subtle enhancement of the posterior and anterior edges compared to the left and right edges. This could be head movement - movement in the chin-to-chest axis is most easily achieved by a subject for many reasons, including (especially) swallowing - pulsatile flow of CSF in the meninges, or it could be due to a small amount of gradient heating during the acquisition, because such heating tends to look like motion in the phase encoding direction. It's difficult or impossible to tell between these three explanations with these simple statistical views. What's important here is that this sort of edge enhancement is absolutely normal and should be expected. (Realignment algorithms, as typically used in an fMRI post-processing stream, can be expected to eliminate much of the reduced TSNR in the anterior/posterior edges.)

If you can negotiate your gaze around the CSF-filled spaces you may be able to make out that gray matter regions have higher standard deviation than white matter regions, producing correspondingly lower TSNR for gray matter. This is entirely predicted by neurophysiology, because we know that GM is many times more metabolically active than WM. Blood vessel density, to support the higher metabolic rates, is higher in GM and this translates into larger signal fluctuations in our BOLD-contrasted time series. Welcome to resting-state fMRI! It's precisely this "ongoing activity" that we map with that technique.

Now we're down to the lowest of features in these "good" statistical images, so I'll include a labeled composite of the SDEV and TSNR images to help you out:

(Click to enlarge.)

Note the circles of scalp fat around the head, visible most prominently in the TSNR image. Although a fat suppression pulse is applied before each EPI slice, the suppression is imperfect and tends to leave a small residual signal. But the good news is that these fat circles appear only weakly in the SDEV image, suggesting that head movement and/or gradient heating effects were low during the acquisition. Indeed, the fat circles in the SDEV image get brighter the more superior the slice. This is consistent with small movement (head rotation) in the chin-to-chest direction, as can happen from slow compression of foam padding and/or relaxation of neck muscles. Swallowing, fidgeting and other head motion tends to produce much larger deviations, as we'll see in later posts.

Finally, ghosts again. Note that the N/2 ghosts from the majority brain signal are low intensity; down towards the background noise level (where they should be). The most intense ghosts in the TSNR image are due to the residual fat signal I just mentioned. This is expected because fat protons resonate several hundred Hz away from water protons, meaning that there will always be an accrued phase difference between fat and water signals that will manifest as increased ghosting for whichever of the two species is placed off-resonance. Convention (and common sense) suggests placing water on resonance, because that's the signal of interest in the brain tissue, and that means that the residual subcutaneous scalp fat will produce (unavoidable) N/2 ghosting. This is normal; the only way to reduce fat ghosts is to improve fat suppression, which is a separate topic for another day.

Next post: "good" coronal and sagittal EPIs.

__________________________

Notes:

1. At Berkeley we use a modified version of ep2d_bold, called ep2d_neuro. In these tests there is no functional difference between the two sequences. We have our own local default pulse sequence so that we can have, if desired, thinner slices, a user-defined number of dummy scans, 10 microsec precision in TR, and some other relatively minor tweaks. If you want to replicate these tests then simply set up ep2d_bold with the same parameters as used here.

2. Bandwidth = 2790 Hz/pixel, echo spacing = 0.47 ms. These days I don't tend to drive the echo spacing below 0.5 ms without very good reason, but for the purposes of this post - assessing normal variations in the signal and ghosts - these settings may be considered optimal. I don't want to get deep into the effects of echo spacing here, but a short version is that setting the bandwidth and echo spacing (as a pair) is an exercise in avoiding mechanical resonances and excessive ramp sampling that can enhance any fluctuations in the power feeding the gradient amplifiers. I'll deal with these issues in depth in dedicated posts. These are two of the artifacts we want to become intimately familiar with!

3. Want the raw data from this post? You can download a zip file containing all the DICOM images here:

http://dl.dropbox.com/u/26987499/Good_EPI_Axial_64x64.zip

If you don’t already have a DICOM viewer, check out Osirix for Mac OSX (available via a link in the sidebar). ImageJ from NIH - also in the sidebar - has some nice features for ROI analysis, too.

If you want the movies and jpeg images that appear above, these are available here:

http://dl.dropbox.com/u/26987499/Good_EPI_Axial_jpg_movies.zip

You can use the data, movies and jpegs for any educational purpose you like. No need to acknowledge the source unless you really want to, in which case please cite practiCalfMRI.blogspot.com and the Henry H. Wheeler, Jr. Brain Imaging Center at UC Berkeley.

Understanding fMRI artifacts

2011-11-15T20:36:00.000-08:00

Introducing the series

The workhorse sequence for fMRI in most labs is single-shot gradient echo echo planar imaging (EPI). As we saw in the final post of the last series, EPI is selected for fMRI because of its imaging speed (and BOLD contrast), not for its ability to produce accurate, detailed facsimiles of brain anatomy. Our need for speed means we are forced to live with several inherent artifacts associated with the sequence.

However, in addition to the "characteristic three" EPI artifacts of ghosting, distortion and dropout, when we're doing fMRI we are more concerned with changes over time than with the artifact level of an individual image. So, in this series we need to assess the sources of changes between images, even if the images themselves appear to be perfectly acceptable (albeit subject to the "characteristic three").

What's the data supposed to look like?

It would be rather difficult for you to determine when something has gone wrong during your fMRI experiment if you didn't have a solid appreciation of what the images ought to look like when things are going well. Accordingly, I'll begin this series with a review of what EPIs are supposed to look like in a time series. We'll look at typical levels of the undesirable features and assess those parts of an image that vary due to normal physiology. This is what we should expect to see, having taken all reasonable precautions with the subject set up and assuming that the entire suite of hardware (scanner and peripherals) is behaving properly.

Good axial data will be the focus of the first post in the series. (Axial oblique images will exhibit qualitatively similar features to the axial slices I'll show.) In the second post I'll show examples of good sagittal and coronal data. Artifacts may appear quite differently and with dissimilar severity merely by changing the slice prescription, so it's important to keep in mind the anisotropic nature of many EPI defects. Motion sensitivity is also different, of course. Motion that was through-plane for an axial prescription is in-plane for sagittal images, for example.

Ooh, that's bad. Is it...?

With a review of good data under our belts it will be time to look at the appearance of EPI when things go tango uniform. I will group artifacts according to their temporal behavior - either persistent or intermittent - and their origins - either from hardware, from the subject, or from operator error. You should then be able to understand and differentiate the various artifacts and be able to properly diagnose (and fix) them when it counts the most: during the data acquisition. Waiting until the subject has left the building before finding a scanner glitch is a bit like doing a blood test on a corpse. Sure, you might be able to determine that it was the swine flu that finished him off, either way he's dead. Our aim will be to do our “blood tests” while there is still a chance of administering medicine and perhaps achieving a recovery.

Defining what we mean by "good" data

One of the hardest tasks facing someone new to fMRI is the ability to recognize when data is ostensibly ‘good.’ On the face of it such a problem might seem strange, even if we recognize from the outset that ‘good’ and ‘bad’ are subjective assessments. So let’s take a closer look at these definitions.

For the most part what you, as an experimenter, mean when you say that data is good is that the fMRI experiment showed the sensitivity to the stimulus-induced hemodynamic changes you expected, i.e. that a particular threshold of statistical significance was attained. In short, your experiment worked. However, what an MRI physicist usually means by ‘good' data is often subtly and crucially different.

To a physicist, ‘good’ data is obtained whenever the EPI time series yields images that have an appearance and quality that can be reasonably expected on that particular scanner and when using that particular set of parameters. Note that this assessment says nothing about the appropriateness of the time series for your experiment. You might have chosen a bum set of parameters for your intended application. For example, you might have been expecting activity in the inferior prefrontal cortex, but you went and used a TE of 40 ms to ensure lots of functional contrast in occipital lobe. Oops. Now you have a hole in the image where you wanted signal. Bad data, or bad experiment...?

So, to a physicist, ostensibly 'good' data is achieved any time the scanner didn’t introduce any (strong) artifacts such as global signal intensity drifts, when the peripheral equipment didn’t introduce any (large) RF noise, when the parameter selection avoided any settings not commensurate with maximum scanner performance (such as might be needed to avoid strong mechanical vibrations), and when subject motion was in some as yet to be defined ‘normal range.’ In other words, we are interested in the quality and stability of the images, not whether the time series was able to answer your question. See the difference?

There are many, many reasons why your experiment might not have worked as expected even though the images had low artifacts and were temporally stable. Perhaps your voxels were too large for the particular brain regions you are interested in. How is that the fault of the scanner? It's entirely possible to acquire good quality data with a sub-optimal protocol then implicate the data for failing to produce the magical colored blobs. But it is a planning problem, not a data quality problem.

Winning the data lottery

Conversely, you may have used the ideal protocol but been sloppy with your subject's head restraint and obtained a positive result in spite of unnecessarily high motion. (Congratulations! Go buy a lottery ticket!) To a physicist, however, the fact that you managed to dodge lots of bullets doesn’t imply that your data was ‘good.’ When we are discussing the quality of a particular time series acquisition we are saying absolutely nothing about the experiment per se; we are simply going to determine whether that particular time series could, under normal circumstances, have been acquired better.

Post hoc analysis

Now, you might ask why we are even bothering to look at EPIs directly. Why not just run a script on the time series and have it spit out some sort of quantitative assessment? (See Note 1.) That is a great idea! But it is generally only useful to do after the time series has finished, often requiring the data set to have been ported off the scanner. What we would prefer to do is use our expertise and intuition to spot the moment something isn’t right during an EPI run – during the experiment - then quickly determine whether the problem has righted itself or will persist and necessitate some sort of remedial action on the part of the experimenter. And the fastest way to do that is to observe the data as it acquires, in real time. Of course, if you are fortunate to have a scanner or software that can analyze data in real time that is an added bonus to using your skills as an observer. Most labs don’t have such facilities yet, however.

Your analysis scripts will be extremely useful to quantify the data quality after the experiment, perhaps just seconds after the end of one time series acquisition and before the next one. Therefore, we will take a quick look at a couple of simple statistical images that can be produced on most scanners’ software, in the midst of your scan session. Unless you are fortunate to have a scanner that produces such statistical images in real time you won’t be able to produce these maps continuously through the session, but you might be able to assess each task block, say, and assure smooth progress, repeating suspect or definitively bad runs. Remember, the more ways you can assess your data during the scan session, the more likely you are to detect (and be able to remedy) a problem!

Contrasting your images for artifact recognition

We're very nearly ready to get down to business. In the posts to come I will routinely review the same images with different contrast levels. Let this segment serve as a general explanation of what I'm doing and why.

When you're hunting for artifacts it is a good habit to assess your EPIs with at least two different contrast settings. You want to assess the images with contrast as they might be used in a publication - what I'll term anatomically contrasted images - so that you can determine whether the appearance of the brain is acceptable. But you then want to review the images having brought up the background "noise" region to a visible intensity, to reveal the N/2 ghosts and any problems that may be lurking in the crud. Many artifacts in EPI are furtive and are often best detected at the level of the image background, where there should be only noise plus the ubiquitous N/2 ghosts. Indeed, the ghosts themselves can be used as sensors for problems, as you'll see when we encounter several of the artifacts later on in this series.

As a very rough rule of thumb, then, large errors, such as sizable subject motion or intense gradient spiking, will be highly visible in anatomically contrasted EPIs; that's the first place to look. Small subject motion and subtle scanner hardware issues, such as mechanical resonances, may only be visible in background-contrasted images. You could (perhaps should) even make this algorithmic:

Step 1 - contrast the anatomy and check that the brain appears as you expect it to look. Large holes? Stripes? Signal moving all over the place? If the answer to all of these is "no," proceed to step 2.

Step 2 - contrast the background noise to reveal the N/2 ghosts. Inspect the ghosts. Is their intensity acceptable? Are they stable over a few TRs? If all is well here, proceed to step 3.

Step 3 - with the image still contrasted to reveal the ghosts, inspect regions in the image that should be ghost-free, i.e. noise. Is it noise, or does anything with structure pop in and out of existence?

Okay then, that's how we're going to do it. Let's get to work by assessing the fluctuations that you can expect to see in "good" data. See you at the next post!

_________________________

Notes:

1. I'm sure there are dozens of useful post hoc diagnostic tools out there. One suite that I'm partially familiar with is courtesy of Matthew Brett from his Cambridge days. Another suite that looks useful but that I've not tested myself comes from the Gabrieli Lab. (The parent links to these sites are available in the right sidebar of the blog.)

Physics for understanding fMRI artifacts: Part Twelve

2011-11-01T20:25:00.000-07:00

Apologies for the lengthy delay getting this post out. New academic year, teaching, talks, etc. etc. Anyway, I hope that this opus will be the final post in the background physics series for the time being. I reserve the right to append further posts down the road, but with this post I hope you will be in a position to understand the origins of artifacts in real (EPI-based) fMRI data. So, after today we'll change tacks and start reviewing what "good" data should look like. First things first though. Time to put all your k-space knowledge to good use, and review the pulse sequence that the majority of us use for fMRI.

The Echo Planar Imaging (EPI) pulse sequence

In Part Ten we looked at a pulse sequence and its corresponding k-space representation for a gradient-recalled echo (GRE) imaging method. That sequence used conventional, or spin warp, phase encoding to produce the second spatial dimension of the final image. A single row of the k-space matrix was acquired per RF excitation, with successive rows of (frequency-encoded) k-space being sampled after stepping down (or up) in the 2D k-space plane following each new RF pulse.

One feature of the spin warp imaging scheme should have been relatively obvious: it's slow. Frequency encoding along kx is fast but stepping through all the ky (the phase-encoded) values is some two orders of magnitude slower, resulting in an imaging speed from tens of seconds (low resolution) to minutes (high resolution). That's not the sort of speed we need if we are to follow blood dynamics associated with neural events.

Instead of acquiring a single row of k-space per RF excitation - a process that is always going to be limited by the recovery time to allow the spins to relax via T1 processes - we need a way to acquire multiple k-space rows per excitation, in a sort of "magnetization recycling" scheme. Ideally, we would be able to recycle the magnetization so much that we could acquire an entire stack of 2D planes (slices) in just a handful of seconds. That's what echo planar imaging (EPI) achieves.

Gradient echo EPI pulse sequence

The objective with the EPI sequence, as for the GRE (spin warp) imaging sequence we saw in Part Ten, is to completely sample the plane of 2D k-space. That objective is unchanged. All we're going to do differently is sample the k-space plane with improved temporal efficiency. Then, once we have completed the plane we can apply a 2D FT to recover the desired image. Pretty simple, eh?

As before, sampling (data readout) need only happen along the rows of the k-space matrix, i.e. along kx. So we need a way to hop between the rows quickly, spending as much time as possible reading out signals under the frequency encoding gradients, Gx, and as little time as possible getting ready to sample the next row. EPI is the original recycled pulse sequence, so I'll color the readout gradient echoes in green:

The first four (and a half) gradient echoes in a gradient echo EPI pulse sequence.

To keep things simple I've omitted slice selection and indicated a 90 degree RF excitation; this could of course be any flip angle in practice. (See Note 1.) I've also shown just the first four (and a half) gradient echoes in the echo train. The full sequence repeats as many times as there are phase-encoded rows in the k-space matrix. A typical EPI sequence for fMRI might use 64 gradient echoes, corresponding to 63 little blue triangles in the train shown in the figure above. But for the example k-space plane below, the k-space grid is 16x16 so assume for the time being that the full echo train would consist of 15 little blue triangles separating eight positive Gx gradient periods and eight negative Gx gradient periods.

Before we consider the 2D k-space representation, let's take a moment to assess a few features of the pulse sequence. First, note that the negative Gx period (in orange) at the start of the sequence is designed to balance half of the subsequent positive Gx period (in green). (See Note 2.) This is a gradient echo, exactly as we saw in Part Eight (and again in Parts Nine and Ten). Furthermore, after the first positive Gx gradient the sign alternates for subsequent Gx periods, meaning that we will end up with a "train" of gradient echoes, each one being sampled (read out) with frequency information. In the absence of any other gradients, all these gradient-recalled echoes would contain essentially the same one-dimensional spatial information. So let's look at the second dimension via the k-space plane.

Here's the k-space representation for a 16x16 matrix, with the trajectory color scheme corresponding to the gradients in the figure above:

You don't see the individual data samples under the read gradients, these are assumed. And as a general rule the dimensions of the matrix can be set somewhat arbitrarily. I say "somewhat" because as the number of data samples under each read gradient is increased, the time between successive phase encode gradients (the blue triangles) increases and so does the overall duration of the pulse sequence. As we will see below, there are practical limits to the amount of time that can be spent sampling the k-space matrix in EPI and you must be judicious in how you spend that time.

Let's add another parameter to the pulse sequence. In Part Ten when we looked at the regular (spin warp) GRE imaging sequence I didn't label the effective echo time, TE. However, I did indicate the point at which the phase imparted by the read gradient is returned to zero, and I mentioned that this point corresponds with a journey through the k-space origin. It is this journey through the k-space origin that we use as our definition of TE.

In both EPI and GRE, and for MR pulse sequences in general, we define the TE as the time corresponding to the interval between the creation of transverse magnetization, i.e. from the center of the RF excitation pulse, to the moment the k-space trajectory passes through the k-space origin. This isn't an arbitrary designation, even though there are clearly signals (indeed, the majority of signals!) being acquired at times other than at TE. Signal intensity at the k-space origin should be maximal - assuming the imaging gradients dominate any artifactural gradients in the sample. Then, because variations in signal intensity define the image contrast, and it is image contrast that is usually of primary interest/utility (edges and small details tend to be secondary in most imaging applications, not just fMRI) it is a good compromise timing at which to understand the overall image features.

Here is TE indicated on a complete EPI pulse sequence comprising 32 total echoes:

Courtesy: Karla Miller, FMRIB, University of Oxford.

The TE is of vital importance for understanding EPI artifacts and for setting up BOLD contrast. But it is also worth remembering that most of the signals acquired in the echo train are not actually acquired at TE but at differing offsets from it. As we saw in Part Eleven, high spatial frequencies appear at the extremes of k-space, so we can already see that the early and late echoes in EPI - which are the echoes containing edge information in the phase encoding dimension - will have markedly different contrast properties than the signals towards the k-space center. More on the practical manifestations in the artifact recognition posts to come (but see Note 3 if you fancy trying to predict what we might observe experimentally).

Processing 2D k-space for EPI

So I may have lied. Slightly. It transpires that there's an additional step necessary before the 2D FT can work its magic and obtain an image from the 2D k-space presented above. Did you notice how the green arrows point in alternating directions in the k-space trajectory? In effect, the way we have sampled alternating rows of kx-space is akin to time going forward then backward from the spins' perspective. Phase accrual, which as we know from Part Ten is how the spatial encoding actually happens, appears to be running backwards on alternating rows. We need to reverse the direction of all of the odd (or all of the even) lines to make the phase evolution consistent:

Left: original k-space matrix acquired by the single-shot EPI sampling trajectory. Right: the matrix after reversing alternate rows of k-space.

The necessity of reversing all the odd (or even) rows of k-space prior to the 2D FT would be a mere processing formality were it not for something that you may have picked up on in the previous paragraph. I said that the phase accrual under the readout gradient appears to be forwards in one kx row, then backwards in the next, then forwards again, and so on through the matrix. If the only sources of phase shift were the imaging gradients, and if these gradients were ideal, we'd be in good shape. However, any additional (or erroneous) phase shift, e.g. arising from magnetic susceptibility gradients in the sample, or from imperfections in the imaging gradients themselves, will cause a problem. (See Note 4.)

By way of example let's look at one source of spurious phase shift that we can easily represent in a diagram you've seen before. Consider the trusty one-dimensional gradient echo sequence that we first saw in Part Eight. However, in addition to the imaging gradients that we can control, let's add another background gradient (that we can assume arises from magnetic susceptibility differences across the sample) in the same direction, in this case the x axis. With the imaging gradient and background (susceptibility) gradient in the same direction we can instead consider just the effective gradient, Geff:

An illustration of the effect of a background (magnetic susceptibility) gradient on the echo refocusing time for a gradient echo sequence. The (resultant) effective gradient is shown in orange.

For simplicity I have suggested that the background gradient is linear, just like the imaging gradient, Gx. (Real background gradients can and will have complex spatial dependencies.) And, although real background gradients persist indefinitely - as long as the subject is in the magnet they are "on" - I've draw the background x gradient only for the period when it coincides with the Gx imaging gradient pulses, because that's when it has its effect on the k-space problem we're considering. What happens before the RF pulse or after the data sampling (analog-to-digital conversion) period is of no consequence for us here.

Clearly, in the absence of the background gradient the echo would refocus at the time of the dashed line, as we saw in previous posts. But the additional positive term in Gx actually causes the echo to refocus late, at the dotted line, because that's the time at which the effective gradients, in orange, are balanced.

What does this late-arriving echo mean for our k-space representation? Well, we intended the center of k-space to be at the dashed line position, whereas the bulk of the signal and the actual position of zero phase (the echo top) corresponding to the center of k-space happens later, at the dotted line. If we were doing a one-dimensional FT to get a 1D profile this isn't a big deal; we accrued a bit of spurious phase but the echo still occurred inside the period of signal detection. (We nearly always look at magnitude images so a little phase twist in the profile would be eliminated.) Likewise, if we were doing a conventional 2D spin warp phase encoding for our 2D image it's another case of having a constant phase offset in each kx row that essentially "falls out" of the 2D FT with no major consequence for the resulting image. (Again, taking the magnitude of the image removes any phase twist.) But EPI is different.

In EPI we use a train of readout gradient episodes, first positive, then negative, then positive again, etc. We are reversing the sign of the readout gradient, but of course the background gradients are unchanged! Any background gradient that is additive to positive read gradients will be subtractive from negative read gradients, causing the echoes collected each time to be either early or late relative to the expected echo tops. (I hope it's clear that if we reversed the signs of the Gx imaging gradients in the previous figure, leaving the background gradient as it was, and then reconsidered the effective Gx gradient, that the echo top would arrive earlier than the dashed line.) The exact ordering - additive or subtractive, early or late - isn't really an issue. All we need to recognize is that we have an alternating pattern for each successive gradient in the EPI readout echo train. And it is this alternating pattern that causes the problem for the resulting image.

We will consider in detail the effect of the alternation in the section on ghosting, below. For now, let's just complete the picture that we get in our 2D k-space for EPI in the presence of a small background gradient. We're going to have a zigzag offset across the ky dimension of k-space, one that persists after reversal of alternating k-space rows:

Left: original k-space matrix with an offset, in this case a small delay at the start of each sampling period. Right: the same matrix after reversing alternate rows of k-space. A zigzag across the phase encoding dimension is now clearly apparent.

In this figure the green arrows show the effective k-space matrix overlaid on the ideal matrix. But we're not FTing the ideal matrix, we're FTing the zigzagged k-space matrix in green on the right. More on the practical consequences of doing this in the next section. (See Note 5.)

EPI artifacts

You never get something for nothing. In exchange for the speed of EPI we pay a price in terms of image quality. There are essentially three artifacts that might be considered characteristic of the pulse sequence: ghosting, distortion and dropout. Strictly speaking, the latter - dropout - isn't a characteristic of EPI per se but is instead a consequence of using a relatively long TE (compared to the brain's T2*) to acquire the image. Another gradient echo sequence, such as a conventional (spin warp) phase-encoded gradient echo sequence, would suffer from similar signal dropout were it to be acquired at the same TE.

Anyway, let's look at the three artifacts that will plague every EPI you ever acquire for fMRI and consider in brief some of the factors contributing to the level of artifacts thus produced. One point to keep in mind as you read is that the severity of each artifact will nearly always vary across the image; these aren't usually constant, global effects. The regional variation is a function of the magnetic susceptibility gradients across the sample. The larger the heterogeneity of the magnetic field the stronger the (local) artifacts and the more anisotropic will be their spatial pattern. Put another way, some parts of the brain - occipital and parietal cortices, for example - can be imaged relatively free of major artifacts whereas other parts of the brain - mid-brain structures, frontal and temporal lobes - will be plagued by them.

Ghosting

Often called N/2 ghosts or Nyquist ghosts (for reasons that will become clearer later), these artifacts are an unavoidable consequence of the back and forth sampling trajectory in k-space. Let's look at the appearance of the ghosts, in particular their spatial location, before examining their origins further.

Consider this EPI through a spherical phantom:

The position of N/2 ghosts - a pair of faint "shadow" images - on a single-shot echo planar image. The center of the phase encoding dimension is indicated with a dashed line. The ghosts are displaced by half the FOV from the center.

The first thing to note is the dimensionality of N/2 ghosts: ghosting occurs in the phase encoding dimension only, because that's the dimension that suffers from the zigzag errors in 2D k-space, as we've seen above. The frequency encoding dimension (left-right in the image above) doesn't show ghosts.

Next, note how the ghosted images appear displaced by half the field-of-view (FOV), with the bottom half of the (real) image generating a ghost image that appears at the top of the FOV and the top half of the (real) image generating a ghost image that appears at the bottom of the FOV. Single-shot EPI ghosts always appear at these locations, which is why they are often referred to as N/2 ghosts, where N indicates the total number of pixels in the phase-encoded dimension. (Clearly if the FOV is defined by N pixels then N/2 defines half the FOV.)

You should also note that portions of the ghosts overlap with the real image. As we will see when we look in-depth at real data, such overlap can have profound implications for your fMRI statistics because, as a general rule, anything that perturbs the time series data, such as subject motion, will affect the ghost intensity, thereby increasing the variance of the ghosts above that for real signal regions. So, if you have a ghost parked on your most critical anatomical region, expect to get crappy statistical power! More on tactics to reduce and relocate ghosts in future posts. (See Note 6.)

A final comment on ghost appearance: the magnitude of N/2 ghosts can vary considerably, but a rough rule of thumb is that they should be no more than 5% of the intensity of the main signal regions on a well-behaved scanner with good parameter selection and minimal subject motion. It's quite feasible to get ghost levels as low as 1% with a modicum of effort.

Okay, now that you know what they look like, let's go back to the zigzag in k-space and get a basic understanding of why the ghost images arise. In essence, we can consider the final echo planar image that emerges from the 2D FT as the sum of an ideal image with an artifactual (ghost) image, the latter arising purely from the offsets producing the zigzag across k-space. Here's a cartoon of this process:

The k-space for the ghosts (right panel) comprises erroneous "zig" offsets that sum with the odd rows of the ideal k-space (center panel), plus "zag" rows that are just zeros, i.e. no signal or noise, as indicated by white dotted lines. Together, the zigzagged erroneous k-space (right panel) adds to the ideal k-space (center panel) to yield the actual k-space (left panel).

Now, since we know from Part Eleven that an image's FOV is defined by the k-space step size, it should be intuitive that the erroneous k-space step in the right panel is effectively 2*delta-ky (where delta-ky is the actual k-space step size in the phase encoding dimension) and will therefore produce an image that has half the FOV of the ideal image. (Bigger step in k-space generates smaller FOV in image space.) In other words, the ghost signals arising from the zigzag in k-space appear at positions appropriate for an image with half the FOV of the actual FOV. The ghosted signal ends up misplaced by that amount - half the FOV - in the final image.

So, what can we do about these pesky ghosts? We first spend some time minimizing the experimental processes that lead to the zigzags in k-space, then we usually use a correction step that involves a handful (three is typical) of "navigator" gradient echoes that are inserted into the pulse sequence between the RF excitation pulse and the readout echo train. These navigator echoes are designed to capture information on the zigzags that can be expected in the readout echo train, allowing a post-processing phase correction to be applied to the 2D k-space before the Fourier transform. The correction step is most likely built-in to your scanner's software, so it's not something that you have to worry about unless the default correction proves insufficient for some reason. (See Note 7.)

Distortion

Spatial distortion is probably the most insidious of the EPI artifacts. I label it as insidious because it isn't until the effects become grotesque that most fMRIers acknowledge the problem, even though it's present in all EPIs. Want a good example? Go back and look at the EPI of the sphere above, the one showing the position of the N/2 ghosts. Does that bright white thing in the center of the image look like a circle to you? That's what it's supposed to be because it's a slice through the center of a sphere. Doesn't it look more like a section through a prolate spheroid? Bear in mind that this image was obtained from a phantom that can be shimmed better than any head! (Spheres have infinite axes of symmetry so it's relatively straightforward to homogenize the magnetic field with three first order and five second order shims that you find on most modern scanners.) In prosaic terms, the level of distortion in your real brain data is going to be worse than what you've just seen in the phantom.

So where does the distortion come from and what can we do about it? It's a consequence of the "recycling" of magnetization through a gradient echo train, i.e. it is due to the time it takes to produce the entire data readout. As fast as EPI is it's not instantaneous, meaning that background susceptibility gradients (yet again!) can make their unwanted presence felt. Let's look in more detail at the effect by first considering the situation in conventional (anatomical) imaging, where distortion can usually be safely ignored.

In conventional anatomical imaging we acquire one line of (frequency-encoded) k-space per RF excitation, stepping through successive k-space rows with a new RF excitation and a different phase encoding gradient value. Taking the frequency encoding axis first, then, the bandwidth (the spread of the spatially-encoded frequencies imposed by the gradient) of that dimension is given by the inverse of the time between data sampling points (the so-called "dwell time"). This is equivalent to our previous observation in k-space, that small steps in k-space define the image FOV; all that differs here is the conjugate pair under consideration, i.e. time and frequency rather than k-space and real space. So, the time between individual readout data points under the frequency encoding gradient defines the entire bandwidth of that dimension (and also the FOV, except that we're not thinking in terms of k-space at the moment).

What's a typical dwell time, hence bandwidth, for a readout dimension? On a modern scanner with typical gradients we might acquire as fast as 5-10 microseconds per point. Thus, a typical readout dimension bandwidth is the reciprocal of this range, i.e. 100-200 kHz. Now consider the effects of background magnetic field heterogeneities on this dimension. If I tell you that a typical "bad" region in the brain produces errors in the main magnetic field of around 1-3 parts per million (ppm), this corresponds errors of one to three millionths of the operating frequency of the scanner, or about 100-300 Hz at 3 T. Compared to the bandwidth of the readout axis, that's pretty small. Let's suppose we acquired 128 data points under the readout gradient. That means each pixel in the readout dimension is 128th of the bandwidth, or between 800 and 1,600 Hz per pixel. At most the magnetic susceptibility gradient causes a spatial error that's a third of a pixel. Pretty tiny.

Now consider the conventional (spin warp) phase encoded axis in a non-EPI scan, which we saw in Part Ten. Here the situation with respect to distortion is even better! Although the bandwidth in the phase encoding dimension can be defined as the reciprocal of the time underneath the phase encoding gradient episode (given by tp in Part Ten), there's a subtle difference in the way the magnetic susceptibility gradients manifest in the image. Note that if the magnet heterogeneities remain fixed throughout the acquisition of the image - as they should - then each row of k-space will possess precisely the same erroneous phase shift. This constant phase error effectively passes through the Fourier transform, adding a net phase twist to the image that is neatly removed by the expedient of looking at a magnitude image. In effect, we can consider the bandwidth of the conventional phase encoding axis as being infinite, rendering negligible the 100-300 Hz misplacement we saw above. The conventional (spin warp) phase encoding dimension is distortion-free.

And then there's poor old EPI. We have the opposite situation than with conventional phase encoding. Instead of acquiring one phase encode step per RF excitation, we're going to acquire every single phase encoding blip in succession following a single RF excitation pulse. The time between each phase encode blip, which is the echo spacing as shown in the first figure of this post, defines the bandwidth of the phase encoded axis. A typical echo spacing might be 0.5 ms, yielding a bandwidth of 2,000 Hz defined by, say, 64 pixels. That corresponds to around 30 Hz per pixel. If we have a 100-300 Hz error arising from magnetic susceptibility gradients then we can have a spatial distortion in this dimension that could be 3-10 pixels in magnitude. For a typical pixel of 3-4 mm that is quite a lot if displacement!

Okay, let's recap. In EPI we have a frequency encoding dimension that suffers only modestly from distortion - it will generally be at the sub-pixel level - and we have a phase encoding dimension that suffers extensively from distortion, to a level of several pixels. I'd like to emphasize that the distortion of the phase encoding axis isn't global, rather it is a regional effect that is determined by the local magnetic susceptibility gradients. And, because these background gradients aren't usually linear (indeed, they are nearly always very complex with high spatial order) we don't have a neat linear distortion that is easily remedied. The following example should immediately convince you that this is the case:

The effect of a long echo spacing (left) and halving the echo spacing (right). Halving the echo spacing also halves the distortion in the phase encoding dimension (which is anterior-posterior here).

The weirdly distorted frontal lobe regions are stunningly obvious in the image on the left, especially. While we can't determine by inspection precisely what spatial variations gave rise to the spikes and troughs in the brain signal, it's readily apparent they weren't linear. What's more, the distortion can be a compression or a stretch, depending on the relative sign of the susceptibility gradients in the sample and the phase encoding gradient polarity. In the above image the frontal lobe regions are being mostly, but not exclusively, stretched.

The non-linearity and a combination of stretches and compressions makes fixing distortion a tricky proposition. There are methods, typically involving some sort of magnetic field mapping (to determine the susceptibility gradients that gave rise to the distortion) that may, under certain circumstances, be able to provide a rudimentary correction. I'm not going to get into this issue here, except to note two important issues: generally, it is easier to relocate (distortion-correct) pixels that have been stretched from their correct location than it is to relocate pixels that have been compressed. This is because compression causes displaced pixels to coalesce, massively complicating the mathematics of correction. Secondly, all distortion correction schemes produce only approximations to the anatomically-correct pixel locations. There are numerous assumptions in the methods that make the approximations more or less valid. There is as yet no robust, foolproof way to accurately correct distortions, which is presumably why the use of distortion correction isn't as widespread in the fMRI literature as you might expect. Instead, many fMRIers try as best they can to minimize the distortion at the acquisition stage, then let non-linear warping algorithms (spatial normalization) tackle the distortion at the same time as one is trying to morph the subject's brain anatomy into a standard space (such as Talairach or MNI coordinates). The whole process is complex and fraught, so I'll deal with distortion correction in a future post.

A final note on distortion. It may be easier to minimize distortion than to even try and fix it. As the above example images show, by halving the echo spacing duration the level of distortion is also halved. Of course, there are hardware limits to how short we can make the echo spacing, so other alternatives might use parallel imaging methods (e.g. GRAPPA) or segmented (multi-shot) EPI, to reduce the distortion level. But, as you have come to expect (I hope!), these new methods bring with them their own complexities and artifacts, so it's not like we can simply avoid the distortion by turning on another scan option. Some level of distortion is the cost of doing business with EPI. Sorry.

Signal dropout

As already mentioned, strictly speaking this artifact isn't a pure characteristic of EPI because similar regions of dropout will occur for any gradient echo imaging sequence acquired at the same TE. But, in fMRI we have a requirement to generate BOLD contrast by setting the TE within a certain range of values; the optimum BOLD contrast occurs when the TE matches the local T2* of the tissue of interest. Thus, we don't typically minimize the TE in order to minimize signal dropout, we instead set our TE to provide the BOLD contrast we want and we have to accept the concomitant signal dropout. (See Note 8.)

Let's first consider the origins of the effect, then we can consider tactics to reduce it. As with the other artifact sources, dropout arises because of spurious dephasing caused by magnetic susceptibility gradients. As before, the frontal and temporal lobes are particularly at risk because of the nearby presence of sinuses and air-filled bony cavities (ears). Our EPI sequence is T2*-weighted (by design), meaning that the signal contrast across an image will be due primarily to T2* variations across the brain, with T1 and spin density differences contributing secondary contrast. If we set the TE to be, say, 40 ms in order to match the approximate gray matter T2* in occipital cortex, in the frontal cortex where the T2* might be as low as 20 ms there is going to be pronounced signal attenuation (not to mention sub-optimal BOLD contrast).

Here are echo planar images for two different TEs acquired at 3 T. Notice how it's the frontal and temporal lobes, as well as the inferior surface, that pay the biggest penalty in terms of signal loss (dropout) at the longer TE, although if you look very closely you can also see that signal is reduced for the edges of the brain in all slices compared to signal at the shorter TE:

(Click to enlarge.)

What can be done to ameliorate the problem? One tactic is to compromise between the optimum TE for BOLD in regions of the brain that are well shimmed (i.e. that have low magnetic susceptibility gradients), such as occipital and parietal cortices, and a shorter TE that will retain frontal and temporal signals. If you're only interested in occipital cortex at 3 T you might use a TE of 40-50 ms. But if you're interested in other regions, especially inferior frontal and lateral temporal regions, you might opt for a TE of 20-30 ms; the loss of functional contrast in occipital/parietal cortex should be small.

Another simple tactic is to recognize that the majority of dropout is due to phase variation through the imaging plane, so lots of thinner slices will generally provide less dropout than fewer thicker slices. This is illustrated in the following cartoon, where for simplicity I will pretend that the entire slice has a uniform T2*:

In the thick slice (top) the phase variations across the upper and lower portions of the slice tend to cancel, resulting in a small net signal vector, in blue. When the thick slice is divided into two thinner slices (bottom) the partial cancellation of spin phases is reduced, leading to larger net signal vectors for the two thin slices combined. (I've exaggerated the lower blue vectors for emphasis.) Note also that because we nearly always deal with magnitude images - no phase information - the direction of the net vectors is irrelevant; only their magnitude determines image signal level.

Of course, acquiring arbitrarily thin imaging slices reaches a practical limit very quickly, not least because coverage in the slice dimension is decreased. Remember that if it takes (a typical) 60 ms to acquire each EPI slice it will take about the same time whether that slice is 2 mm thick or 4 mm thick! If you need whole brain coverage and you don't want to violate the Nyquist condition for the hemodynamic response, you only have a TR < 3000 ms to utilize. But even if you only need a handful of slices you may run into another limit. Remember that the overall image SNR scales with the volume of the signal being sampled. While there may be some regional benefits (reduced dropout) from acquiring two 2 mm slices instead of one 4 mm slice, you should note that the individual image SNR will be reduced by 50% on purely volumetric grounds. The regional benefit of lower dropout might come at the expense of worse global statistical power. It all depends on your application.

Final thoughts

Variability of artifacts

There is a fair degree of anisotropy in the way EPI artifacts behave, their relative severity, and so on. As with nearly all aspects of fMRI acquisitions, it means that you have (limited) choices that may be able to produce better data for your application. These include the slice direction, the quality of the shim, the slice thickness, the sign of the imaging gradients, TE, echo spacing and so on. What you actually observe in any situation will result from the interplay of all of these factors, implying that no single parameter should ever be modified in isolation; you need to consider concomitant effects. We will look at some of these issues as we go through real data and real artifacts.

Alternative EPI pulse sequences

There are several EPI variants in common use that don't acquire an entire 2D k-space plane in a single shot. These include partial Fourier EPI, segmented (multi-shot) EPI and accelerated EPI using parallel imaging methods such as GRAPPA. The methods can be differentiated based on the way they acquire the phase encoding dimension of k-space, and you should find that the background physics you've learned in this series of posts will enable you to understand each one. That's why I'm stopping at this point; you have the background knowledge that you need to comprehend all the EPI variants as and when we encounter them in practical situations.

________________________________

Notes:

1. Some recent work at NIH suggests that we should be using lower flip angles for fMRI than many of us are, and certainly lower than the Ernst angle (which produces maximum SNR per unit time, i.e. for the TR being used). This is because in fMRI we are usually operating in a regime limited by physiologic noise in our time series acquisitions, we generally don't have to worry nearly as much about thermal noise (unless the scanner has developed a problem). BOLD contrast tends to be invariant to the particular flip angle being used, whereas it is highly dependent on TE. Physiologic noise, on the other hand, tends to scale with raw image SNR, i.e. with flip angle, and it also has a TE dependence. Overall, the NIH group found that one could operate in an optimum regime by reducing the physiologic noise in the time series, even though the raw image SNR was low. This all makes sense to me. I haven't quite got around to recommending that everyone use single digit flip angles, as tested in the paper, but I have started suggesting that people scale down the RF to between 30-50 degrees while I investigate further.

"Physiological noise effects on the flip angle selection in BOLD fMRI."
J Gonzalez-Castillo, V Roopchansingh, PA Bandettini and J Bodurka.
Neuroimage 14; 54(4):2764-78 (2011).

On a related note, I was talking to a member of Peter Bandettini's group at a conference earlier this year, and I hear that this reduced flip angle suggestion holds for resting state fMRI as well as task-based experiments as tested in the above reference. This also makes sense to me.

2. While there are practical consequences arising from the sign of the so-called "dephasing" gradients (colored orange), for the purposes of understanding the EPI sequence and its k-space trajectory the signs of the orange and green gradients don't matter at all, provided they are balanced in the manner I've drawn them. I just happened to start with a negative dephasing gradient for Gx and a positive dephasing gradient in Gy, but I could just as easily have reversed one or both of them. In later posts we will look at the consequences of the read and phase encode gradient signs because our spatial encoding gradients tend to interfere with the intrinsic and anisotropic (and unwanted) gradients arising in your subject's head. So, while there are often practical consequences, at this stage of the game it's not important whether we start with negative or positive gradient lobes.

3. Remember convolution from Part Six? Recall how an exponential decay function applied to a time domain signal causes a broadening of the frequency domain representation? Well, T2* decay during the echo train is going to cause a certain amount of smoothing to the image. If you simply must know more about this subject now, take a look at pages 262-5 of Introduction to Functional Magnetic Resonance Imaging (2nd Edition) by Rick Buxton. I've not yet managed to find a good (free) online description of T2* broadening I'm afraid, but I'll keep searching.

4. There are several practical sources of erroneous phase shifts that can cause zigzag patterns in EPI k-space. Some of them are discussed in my user training doc/FAQ available from an administrative post in April. I'm not going to get into the practical causes in this post, except to mention examples in order to illustrate the consequences for the resulting images. Instead, we will consider each artifact source independently as we look at real data in subsequent posts because it's important to get a sense of relative importance of each source, and whether or not you can do anything to improve the situation at the time of the experiment.

5. These zigzags in k-space may be localized or global, depending on the origin of the particular offset. In other words, for many parts of the brain the acquired k-space trajectory may match quite well the ideal (target) trajectory, whereas in other parts (e.g. frontal lobes, temporal lobes) there may be appreciable mismatch. Remember that there is no single k-space trajectory for the entire brain, only a global target (ideal) trajectory. Localized variations in k-space (i.e. zigzags) may be due to magnetic field heterogeneity (i.e. imperfect "shim"), but global variations may result from a mismatch between the ADC sampling periods and the readout gradient waveforms (e.g. because of the time it takes to ramp up the readout gradients). The important point to remember, as you read about the three typical artifacts in EPI, is that there will likely be considerable spatial heterogeneity in the severity of artifacts across the brain. Not all regions of the brain are created equally when it comes to EPI!

6. The obvious tactic of simply increasing the FOV until the ghosts fall into noise regions, thus not overlapping valuable signals, can be employed but tends to penalize severely spatial resolution in the phase encoding dimension. In most fMRI studies, therefore, the N/2 ghosts are allowed to overlap signal regions and we do our best to minimize ghost intensity.

7. There are numerous published schemes for correcting N/2 ghosts, each with pros and cons. But you probably have limited scope to change the correction scheme delivered with your scanner software. Each method typically requires a specific pulse sequence as well as online software to generate and apply the correction terms to the k-space. Put another way, establishing new methods for ghost correction is likely going to be a lot of work for someone at your facility!

You might also be interested to learn that there is one version of EPI that is entirely ghost-free. It's usually referred to as "fly-back" EPI on account of using all positive (or all negative) readout gradient periods:

With fly-back EPI the k-space rows are all read with the same polarity, in the same direction, eliminating the need to reverse alternate rows and also eliminating the zigzag phase errors across the phase encoding dimension. Thus, no ghosts.

Wondering why we don't all default to using fly-back EPI for fMRI, negating the ghost problem entirely? It's because the demands on the gradients are far higher for fly-back EPI than for regular EPI. You can only slew (switch) the gradients so fast before you start causing peripheral nerve stimulation in your subject. Furthermore, the time between readout periods is longer for fly-back EPI than for regular EPI, making the distortion problem worse in exchange for eliminating the ghosts. Doh!

8. There are customized EPI pulse sequences that can reduce dropout, such as so-called "z-shimming" methods, but these usually come at the expense of lengthening the pulse sequence (through the incorporation of the compensation scheme) and thereby decreasing the number of slices that can be acquired for a given TR. Few of these custom sequences are offered on commercial scanners anyway. I may deal with z-shimming methods in a future post, but because they're not in widespread use I won't make it a high priority.

Light relief (to buy me time).... This year's IgNobel in Medicine

2011-10-11T11:08:00.000-07:00

Anyone who has ever experienced an fMRI scan knows two things about the effects of the method on a subject: (1) it's soporific, and (2) like a long car journey, you don't need to pee until five minutes after you've started. So this year's IgNobel Prize in Medicine, awarded jointly to two groups, caught my attention. Their work shows how the need to urinate can affect performance on some simple mental tests - just the sort of tests that we use in our fMRI experiments.

Implications for fMRI?

An enjoyable summary of the winning researchers' work is available on this Scientific American blog. According to this summary (Yeah, I haven't got around to reading the papers themselves yet. I'm training to be a mainstream science journalist ;-), needing to pee could have your subjects performing better (yes, better) on delayed gratification tasks, but worse on cognitive tasks. I take these results at face value - I have to, I've not read the papers - but I do want to think a little more about the implications for fMRI studies. It's hard enough keeping people awake, let alone motivated to do a task. And as for providing *additional* motivation for a task... The mind boggles!

"I feel the need, the need to pee!"

So, short of rejecting subjects who rush to the toilet the moment they get out of the scanner, what else could we do to control for the effects? Perhaps we could insist that subjects must be able to sit in a waiting room for 20 minutes post-scan - no pee - and only then opt to retain their scan data.

What else might produce a similar effects in subjects? General discomfort? You have to wonder, given the "need-to-pee" effect, whether a subject's general state of (un)happiness in the scanner might well be interfering with his mental performance. If so, having pressure points in a subject's lower back or whatever could have him significantly altering his task ability.

Alternatively, perhaps the need to pee and general discomfort merely increases a subject's propensity to move. We all know that this is one of the Stages of Having to Pee:

So perhaps this amusing research has some important ramifications for fMRI studies after all. As with so many other state factors - caffeine use, stress, menstrual cycle, etc. - it could just be another in a long litany of issues that contribute to our relatively poor inter-subject variability. You know my feeling on the matter: if you can control for it, control it. And if you can't control for it but you can measure it, MEASURE IT! Would it really be the end of the world if you were to ask your subject to rate her "need to pee" as she exits the scanner? How amusing would it be to see your effect disappear having regressed out the "need to pee" score?

PS I really will have a last post on EPI k-space along very soon! I promise!

More on decoding - Jack Gallant radio interview

2011-09-27T18:11:00.000-07:00

KQED radio had a half hour segment with Jack Gallant this morning, discussing the study published by Shinji Nishimoto et al. last week.

Here's a link to the KQED archive. An MP3 is also available.

"Reconstructing visual experiences from brain activity evoked by natural movies."

2011-09-22T23:02:00.000-07:00

(Waiting for the next post in the physics series? Apologies. New academic year = teaching. Next post, on EPI, will be along within a few days...)

Gallant Lab strikes again!

Another shameless plug for fMRI at Berkeley! A study with the title of this post was published online today in Current Biology. Seems like it's generating as much media buzz as their previous study (Nature, 2008) using still images. Congratulations to Shinji, Joseph (An), Thomas and co. on another fine study. (And people said the Varian didn't work. Ha!)

Some links for more information:

UC Berkeley News Center.

Gallant Lab website. Best place to start, see the FAQ.

YouTube videos of some of the results:

Physics for understanding fMRI artifacts: Part Eleven

2011-08-15T09:26:00.000-07:00

Resolution and the field-of-view as seen in k-space

Understanding how distances in k-space manifest as distances in image space is quite straightforward. All you really need to remember is that the relationships are reciprocal. The discrete steps in k-space define the image field-of-view (FOV), whereas the maximum extents of k-space define the image resolution. In other words, small in k-space determines big in image space, and vice versa. In this post we will look first at the implications of the reciprocal relationship as it affects image appearance. Then we'll look at the simple mathematical relationships between lengths in k-space and their reciprocal lengths in image space.

Spatial frequencies in k-space: what lives where?

I mentioned in the previous post that there's no direct correspondence between any single point in k-space and any single point in real space. Instead, in k-space the spatial properties of the object are "turned inside out and sorted according to type" (kinda) in a symmetric and predictable fashion that leads to some intuitive relationships between particular regions of k-space and certain features of the image.

Here is what happens if you have just the inner (left column) or just the outer (right column) portions of k-space, compared to the full k-space matrix arising from 2D FT of a digital photograph (central column):

An illustration of the effect of nulling different regions of k-space from a full k-space matrix, applied to a digital picture of a Hawker Hurricane aircraft. The full k-space matrix and corresponding image are shown in the central column.

Inner k-space only:

The inner portion of k-space (top-left) possesses most of the signal but little detail, leading to a bright but blurry image (bottom-left). (See Note 1.) Most features remain readily apparent in the blurry image, however, because most contrast is preserved; image contrast is due primarily to signal intensity differences, not edges. If this weren't true we would always go for the highest signal-to-noise MRIs we could get, when in practice what we want is the highest contrast-to-noise images we can get! Imagine an MRI that had a million-to-one SNR but no contrast. How would you tell where the gray matter ends and the white matter begins? Without contrast no amount of signal or spatial resolution would help. So much for SNR alone!

Outer k-space only:

If we instead remove the central portion of k-space (top-right) then we remove most of the signal and the signal-based contrast to leave only the fine detail of the image (bottom-right). Strangely, though, it's still possible for us to make out the main image features because our brains are able to interpret entire objects from just edges. In actuality, however, there is very little contrast between the dark fuselage of the Hurricane, the dark shadow underneath it and the dark sky. Our brain infers contrast because we know what we should be seeing! If we were to try doing fMRI, say, on a series of edges-only images we would run into difficulties because we process the time series pixelwise. With a relatively low and homogeneous signal level you can bet good money the statistics would be grim.

Whole k-space:

The central portion of k-space is important because it provides the bulk of the image signal as well as the signal-based contrast, while the outer portions of k-space provide image detail, in particular establishing the boundaries in image contrast. Having only one or the other might not prevent us, by inspection, from being able to recognize an object in an image, but it may not suffice for pixelwise processing. The objective in fMRI isn't simply to be able to recognize an image as that of a brain!

So why bother to categorize k-space in this manner? Well, for starters, in several future posts we will need to consider the effective k-space matrix to understand many properties of an EPI time series as used for fMRI. When we look at spatial smoothing, for example, it will be imperative for you to understand where in k-space the primary effects of a smoothing function are manifest. A second reason concerns artifact recognition. This simple, intuitive partitioning of k-space regions can be extremely useful when it comes to diagnosing certain data artifacts. Because of the reciprocal relationship, a feature that is widespread (spatially) in image space will likely be focal in k-space. Tracking down a focal artifact source can be considerably easier to do. But I digress. We will look at artifact recognition in the next series of posts. For today I am going to focus on clean data and restrict the topic to features in an ideal image.

Why does the signal level change across k-space?

Here, I am going to offer you two alternative explanations. First, an MRI explanation considering the action of the imaging gradients. We know that whenever a phase is imparted across a sample there will be some partial signal cancellation (as we saw in Part Eight). Near to the center of k-space the signal is high because the amount of phase applied by the imaging gradients to the sample magnetization is low; the degree of signal cancellation is low. The more spatial information (detail) we try to encode, the more phase we need to impart to the signal, the farther the signal level will be reduced. In the outer regions of k-space, where the imaging gradients are comparatively large and the concomitant dephasing is relatively large, the signal level will be diminished.

But the pictures I've presented aren't MRIs, they're digital photographs. Thus, an alternative (and technically more correct) explanation is to consider the spatial frequency content of the image. The image of the Hurricane contains more broad areas of relatively uniform intensity - clouds in the sky, the grass, large blobs of camouflage painted on the wings and fuselage, etc. - than it does edges and other fine details. And since we now know that edges live in peripheral (high) k-space regions whereas spatially broad features live towards the center of k-space, we can consider the k-space plot as a kind of "spatial content map." There is simply more image content to map that changes slowly with distance than there is content that changes rapidly with distance (i.e. detail).

In physical terms, going into high k-space regions means we are encoding high spatial frequencies. And an edge is something with a high spatial frequency; the feature changes rapidly over a short distance. At this point we can even make a prediction. It's reasonable to predict that to get more resolution - the ability to resolve finer structure - we will have to push out to higher k values. We'll deal with this last point in the image resolution section, below.

Defining parameters in k-space to yield the image you want

Okay, now you have a rough idea of what features live where in k-space, it's time to return to entire k-space matrices and learn how to establish k-space to yield an image having the spatial properties that you want. Here, in illustrative form, are the spatial parameters we need to consider:

Delta-kx and delta-ky are the steps in each k-space dimension, while 2kxmax and 2kymax are the spans of k-space. In this example the k-space steps and spans are equal so the resulting image is a uniformly sampled square, but that doesn't have to be the case. FOVx and FOVy define the image size while delta-x and delta-y define the pixel size, i.e. the spatial resolution. (See Note 2.)

Image field-of-view

The relationship between k-space and the image FOV is straightforward. The reciprocal of the k-space step, delta-k, defines the image space extent, the FOV, i.e.

Here is a k-space matrix with small delta-k and its corresponding image:

If the span of k-space (i.e. the maximum k value) is left constant but the step size is changed, the effect on the image is to alter the FOV. If delta-k is increased the result is a reduced FOV, i.e. a zoomed image. Here is the same image with delta-k doubled in the y direction only, resulting in an image that is unchanged in x but that has half the FOV in y (see Note 3):

You might be wondering why this inverse relationship holds. What is it about the delta-k value that sets (or restricts) the image to a particular size? The relationship arises because of a restricted ability to interpret the phase changes imparted across the magnetization in the sample by the imaging gradients. Between gradient increments - that is, between successive sampling points under the frequency encoding gradient for the x dimension, or for each phase encoding step in y - we can't impart more than 360 degrees of phase because we cannot discriminate between 360 degrees and 0 degrees. Phase changes greater than 360 degrees "alias," so that a change of 450 degrees would be measured as only 90 degrees, and so on. It's the Nyquist sampling theorem, which we saw in Part Six, in another guise. (The algebra demonstrating the 360 degree phase discrimination limit is in Note 4.) All you need to remember is the simple inverse relationships given in the blue box, above.

Since we know that k is the time integral of the (readout or phase encode) gradient being applied to encode spatial information, this inverse relationship produces an interesting observation: big images are easy to produce, smaller images are more difficult to produce. Obtaining a small delta-k requires just a low amplitude gradient or a small amount of time under a gradient, which is obviously easier to achieve experimentally than either a large amplitude gradient or a protracted gradient period.

At first glance this seems a little counter-intuitive, but that's because it's not the whole story. The image is comprised of a fixed number of pixels (arising from the same number of k-space samples), so getting large images is not the freebie it might appear at first blush. If your image is, say, 64x64 pixels total then an image a meter on a side isn't going to be very useful for brain imaging! We need to know the resolution of the image - the size of the pixels - before we determine whether the k-space matrix is in fact appropriate.

Image resolution

You saw above how restricting the k-space coverage to a small, central region from a larger matrix has the effect of blurring the image. So it should come as no surprise that simply extending the size of the k-space matrix will have the effect of increasing the image resolution.

Defining resolution is straightforward now that we have already got the FOV relationships established. All we need do is divide the FOV in each dimension by the number of pixels defining that dimension to see the k-space relationship:

Remember that for an Nx x Ny image we acquire Nx/2 and Ny/2 values of k-space either side of zero, making the total span of k-space equal to 2kxmax and 2kymax.

Here's an example of a k-space matrix having a large extent, yielding a high-resolution image:

If delta-k is maintained, to keep the FOV constant, but the extent of k-space values is restricted then the image resolution decreases:

Here the image is blurred because the unsampled white area around the central (sampled) square of k-space has been "zero-filled." (See Note 5.) This produces a smoothing effect in the image. If the zero filling were not performed then the image would appear pixellated instead of smooth. Either way, the actual resolution of the image is reduced from the previous, high maximum-k situation.

What does this reciprocal relationship between total k-space extent and image resolution mean experimentally? Getting more resolution in the image requires larger k values, requiring either larger amplitude or longer gradient episodes (or some combination of the two). Thus, we can now see that while it is easy to get a large image FOV, it is difficult to get high image resolution! We will have to drive the gradients harder (larger amplitude) or leave them on for longer to attain smaller pixels. Indeed, this is probably the biggest single limit to MRI performance. Gradients can't be made arbitrarily high amplitude for engineering and safety reasons; large, rapidly switched gradients tend to cause peripheral nerve stimulation in the subject, as well as unwanted residual effects (eddy currents) in the magnet hardware, for example. And gradients can't be enabled for arbitrarily long durations or the signal that we're using to encode spatial information is likely to have died away to (near) zero, making our desired high-resolution image very low signal indeed.

In physical terms, the reason why k must be pushed very high to get high spatial resolution in the image is due to the ability of the imaging gradient to impart a significant phase difference between two nearby spatial positions. For a position, x and a nearby position x' we must be able to distinguish the phase difference in order to resolve these two positions as unique, and not the sum of the two. Only as the imaging gradient's time integral gets very large (leading to high k values) does there start to manifest a measurable phase difference between each increment in k-space. Sadly, as far as we know today, there is no way around this resolving power limitation. It's not the fault of the pulse sequences per se but a fundamental limitation in the way we encode spatial information with magnetic field gradients.

Okay, that's enough of the general properties of k-space for the time being. In the next post we will return to k-space trajectories and pulse sequences. We're finally ready to see the workhorse of the majority of fMRI experiments: the echo planar imaging (EPI) sequence.

---------------

Notes:

1. The overall signal level and the image contrast are predominantly established by other aspects of the pulse sequence, such as the excitation flip angle and repetition time, not the k-space coverage scheme. So here all we're considering is where that signal, established prior to spatial encoding, ends up residing in k-space.

2. Some of you may already be aware that the real shape of pixels in MRI isn't actually rectangular. They are in fact defined by a "point spread function." Without any sort of filtering the PSF is sinc-like because the sampling window is a square; we saw in Part Six that the FT of a square function is a sinc. In EPI we also have some smoothing in at least one dimension arising from T2* relaxation. I don't want to get sidetracked with these issues at this point, instead we will deal with the true shape of pixels when we consider the Gibbs artifact (or ringing) in a future post, because these are two sides of the same coin.

3. There's one more FOV consideration: is the image FOV big enough? In Part Six we saw the effects of aliasing when there were insufficient data points to properly sample a waveform. This was encapsulated by the Nyquist theorem. And here's where that data sampling restriction enters into imaging. We will look in detail at aliasing and the image FOV in a later post, early in the series on artifact recognition.

4. In order to see why the FOV should be inversely proportional to delta-k we need to consider the phase evolution between successive k-space points. For frequency encoding along x we have:

For phase encoding along y we have:

There are two observations we can make from the two relationships just derived. Firstly, these relationships reinforce the fact that MRI axes are frequency axes - delta-wx and delta-wy - and not really spatial axes at all. However, the spatial labels can be made appropriate after a little bit of algebra. (That's rather the point with this whole k-space formalism. We use it because it's a more intuitive, convenient way to think about imaging than the time/frequency relationships.)

Secondly, we can now see that the maximum phase shift imparted by the imaging gradients to the extreme spatial positions along each image dimension (i.e. to the edges of FOVx and FOVy, which are defined as the equivalent frequency ranges of delta-wx and delta-wy, respectively) is exactly 2pi radians, or 360 degrees, during each k-space increment. The result is identical for the frequency encoding axis (x) and the phase encoding axis (y). This isn't a coincidence.

Noting that the FOV is really a frequency range with a central "carrier frequency," the 360 degree phase range is really +/- 180 degrees relative to the carrier frequency's phase (which is the nominal zero position). At spatial positions less than the FOV (i.e. within the image, but not at or beyond the limits of the FOV) the phase imparted by each k-space increment is somewhere in this +/- 180 degree range. Outside of the image FOV is where things get interesting. Remember aliasing from Part Six? We cannot meaningfully encode any more phase than 360 degrees because 360 degrees and 0 degrees are indistinguishable; phase is modulo(360).

We will consider FOV and aliasing again in one of the first posts on artifact recognition, because aliasing is one of the issues that can easily affect fMRI data if the scanner operator isn't careful. (Yeah, you can't blame aliasing on your subject! It's pure pilot error!)

5. I don't have plans for a post on image smoothing or zero-filling at the moment, but the issues will be covered in part when I do a post on the Gibbs ringing artifact in the artifact recognition series. For fMRI most intentional image smoothing is done offline, in the image domain, as part of a processing pipeline rather than on the scanner (before 2D FT), so it's not really a scanner acquisition issue.

Physics for understanding fMRI artifacts: Part Ten

2011-08-06T01:40:00.000-07:00

(For the answer to the homework k-space diagram given at the end of Part Nine, see Note 1.)

K-space in two dimensions

As anyone knows who has encountered MRI professionally, whether in research or medicine, there seems to be an endless array of pulse sequences to choose between. The variety can be overwhelming at first. Nor is the situation helped by different vendors using different acronyms - we always use acronyms in MRI! - for what are essentially the same sequence.

It's little wonder, then, that most neophytes' eyes glaze over when it comes to comparing and contrasting any two pulse sequences if the taxonomy appears to be ad hoc. Where on earth to start? But it turns out that most pulse sequences can be categorized fairly easily, and their heritage traced, by separating the part(s) of the sequence that is responsible for spatial encoding, from the part(s) of the sequence that will provide the tissue or functional contrast. Occasionally there is overlap within the sequence of these two missions, but even then it's usually straightforward to understand the spatial encoding and interpret its genesis.

A useful pictorial representation of imaging pulse sequences

It turns out that there are only a handful of spatial encoding methods in common use these days, almost all with roots in the late 1970s or early 1980s. While new pulse sequences appear in the literature all the time, when you look at their k-space representations you'll be able to see how each new method has developed from a small number of key ideas from those early years. It's possible to categorize the encoding methods without k-space, but the k-space formalism makes comparisons trivial (in MR terms).

Spatial encoding methods can be separated into families derived from a central idea. For instance, following Lauterbur's original imaging paper in 1973 (which led to the family of projection reconstruction methods), in 1975 Richard Ernst's group came up with a sequence that utilized a 2D Fourier transform to yield the final image. (See Note 2.) It was a remarkable breakthrough and is the grandparent of nearly all medical/biological sequences still in common use today.

Still, even geniuses miss opportunities every now and then. And in 1980 a group at Aberdeen came up with a far more practical implementation of Fourier imaging, using amplitude-modulated gradients in a "constant time" pulse sequence, rather than the fixed amplitude, variable time scheme of Kumar, Welti and Ernst. It is this constant time scheme, which the Aberdeen group termed "spin warp" phase encoding, that provides the basis for most clinical (anatomical) scanning used today. It's also a good scheme to look at when first encountering 2D k-space, so we'll consider it in detail in this post.

The goal revisited

In the first part of the last post (see Part Nine) I used two examples of digital images to illustrate how the information content in a 2D plane of image pixels can be equivalently represented in reciprocal 2D space, or k-space. I mentioned that both the images and the k-space comprised 512x512 points, but later on when I started to draw (one-dimensional) k-space trajectories I did so on a k-space plane that was represented by just a set of axes, not discrete points. In case you think that image space and k-space in MRI are continuous, I'm going to spend a moment considering the digital k-space plane explicitly. (Like real space, k-space can also be continuous rather than digital, but that's not how MRI works.)

Here is a 16x16 plane of k-space points (see Note 3) overlaid on some actual signals to reinforce the point that we're digitizing a continuous process:

Courtesy: Karla Miller, FMRIB, University of Oxford.

The goal is to traverse the entire k-space plane, i.e. to use our gradients to follow a trajectory that crosses every single point (as defined by the white grid itself), acquiring data (with our receiver coil), one point for each grid coordinate, as we go. Once we have traversed the entire 2D plane (and assuming a suitable data acquisition scheme) we will have 16x16 k-space data points and will then be in a position to apply a 2D FT and get a 16x16 image out. (See Note 4.)

It is worth emphasizing here an important point that is often overlooked when people are thinking about k-space. We aren't actually doing a 2D FT to achieve the k-space representation and the pictorial analysis of the gradient actions. Rather, we are simply recognizing that a 2D FT of the image plane is its 2D k-space representation; hence, the action of the imaging gradients is to trace through each point in that 2D k-space. Semantics? I don't think so. The 2D FT is ultimately required to recover the actual image, but it is performed on the completed 2D k-space plane, i.e. only after the k-space plane has been properly sampled by the action of the imaging gradients. We don't need to do the 2D FT in order to understand how the pulse sequence is encoding spatial information! And that's what makes the k-space picture so valuable as a sequence comparison tool. All (2D) imaging sequences, from EPI to RARE/FSE, must achieve the same completed plane of k-space before the final image can be recovered (by 2D FT). Even spiral isn't immune to this requirement! (One of the processing steps for spiral is to "re-grid" the k-space trajectory so that it is rectilinear and can then be fed into a regular 2D FT algorithm.)

With this appreciation of the intuitive meaning of k-space under your belt, it's time to see the action of the imaging gradients as they trace through the entire 2D k-space plane.

Gradients along the x direction (again)

Last time out we didn't consider the data acquisition at all, but in this post it's going to be reintroduced. We know we need to acquire a signal corresponding to each point on the k-space grid and an easy way to achieve that is to acquire one line of kx information at a time, then figure out how to move across ky to hit each row of kx points in turn. The sampling process is therefore essentially the same as we saw previously for the gradient echo in Part Eight:

Except that now we can see that the period of signal acquisition (analog-to-digital conversion) occurs coincident with a journey from maximum -kx to maximum +kx, as we saw in Part Nine:

Eight data points are sampled during period 2, then the central value at kx=0, then a further seven data points during period 3. (Yeah, I know, the green arrow should have stopped one square earlier than I drew it. Sorry!) In the next post we will look at the effects of the number of k values, the space between them (delta-k) and the extreme k values attained, because these parameters determine the image resolution and field-of-view (FOV). For now don't worry about them, let's just get an image to be going on with.

Gradients along the y axis

Clearly, to get off the ky=0 axis we're going to have to apply a gradient, Gy (just as was implied by the homework problem at the end of Part Nine). What's more, if we are interested in sampling a rectilinear grid then there is no point having Gy turned on when readout under Gx is happening, otherwise we will end up with a diagonal trace in k-space (as per the homework example last post) and we will spend a lot of time "missing" the grid points. (I won't discuss non-rectangular sampling in this series of posts. Perhaps I'll do a separate post on spiral scanning and its ilk at some point in the future.)

However, the -Gx period that precedes data acquisition is of no value to the data either, it just gets us to one side of kx space so that we can zip along one entire line of kx points in one go. Why not put the Gy gradient coincident with that? All that matters is that we have moved in ky prior to the readout line along kx. Handy! In this situation the actual vector in k-space - diagonal or otherwise - is of no consequence. We simply want to have arrived at the target kx,y coordinate as quickly as possible prior to the start of data acquisition under +Gx.

This can be achieved with the following pulse sequence (see Note 5):

Here's the corresponding k-space trajectory:

Now let's repeat the previous pulse sequence, starting from scratch, but this time we will reduce the amplitude of Gy by 1/8th. Look closely for the slight reduction of Gy; the previous value is indicated by a dashed line:

This time through we only get 7/8ths as far in the +ky direction before the data sampling commences along the kx points:

If we keep on repeating this acquisition process, stepping down each time by another 1/8th increment reduction of ky, then after sixteen total experiments we will have traversed every line of ky as well as every kx point along the sixteen ky rows. Here are the sixteen k-space trajectories on a single diagram (see Note 6):

Let's recap for a moment. Even though the data sampling only happens as the rows of k-space are being traced out, i.e. under the x gradient, it's clear that in stepping down a row each time through the pulse sequence we also hit every column and thus completely sample the 2D plane. So, although the green traces above all appear to be identical, in actuality each row of kx values is slightly different. The incremented Gy gradient that precedes the Gx readout gradient (and data acquisition) has imparted a phase increment to the sampled magnetization. I'll come back to this point below.

One dimension's just like the other one

The mathematics of the second, phase-encoded k-space dimension is quite straightforward. Indeed, to have followed the trajectories above you must have mentally integrated the changing area under Gy, first by recognizing that the initial Gy gradient had the same amplitude and time as the +Gx episode in the sequence. Thus you deduced (whether you realized it or not) that the areas under Gy and Gx were initially the same, making the first trajectory a diagonal journey from k0,0 to the top-left corner of the k-space matrix. You did the math without even realizing it! (That's the beauty of the k-space formalism. No math required!)

But for those of you who would like to see the equivalence of the two k dimensions, the definition of ky follows naturally from what we saw in Part Nine for kx. Recognizing that all the y gradient is doing physically is adding another phase shift to the signals, we can simply add one more phase term to the time-dependent signal equation and then recast the signal in terms of space and reciprocal space as the conjugate variables. Adding terms for y and ky to what we had in Part Nine for x and kx, we get:

While the definition of ky appears different to kx it's actually just because in this particular pulse sequence the y gradient amplitude was stepped in equal increments. If we define the maximum value used as Gy then we can write Gy = (N.gy)/2, where there are N total steps (we had N=16) and the increment is gy, and the two k variables are seen to correspond. In both cases the definition of k is simply the area under the gradient (which is G.t for a square gradient) multiplied by some constants (gamma/2pi) that we can ignore.

Note that whereas the kx value evolves with the x gradient's duration, for the ky (or phase-encoded) dimension the gradient time, tp is constant and the amplitude is changed. In both cases, however, the areas under the gradients are changing in equal amounts between each point in the k-space matrix, making delta-kx and delta-ky equal. (See Note 7.)

With our k-space plane fully sampled we can now do a 2D FT to get an image. Recalling that a 2D FT is actually two 1D FTs performed in succession, we can see that Fourier transformation along ky gives us the y dimension of an image, I(y):

Gaining an intuitive understanding of phase encoding

That's it! You've covered a 2D k-space plane completely and have arrived at the point where a 2D FT will yield an image. Job done, right? However, some of you - quite possibly those of you who have attended a class on MRI and have been introduced to phase encoding from a different perspective - might still be wondering what it means to do phase encoding, and might be struggling to understand how it is different from frequency encoding. How does phase encoding work?

Often it is much easier to see how frequency encoding works - as a projection of the spins in one dimension - because it seems a bit like a shadow cast from an object by a bright light. Well, phase encoding is actually the same! If you took the kx=0 column from the 2D k-space plane and FT'd it you would obtain a 1D projection of the object along the y axis. Sure, it took a lot of steps - sixteen separate acquisitions, one ky value on the kx=0 line per acquisition - to get to the point where you could do this, but the 1D profile of the kx=0 line from our 2D plot would look very similar indeed (discounting experimental issues) to a separate experiment where one used a single frequency encoding gradient to get a 1D profile along y.

How is it that these processes - frequency and phase encoding - are seen to be equivalent in the final analysis? Because, as we saw in Part Four, frequency is just the rate of change of phase! In this case, what we're doing is imparting a stepped phase into a succession of waveforms, so that when we project in a dimension orthogonal to those waveforms (by "joining the dots," or data points) we are actually constructing new waveforms in the orthogonal dimension. The only true differences between the frequency-encoded dimension and the phase-encoded dimension are practical ones, as we will see in the artifact recognition posts to come.

Final thoughts to help you understand phase encoding, and to reinforce the mathematical relationship between frequency and phase:

A frequency encoding gradient imparts a phase that changes quickly (and continuously) with time, because the frequency encoding gradient is turned on and left on.

Conversely, the phase encoding gradient is used discretely. A stepped phase shift is imparted to the magnetization and then the phase encoding process is suspended, leaving a "phase memory" in the magnetization for when it is ultimately read out under the frequency encoding gradient.

In essence, then, frequency encoding could be considered "fast, continuous phase encoding." Phase encoding is slow and discrete. Want an analogy? Frequency encoding is like driving uninterrupted down the freeway from SF to LA without traffic to impede your progress. Phase encoding is like being in a stop-start traffic jam all the way. The total distance is the same, of course, but in the case of stop-start traffic the journey takes a lot longer.

Next post, a quick review of image field-of-view and resolution as they appear in k-space. Until then I would strongly urge you to consider the Further Reading, below.

-------------------

Further Reading:

Introduction to Functional Magnetic Resonance Imaging: Principles and Techniques (2nd edition) by R.B. Buxton.

For a description of how phase encoding works, read the section entitled "Phase encoding" on pages 214-6. Rick shows how frequency and phase encoding are equivalent in gradient terms and explains, with the assistance of some simple figures, how one can consider little chunks of gradient history to understand the similarities of the two processes. If you understand this section you will have developed an intuitive feel for gradient spatial encoding in MRI! Then keep reading for a nice description of k-space on pages 216-9.

Functional Magnetic Resonance Imaging (2nd edition) by S.A. Huettel, A.W. Song & G. McCarthy.

Another description of 2D spatial encoding (frequency and phase encoding) using the k-space formalism appears on pages 109-117. You will also be introduced to spatial frequencies in k-space, should this whole consideration of k-space have left you slightly numb and wondering what reciprocal space is all about.

Notes:

1. One possible answer to the k-space trajectory problem at the end of Part Nine:

It's a junk sequence purely to illustrate the convenience of the k-space representation, it doesn't have any specific role or name.

Note that in many k-space diagrams there is no information on the speed through any segment. In my quiz diagram I implied that each arrow represented a similar time for each gradient episode, but that doesn't have to be the case. Thus, another correct pulse sequence diagram could have used these relative magnitudes to keep the time integrals of the gradients consistent with the k-space representation:

In this case the speed through the diagonal trajectory would be half the prior version. Does it make a difference? Experimentally, yes. (There would be more relaxation effects in the longer duration variant.) But the k-space values achieved in each of these two pulse sequences would be the same, leading to similar spatial encoding until the effects of relaxation are considered.

2. Interestingly, in the first Fourier imaging paper the authors presented both 2D and 3D imaging variants. At the time, slice selection hadn't been invented so the 2D method was a bit of a kludge in some respects, at least when compared to today's sequences. The authors used the same type of sample as Lauterbur had used - a parallel arrangement of cylindrical tubes of water- to eliminate the requirement of some sort of slice for the third dimension; symmetry of the sample did the trick instead.

Here is one of the very first 2D Fourier images of water-filled tubes, with "pixel" intensities increasing in the order (blank), ., *, A, B, C, D, E:

From: Kumar, Welti & Ernst, NMR Fourier Zeugmatography. J. Magn. Reson.18, 69-83 (1975).

Many of the early imaging papers (1970s and early 1980s) used alphanumeric plots for images because gray scale printers were hard to come by back then. Don't laugh. One day your children will think your iPhone 4 looks like a flint hand axe.

Given that most modern MR imaging methods derive from this first use of a multi-dimensional Fourier transform to convert from signals to images, this seminal work was probably sufficient to merit a Nobel prize by itself. But for largely political reasons, it appears, the Nobel committee waited until 2003 to give a prize for MRI, awarding it jointly to Sir Peter Mansfield and Paul Lauterbur in the category of Phsyiology or Medicine. Was Richard Ernst miffed? Probably not. He'd already snagged the 1991 Nobel Prize in Chemistry for his contributions to MR spectroscopy, most notably the first use of the Fourier transform (in one dimension) to interpret the signals in an NMR spectrum. (Intriguingly, there is mention of his work in MR imaging in the press release announcing the 1991 Nobel award.)

For those of you interested in the history of MRI, there's a short essay courtesy of the European Magnetic Resonance Forum.

3. While there are 16x16 "pixels" of this k-space plane, we actually sample only at discrete points, as represented by the white grid lines and the yellow border. Thus, there are actually 17x17 points in the grid. In practice we tend to sample the kx=0 and ky=0 axes and then acquire 8 points on one half and 7 on the other, for a total of 16x16. The small asymmetry isn't usually a problem, nor does it usually matter which halves of k-space have the 8 or 7 samples. Why don't we just acquire 17x17 points? These days it's mainly for historical reasons. Back when computers were slow it was generally desirable to apply "fast" FT (FFT) algorithms that require powers of two input to work. Honestly, don't worry about it. It's a tertiary effect when it comes to creating artifacts in images!

4. There isn't a one-to-one correspondence between any one point in image space and any one point in k-space. It's a more nuanced reciprocal relationship. I'll mention some of these issues in passing in this and subsequent posts, but otherwise I leave these issues to further reading.

5. In practice it is common to start with -Gy values and cover the negative portion of ky-space first. It shouldn't matter which is performed for spin warp imaging. The direction can have important consequences for EPI, but I'll mention those when we get to that point. For now, top-down or bottom-up can be considered equivalent. We just want to fill the plane.

6. As already noted, there will be a small asymmetry in k-space when utilizing the central lines of k-space and an even number of data points per dimension. Thus, in my example I have eight points on one half of k-space and seven on the other, plus the central lines. Note that in drawing the green kx lines I actually overshot by one data point in the figures! In practice that last (right-hand) column of data points wouldn't be acquired.

7. In this example, delta-kx and delta-ky were equal and so were the maximum extents of k-space, i.e. the k-matrix was square. That yields a square (16x16) image with the same field-of-view and resolution in the two dimensions. This doesn't have to be the case, of course. In the next post I will cover the FOV and resolution issues in detail, albeit using square k-space examples. Just remember that images don't have to be square and there can be occasions when a rectangular image is preferable.

Lessons from epidemiology

2011-08-05T01:18:00.000-07:00

Ben Goldacre, psychiatrist, occasional fMRIer and critic of rubbish medical research over at BadScience.net, has produced a radio documentary that covers many of the pitfalls of modern medical science:

Science: From Cradle to Grave

It's aimed at a general audience but there are important reminders for us in fMRI-land.

Confounds abound

Epidemiology is a lot like fMRI when it comes to discriminating correlation from causation. As with many areas of research using human subjects, there are usually limits to the factors that can be controlled between groups, or even across time for an individual subject.

But there are often some simple things that we can measure - like heart and respiration rates during fMRI - and thus control for. Surely we should be measuring (and ideally controlling for) as many parameters as we can get our hands on, especially when the time and expense are comparatively minor. Get as much data as you can!

Resting state fMRI: a motion confound in connectivity studies?

2011-08-05T01:14:00.000-07:00

Neuroskeptic has done us a favor and covered a recently accepted paper from Randy Buckner's lab concerning the role of motion when determining connectivity from resting state fMRI. Not only was the amount of motion found to differ systematically between male and female subjects, but this systematic difference was preserved across sessions, suggesting that it is a stable trait. The implications for group studies are discussed in the paper, and Neuroskeptic adds further perspective. It's a warning that all resting state fMRIers should heed.

Non-neural physiology.... again

There are some important limitations to consider, however. While ventricular and white matter regions were used as ways to remove some effects of heart rate and motion, the study did not acquire breathing or heart rate data and so the authors were unable to perform the more advanced BOLD-based model corrections developed by Rasmus Birn and Catie Chang (references below). Instead, they followed what might be considered the "typical" post-processing steps, including global mean signal removal. The methods are fine, my point is to highlight the limitations of the "typical" processing stream in the absence of independent physiological data.

So, could the gender differences be explained with improved physiological corrections? What about the motion correction methods in current use: might they not be up to the job we give them? We'll have to wait for further studies to find out. In the mean time, surely it only makes sense to acquire physiological data with resting state fMRI - heart rate and respiration at the very least, although there are suggestions that time course blood pressure might also be useful - and to try to explain as many confounds as possible before concluding there's a group difference due to brain activity.

References for physiological corrections:

Birn et al., Neuroimage 31: 1536 –1548, 2006.
Birn et al., Neuroimage 40: 644-654, 2008.
Chang & Glover, Neuroimage 47: 1381–1393, 2009. Also 1448 –1459 in the same issue.

Physics for understanding fMRI artifacts: Part Nine

2011-07-29T03:09:00.000-07:00

Conjugate variables redefined

In this post I'm going to provide the first part of a recipe for generating 2D images. It's going to be somewhat algorithmic. I may occasionally mention what a particular step implies, but for the most part I'm going to step through a sequence of events, produce a final recipe for you to follow, then go back and explain what some of the parts mean physically. This isn't the traditional approach to learning about k-space; most text books assume that you need to understand what it all means before you get to learn "the rules of the game." As is my wont, I'm coming at it backwards. My hope is that you will then be able to go back to your text books - I'll tell you where to look for subsequent explanations - and cement a decent understanding of the "why" of k-space, not just the "how."

Conjugate variables revisited

In Part Five of this series I introduced the Fourier transform and conjugate variables. The post focused on the most common pair of conjugate variables: frequency and time. If we have the time domain representation and we want to transform it into its frequency domain equivalent, we apply a (one-dimensional) FT, and vice versa.

But there is another pair of conjugate variables that is more useful and intuitive for imaging applications. (In this case your intuition for one of the variables may not develop until the end of this post, or later! Bear with me.) Whether it's maps, MRIs or architectural plans, the axes of an image are best described in terms of length. If we choose the centimeter as our unit of length, then FTing an axis in cm will yield an axis in 1/cm. You happen to have an intuitive notion of time, frequency and space from everyday life. Don't worry about what the reciprocal of real space means, just accept for now that it exists. We call this reciprocal space k-space because another term for 1/cm is the wavenumber, and the wavenumber is given the symbol k.

Representing pictures in reciprocal space

Let's take a random picture, in this case it's a digital photograph of a Hawker Hurricane plane. It's clearly a 2D picture. We have a digital version of it so we can do mathematical operations on it with a computer. If we do a 2D (digital) FT of the picture we get its representation in 2D k-space:

Note that a 2D FT is simply two 1D FTs performed in succession, once along each orthogonal direction. Furthermore, the image on the left is comprised of 512x512 pixels. Thus, the k-space representation on the right also comprises a matrix of 512x512 k-space values.

One can clearly see in the gray fuzz of the k-space plot on the right, above that it's a Hurricane and not the plane with which it is often confused, the Supermarine Spitfire. What, you don't see it? You don't see the squarer planform of the Hurricane's wing, or the raised section of the fuselage behind the cockpit in all those gray dots...? This is a Spitfire in k-space:

See the difference now...? Still not seeing it? Well, that is rather my point. We don't interpret k-space very naturally, our brains have evolved to interpret real space - as displayed in the conventional photos. But this silly illustration brings us to a fundamentally important concept: if we are ever given a 2D k-space representation of something and we want to reconstruct the (real space) image, all we have to do is run it through a 2D FT and we get a picture that we can interpret easily.

What does it mean to have a Hurricane or a Spitfire represented in k-space? We actually don't need to know at this point. Some of you - vision scientists and physicists, most likely - may already know a lot about the concept of spatial frequencies. If you don't, don't worry about it yet. Let's complete an exercise and then reconsider what k-space is.

We can write down the mathematics of the above transforms using what we already know about the FT (from Part Five). Considering the x dimension first, we have:

The y dimension is transformed analogously, remembering that a 2D FT is equivalent to a 1D FT along x then another 1D FT along y (or vice versa). These equations are the same as the FT that was introduced in Part Five, the only difference being a substitution of space and reciprocal space (or k-space) for time and frequency as the conjugate variables being transformed. One dimension of the Hurricane image is represented by I(x), and after Fourier transformation a rather strange pattern of gray dots is produced in the conjugate domain, S(kx). I'm using I and S intentionally; I to mean image, S to mean signal.

It's time to return to MR imaging, but with a new purpose. We know what our target is: a complete 2D plane of k-space. We also know why we want that target: so that a 2D FT will produce our final image. Thus, for MRI all we need to do is obtain a completed 2D plane of k-space (signals) and we will then be a single step - 2D Fourier transformation - away from a final image. It's time to see how we do it.

One-dimensional imaging as seen in k-space

Before jumping into 2D let's first see what happens with 1D k-space, since ultimately the two dimensions will be transformed independently (albeit simultaneously) in the 2D FT. We will first analyze the frequency-encoding gradient experiment that we saw in Part Seven:

Remember, in this exercise we aren't going to detect the signal and FT it. All we want to do is visualize the way that the spatial information is being encoded into the signal (were we to detect it). The goal is to analyze the evolving magnetization in a way that it contains spatial information in terms of kx, such that we could, if we so choose, FT it and get a representation of the x dimension of an image. It's important to make this seemingly trivial distinction because as we will see in the next post, there are times (lots of them!) when spatial encoding happens before signal detection, the spatial information being somehow "preserved" in the magnetization for when it is eventually read out.

Okay, consider again the 1D profile of a water-filled cylinder that we saw in Part Seven, as would be obtained by the above pulse sequence were we to read out and FT the signal during the gradient, Gx:

Previously we were interested in generating the green profile of the object, so we considered the entire signal that was measured in the presence of a read (or frequency encoding) gradient. Now, though, we're going to analyze the situation slightly differently by considering not the entire detected signal (that can be 1D FT'd to produce the profile) but just the signal at any instant in time as the readout gradient "makes its presence felt" by causing phase evolution of the spins depending on their position along x. In other words, we're going to look at a snapshot of this phase evolution that arises from the presence of the gradient.

First of all we need to recognize that until the signal is processed (with an FT) there is no way to distinguish between all the signals arising out of the sample. Even though the spins at position x5 are precessing at a different rate than at x8, etc., all the magnetization is precessing together and being detected simultaneously in the receiver coil. Only the FT can sort the position-dependent frequencies. So we have an interferogram - a distribution of frequencies from across the sample, all piled on top of each other. (See Part Five if you need reminding what an interferogram is.)

Still, this "all piled on top of each other" situation doesn't prevent us from understanding what's happening for any arbitrary restricted "chunk" of spins residing somewhere in the cylinder. After all, we know we will ultimately be able to sort all the signals out - with our FT - so we can consider the encoding that is happening in the time domain as a sort of proxy for doing the full analysis; a win for human insight!

The red line of spins in the above figure is as good a place as any to understand the position encoding. Noting that the red line really represents a very thin 3D slice of spins, the contribution of signal from this slice to the overall signal (the interferogram) can be described mathematically as the product of the magnetization in the slice (i.e. the number of spins present in the slice) with the amount of phase shift that the gradient has generated to this point in time.

Clearly, before the gradient is turned on there is no phase shift at all. The signal from the red slice would be M(x).dx, where M(x) is the total sample magnetization and dx is the fraction arising from just the small volume between position x and x + dx (i.e. the contribution from the red slice).

What's more, the total of all small fractions of signal - the total signal which we can detect - is now just the sum of all of the different positions along x. This summing up is the integral along x, i.e.

Once a gradient is turned on the signal accrues a phase that depends on the strength of the gradient and how long it's been enabled. In that case the signal from the each fraction (such as the red slice) becomes modulated with a phase factor that depends upon the strength of the gradient, position along x and the time for which the gradient is active:

Note that this assumes the rotating frame of reference which is why the main magnetic field strength, B0 doesn't feature; only departures from the on-resonance condition (at exactly B0) need be considered. I've also ignored 2.pi radians for simplicity.

Now the magic step. We're going to define a new variable. Let's call it k. And we will define it as:

(If you don't believe the units come out to be 1/cm, see Note 1.)

Now this is a very handy "arbitrary" definition for kx, and for a succession of reasons. Whereas in the equation for S(t) above the conjugate variables are frequency (w) and time, we are now in a position to redefine the conjugate variables as x and kx:

Nothing has really changed for the signal, we've just made a couple of simple substitutions, and just like that we have a signal equation for which the conjugate variables can be considered as x and kx instead of frequency and time. All to help us visualize what's going on, as we'll see in a moment.

But wait! The equation here for S(kx) is identical to the Fourier transforms we saw above for the Hurricane and the Spitfire (except that here I neglected the 2pi radians). Which implies, remarkably, that this "arbitrary" definition of k actually has meaning after all! When all is said and done we don't need to think about precession frequencies or magnetization phase or any of that stuff to understand spatial encoding in MRI. In fact, all we need to recognize is that when spatial encoding gradients are turned on they cause changes in k, and that when k changes it is actually tracing through the Fourier transform of the object we are seeking to image! Time to see this in pictures.

Tracing kx through time

So kx = gamma.Gx.t, where gamma is the gyromagnetic ratio for protons, the spins we're making images with. Since the gyromagnetic ratio is constant no matter what we do with the gradient strength or its duration, we can effectively ignore it (just as I already did for 2pi radians). That makes kx a remarkably simple variable to compute indeed. It is now just kx = Gx.t, which is clearly equivalent to the area of the gradient being applied, because for a rectangular-shaped gradient the area is given by base times height. (See Note 2.)

There's another important property of kx, too. Note that kx isn't constant during the gradient. In fact, before the gradient starts the value of kx is zero. Then, as the gradient is played out for increasing time, t the value of kx steadily increases. It continues to increase for as long as the gradient is on. Once the gradient is switched off the value of kx remains fixed at the last value it had under the gradient's influence; it doesn't instantly return to zero, nor does it continue to grow. Thus, we can easily see that the value of kx at t/2 is half of its eventual value after time t. In other words, kx evolves.

Now that we know how kx changes with time we can analyze the gradient echo pulse sequence we saw in the last post (Part Eight). Not only are we now in a position to dispense with the sampling period to understand what's happening to the spatial encoding, we can actually ignore the magnetization altogether! All that really matters from an encoding perspective is what is being done to the gradients. So, we can now look at the gradient echo sequence with a different perspective. Here is what happens in terms of kx during the first gradient episode:

By doing a "mental integration" during period 1 we see an increasing area under the square gradient as time evolves. Thus, in a k-space plot (which I have made 2D for obvious reasons - we will ultimately need the ky dimension to get a 2D image) the gradient traces the trajectory from the origin (at time zero) out to a maximum +kx value. The speed along the trajectory is constant (for a square gradient), but that turns out to be a bit less important than the route taken and the ultimate destination.

Next, the gradient's sign is reversed and two further periods of spatial encoding evolve:

We were already at maximum +kx so period 2 simply takes us back the way we came, and we end up back at kx=0. This position corresponds to the gradient echo top we saw in the last post, where all the dephasing has been unwound. As the gradient is left on during period 3 the trajectory takes us into negative k-space, out to maximum -kx. Again, the speed along the trajectories happened to be constant, but as before it's not the speed of the journey that matters so much as the route taken.

What's so remarkable about k-space trajectories is that they represent the gradients being applied to the point where you might reasonably ask why we bother with pulse sequence diagrams at all! Well, there are many uses for each representation. In the coming posts you will see that the k-space representation is most useful for comparing between different imaging sequences, and it can also be very useful for determining the source of some (but not all) image artifacts. In practice, we tend to use the pulse sequence and the k-space representations as a pair.

Getting off axis

I'm going to stop here for today. As an exercise, think about what sort of pulse sequence diagram might get us off the ky=0 axis and start to trace through an entire square of kx,ky space. As we saw at the start of this post, with photos of planes, the moment we have a complete square of kx,ky we are ready to apply a 2D FT and get our image out. That will be the focus of the next post.

I'll leave you with a simple homework exercise. Deduce the pulse sequence timing diagram that generates this k-space trajectory:

Answer next post! Also next post, suggested reading to reinforce the concepts you're seeing here. It usually takes several exposures to "get it." Hang in there, it really is worth it! (See Note 3.)

-------------------------

Notes

1. Gamma, the gyromagnetic ratio, for proton spin is given as 4258 Hz/Gauss, where Gauss is an old measure of magnetic field strength. (1 Tesla = 10,000 Gauss.) Noting that Hz is 1/sec and that a magnetic field gradient strength can be given in Gauss/cm, you can quickly determine that for k the units are (Gauss.sec)/(Gauss.cm.sec), leaving just 1/cm after canceling terms.

2. The shape of the gradient doesn't need to be square, it can be any shape whatsoever. The value of kx is always determined by the area under the shape being defined by the gradient. If it were a trapezoidal gradient instead of a square, we'd calculate the area under the trapezoid. Indeed, trapezoidal, triangular and even sinusoidally modulated gradients are commonly used in MRI. The general expression for determining k is:

3. I think the k-space formalism is the single most useful device you can have at your disposal to understand fMRI acquisitions. I go for years without needing to use angular momentum or product operators or the Bloch equations or relaxation theory, but k-space is something I use daily. As I mentioned in an earlier post, MR physicists think in and talk about k-space routinely. Learn the code, it's not hard. Once you do you will have a perspective on imaging and on artifacts that is invaluable. Hence this series of posts to deposit you at that point! Nearly there. One more post on 2D k-space, then another post on the particular k-space trajectory for EPI. Once we have a k-space interpretation for EPI at our disposal we can understand and differentiate artifacts in rapid fashion.

Physics for understanding fMRI artifacts: Part Eight

2011-07-16T20:18:00.000-07:00

I had initially planned to go into 2D imaging next, but after some consideration I've decided instead to tidy up a few loose ends that follow more naturally from the last post: gradient-recalled echoes and slice selection. Then, in Part Nine I promise to introduce the second in-plane dimension. This route should better allow me to bring everything together at the end of the next handful of posts and permit you to see, and understand, the EPI pulse sequence at a glance. That's the plan. Let's see if we can make it work! (See Note 1.)

Gradient-recalled echoes

In the last post I used a frequency encoding gradient, also called a readout gradient (because it's on while the signal is being recorded, or read out), to produce one-dimensional images - profiles - of water-filled objects. This isn't the typical way that the signal is acquired, however. Instead, it is typical to acquire a refocused, or echoed, signal that has a certain symmetry in time in order to obtain some experimental benefits. I'll mention these benefits later. First, let's see how the gradient echo works.

Here is a simple gradient echo pulse sequence that is adapted from the simple readout gradient-only sequence that was considered in Part Seven:

The first thing to note is that the period of data acquisition (analog-to-digital conversion) has been delayed and now occurs in concert with a readout gradient having a negative sign, rather than being coincident with the positive gradient period labeled 1 in the figure. Also, the duration of data acquisition has been doubled. So instead of acquiring a free induction decay (FID) almost immediately after the 90 degree excitation pulse, we are now acquiring an echo signal at a later time. How and why does this echo form?

To understand the gradient echo let's consider just three spatial positions across the sample, as we did in the last post. As a reminder, here's the situation under the +Gx gradient:

Now let's consider what happens to magnetization residing at these three spatial positions during the time periods labeled 1, 2 and 3 in the pulse sequence shown above. During period 1 the situation (in the rotating frame) is the same as we saw in the last post: magnetization at position x1 begins to lag whereas magnetization at position x2 increases its precessional frequency, and so it moves ahead (below, left panel):

During period 2 in the pulse sequence the gradient's sign is reversed. Thus, magnetization at position x1 starts to precess at a faster rate than that at isocenter, whereas magnetization at x2 slows down. At the end of period 2 the three chunks of magnetization have realigned. This corresponds to the top of the echo signal. It's as if no gradient whatsoever had been applied to this point in time! All the phase shifts induced by the gradients (of either sign) have been cancelled out and the phase of the entire magnetization is zero at the end of period 2 (if we assume no other sources of phase shifts arising out of the sample itself).

By leaving the gradient -Gx enabled past the echo top, however, we see that magnetization at x1 continues to outpace the rest, placing it ahead of the pack during period 3. Likewise, magnetization at x2 is now the laggard. Throughout the process the magnetization at isocenter is unaffected.

Finally, it's important to note that gradient echoes can be applied in either order, i.e. negative or positive gradient episode first. All that matters, in fact, is that the net gradient over time becomes zero at some point. The negative and positive gradients don't have to happen in succession - they can be separated by some time period - and they don't even have to be the same shape! All that matters is that their areas are equal and opposite, to achieve cancellation of accrued phases. The point at which the negative and positive gradients cancel each other is defined as the top of the echo regardless of the timing or shape of the constituent gradients.

The benefits of acquiring a gradient echo

In principle, it is possible to use either half of the signal comprising the gradient echo, or to use the frequency-encoded FID signal following the RF excitation pulse (as in Part Seven), to produce a one-dimensional image of an object. This arises because of what's called the "Hermitian symmetry" of the signal. The two halves of the gradient echo have opposite phases and, ideally, contain equivalent spatial information, plus noise.

So, why do we bother to acquire a gradient echo rather than using a single Gx episode following the RF pulse? Increased signal and reduction of artifacts. The signal boost has two components. First, delaying the ADC period allows some time to switch the scanner electronics between 'transmit' mode and 'receive' mode. This isn't a huge problem with today's modern scanner electronics but in the early days of NMR and MRI these switches could take a long time (milliseconds). But even if the switch on the RF electronics takes only tens of microseconds we need to wait for a few hundred microseconds for the frequency-encoding gradient to ramp up, meaning that we lose the ability to acquire the first part of a FID. These omitted data points would have the largest signal amplitude, meaning we could be paying a disproportionately large signal penalty. The second boost to SNR comes from the redundancy of spatial information in the two halves of the gradient echo, whereas the noise in the two halves is different. Thus, we have a signal averaging effect from FTing the two halves together, netting (ideally) a sqrt(2) ~ 40% benefit in SNR, simply by extending the duration of signal acquisition by a little bit. (A typical readout period might last only a millisecond or two, so getting a 40% boost for such a tiny amount of time is a very good deal indeed!)

The other reason for using both halves of a gradient echo is to reduce artifacts. The Hermitian symmetry in the frequency-encoded information applies only when other phase shifts aren't hampering the experiment; we want just the phase shifts being induced by the readout gradient. But the presence of magnetic field imperfections arising from magnetic susceptibility differences, e.g. between sinuses and brain, motion, and other features of real samples means that in practice using just half of a gradient echo (or a FID signal directly) will produce images having more imperfections than is necessary. Again, for such a comparatively small investment of time - often less than a millisecond - we can simply acquire both halves of a gradient echo and do a much better job.

If you're interested in pulse sequence design and want to know when it might be useful to acquire less than the whole gradient echo, see these handy pages on Partial Echo from MR-TIP, a useful website for technical descriptions. Otherwise, don't worry about it. For 99.9% of fMRI experiments we acquire entire gradient echo signals. Furthermore, gradient echoes don't just appear in frequency encoding either. They have a wide range of uses, one of which appears in the next section. And when we get to the EPI pulse sequence you'll see that balancing positive and negative gradient episodes in a train of (gradient) echoes is the principal workhorse of the total spatial encoding scheme.

Slice selection

In Part Seven we made profiles of three-dimensional objects by conveniently orienting the frequency encoding gradient orthogonal to the long axis of axially symmetric shapes, thereby rendering a one-dimensional profile interpretable. That may be acceptable when one dimension is constant, such as the axis of a cylinder, but what if the cylinder had a bulb at the end, or a narrow waist? And what if it's not a cylinder at all, but a brain? It would greatly complicate interpretation of fMRI data if signal arose from the neck and shoulders in addition to the head, wouldn't you think? (fMRI skeptics have enough ammunition as it is.)

A solution to this issue is to restrict the extent of sampled magnetization in the third dimension via slice selection. This process is akin to electronically slicing up the subject into a near arbitrary number and thickness of slices. I say near arbitrary because there are practical limits, as you'd expect. (See Note 2.)

To achieve slice selection, we first need to select a suitable RF pulse that will excite a band, or notch, of frequencies. If we apply such an RF pulse in the presence of a magnetic field gradient, the combination of RF pulse and gradient will excite just a slab of magnetization rather than the entire object. That's because, as we have seen in this post and the last one, linear magnetic field gradients have the effect of spreading out frequencies with a (linear) spatial dependence, providing an easy correspondence between frequency and position along the gradient direction.

Before considering further the effect of the RF pulse in concert with a gradient, let's look in a little more detail at the RF pulse shapes that might be useful for slice selection. Fourier pairs are a good place to start. In Part Six we saw the sinc function and its Fourier pair, the boxcar function:

If the goal is a slab of excited frequencies with little to no effect outside of the slab, then the boxcar in the bottom-right corner of the above figure looks like it might fit the bill, eh? To get a frequency domain boxcar in 1D - which would correspond to a slab of magnetization in 3D - then all we need to do is modulate the amplitude of the RF excitation pulse so that it is sinc-like as it's played out in time. (See Note 3 if you need a refresher on RF pulses.) The frequency of the oscillating magnetic field pulse is still in the radiofrequency range, all we're doing is modulating its amplitude with the shape that we want, rather than simply turning the RF pulse on and off, as was done to achieve a "square" RF pulse earlier in this post and in Part Seven. (See Note 4.)

Now all that remains to be done is to spread out the magnetization across the sample in such a way that the boxcar's frequency bandwidth corresponds to a range of frequencies across the target region of interest. Let's say that we want to slice a cylinder into a succession of transverse slices. Then all we need to do is turn on the sinc-shaped RF pulse in the presence of a gradient along Z, the cylinder's axis, and then place the resulting boxcar at the frequency range corresponding to the location of the target slice (along Z):

If we want to change the slice thickness we have two options. Either we can change the duration of the sinc-shaped RF pulse (because it has an inverse relationship in the frequency domain, so that a longer sinc pulse excites a narrower frequency notch and vice versa), or we can change the gradient strength. If we double the slice select gradient (Gss) amplitude we will halve the slice thickness (keeping the RF pulse constant):

There is of course a limit to how far Gss can be increased. For typical sinc RF pulses of 4-8 ms duration and maximum Gss of 40 mT/m, the thinnest slices that can be achieved for EPI are about 1 mm on most scanners. (See Note 5.)

In order to move the slice along the Z direction we adjust the central frequency - the carrier frequency - of the sinc RF pulse. In the hand-drawn figure of the sinc-boxcar Fourier pairs above you'll notice the time domain sinc (bottom-left) is centered at zero time. In the rotating reference frame (for 123 MHz) this achieves slice selection at the magnet isocenter (i.e. exactly at 123 MHz, which is 0 Hz in the rotating frame), as drawn in the two cartoons of the cylinder, above. If the sinc's carrier frequency is increased, the resulting slice will slide towards the +Z direction, which is towards the top of the cylinder as shown. (It's also towards the subject's feet on my magnet. See Note 6 in Part Seven for the gradient directions.) Note that in moving the slice by changing its central frequency the slice thickness isn't altered.

Remember: gradients cause dephasing!

There's one more important issue to consider before we've achieved practical slice selection. The RF excitation pulse is being played out in the presence of Gss in order to spread frequencies along one dimension of the object. But as we've seen above and in Part Seven, gradients produce dephasing that can lead to signal attenuation, as well as their intended frequency encoding properties. To reverse the dephasing we now simply apply Gss as a gradient echo, akin to the frequency encoding echo in the previous section (except that we're not acquiring signals with this echo, just selecting magnetization):

I've shown a three-lobed sinc having an arbitrary flip angle, alpha. And it's concurrent with a slice select gradient, Gz. The inclusion of a reversed sign gradient (-Gz above) means that we will have achieved the desired slice selectivity during the RF pulse (and positive Gz lobe, above), then we negate the concomitant dephasing with a separate, reversed gradient episode after the RF pulse. In this way we get what we want (a slice) and eliminate what we don't want (signal attenuation).

You'll see this form of slice selection scheme in nearly all MRI pulse sequences, not just EPI. The area under the refocusing gradient (which is often called Gssr) needs to match the area under the Gss gradient from the mid-point of the RF pulse. In my figure I've drawn the Gssr episode with equal amplitude but opposite sign to the Gss episode, meaning that the duration of Gssr needs to be half that of Gss. Often, however, a little bit of efficiency can be achieved by making the amplitude of Gssr larger, then decreasing its duration in proportion to maintain its area. Very often the Gssr duration is matched to a gradient with another task, and the amplitude is adjusted accordingly. However it's done, the only rule is that the area under the refocusing gradient matches the gradient area from the mid-point of the RF pulse. (See Note 6.)

In the slice selection pulse sequence above I showed the RF pulse flip angle as alpha to reinforce the point that this gradient echo phase reversal scheme is necessary whether you are using a 30 degree flip angle or a 90. The only time the gradient echo isn't required is for RF pulses used for refocusing themselves, typically 180 degree RF pulses. But since these aren't commonly used for fMRI I won't go into more detail here. (See Note 7 if you're curious.)

Next post: two-dimensional imaging and k-space. No, it's not a cosmic phenomenon, it's a handy and intuitive way to visualize pulse sequences.

----------------

Notes:

1. I was wondering whether I should go through the spin echo as a precursor for the gradient echo, but decided against it on the basis that spin echoes aren't crucial for understanding what passes for standard fMRI experiments these days. I am also working from the assumption that most readers will have been exposed to an introductory course on NMR or MRI and so have a rudimentary understanding of - at least have seen at some point - spin echoes. Turns out I'm not the first to think about this quandary. ReviseMRI deals with the question of whether to teach gradient or spin echoes first, here. I concur with ReviseMRI: gradient echo first. Not only that, but when I redo my graduate lectures for later this year I'm going to relegate spin echoes into the 'bonus' material. With apologies to my Berkeley colleague Erwin Hahn, there are too many other concepts for fMRI that should be covered first, especially when lecture time is tight.

2. For typical slice-select gradient strengths on a typical 3 T scanner, with RF pulses that are between about 4-8 ms duration, it's possible to achieve 2D slices as thin as 1 mm or so. A more usual slice thickness for fMRI is in the range 2.5-4 mm, and there are two main reasons for this. Firstly, the volume of tissue scales with the voxel dimensions, and making any voxel dimension (including the slice thickness) too small can render the overall signal-to-noise ratio (SNR) unacceptably low to detect functional signal changes. Secondly, in order to sample the hemodynamic response following an event, we need to sample the entire brain every 2-3 seconds. We cannot, therefore, spend arbitrarily long acquiring loads of very thin slices. Instead, we have to compromise and acquire fewer thicker slices.

3. If you can't remember how an RF pulse "works," by causing additional precession about an effective magnetic field, then you might want to break off and re-read an introductory text. For an fMRI lab I'd recommend Rick Buxton's book, Introduction to Functional Magnetic Resonance Imaging (see pages 134-9 in the 2nd edition) and Functional Magnetic Resonance Imaging by Huettel, Song and McCarthy (see pages 106-9 in the 2nd edition). Both of these excellent books are eminently readable. Not quite beach reading, but as close as you're likely to come!

4. In reality the sinc profile must be truncated, typically after just three, four or five side lobes on either side. This causes the slice profile to be very slightly trapezoidal rather than perfectly square. Furthermore, there are some small ripples outside of the trapezoid, on either side, arising from the truncation. (Recall that squares in the time domain cause sinc-shaped ripples in the frequency domain!) For this reason it is common to acquire interleaved slices (odds then evens, say) to allow half of TR relaxation between adjacent slice excitation. However, some labs prefer to acquire sequential slices because there is less effect of through-plane motion perturbing the T1 steady state. If you want more information on interleaved versus sequential slices, see the section "Should I use interleaved or sequential slices for fMRI?" in my user training/FAQ notes, available via this post. I may do a separate post here on RF excitation profiles at a much later date, but it's worth noting that imperfect RF profiles come a long way down list of concerns and artifacts that hamper EPI and fMRI.

5. For EPI it's impractical to increase the sinc RF pulse much past 8 ms for a number of reasons, one of which is that it increases the minimum echo time (TE), thereby placing a constraint on the number of slices that can be acquired for a given TR. There may also be hardware limitations.

6. There is a method, called z-shimming, which seeks to tweak the slice select refocusing gradient, and maybe some other gradients, in order to combat the dephasing effects of susceptibility gradients arising from the sample itself, not just the effects of Gss. There are, naturally, pros and cons to z-shimming. It's not a technique in widespread use but I'll find the time to cover it in a later post since I can't find a good review article on it, only technical primary references.

7. Refocusing RF pulses that are used to produce spin echoes, and which are usually (but not always) 180 degree flips, do not require the reversed sign rephasing gradient episode that excitation RF pulses do. That's because the dephasing that happens during the first half of the 180 degree RF pulse is matched by the rephasing that happens during the second half; refocusing pulses are inherently balanced.

Physics for understanding fMRI artifacts: Part Seven

2011-07-11T13:12:00.000-07:00

Magnetic field gradients and one-dimensional MRI

Now that you have a basic understanding of the Fourier transform and some of the practical matters that arise from digital signals, it's time to look at a basic imaging pulse sequence and even make some simple images. We're going to use frequency encoding only for the time being, and for now we're going to make one-dimensional images (also called profiles) so that we can introduce an alternative form of timing diagram to represent a pulse sequence.

A magnetic field gradient alters the local resonance frequency

When a sample is placed into the magnet, all the protons (1-H nuclei) resonate at a near-identical frequency. At 3 T that resonance frequency is approximately 123 MHz, as given by the Larmor equation. If we then impose a magnetic field gradient across the sample - your subject's head, say - instead of having the same resonance frequency uniformly across the brain, there will now be a linear dependence in space (see Note 1):

In a real image we might consider 64 different positions along x. These would define the voxels in one (in-plane) dimension of the image. But for the time being we'll consider just three points in the x direction: the central point, and one point either side.

At the center of the magnet the gradient has no net effect, so the resonance frequency at that point is still 123 MHz. We call this point the null crossing, because all three linear gradients, X, Y and Z, are engineered to have no effect here. (See Note 2.) And to keep things symmetric, the gradient null crossing is placed in the geometric center of the magnet - the isocenter - because that's where the main magnetic field has been engineered to be most homogeneous, and we want to do all our imaging in that location to get the best scanner performance.

At the magnet isocenter, then, the resonance frequency is left unchanged by Gx. Otherwise, the resonance frequency is now x-dependent according to:

At isocenter the sample experiences exactly 3 T whereas to the left of center (at position x1) the sample is at a field slightly less than 3 T (there's partial cancellation where the gradient field opposes the main magnetic field), and to the right of center (at x2) it's at a magnetic field slightly greater than 3 T. (See Note 3.) The spatial dependence of the resonance frequency is linear. (If you think I'm doing an obvious point to death, take a look at Note 4.)

We can simplify the representation a little by considering a rotating reference frame, rather than the static reference frame of the magnet (which is known as the laboratory frame). We deal with a rotating frame all the time. Hands up who navigates relative to static positions on the earth? We conveniently discount the fact that we're rotating at a thousand miles an hour about the earth's axis, or that the earth is rotating about the sun, etc. Likewise, in MRI we take a reference rotational frequency - the nominal resonance frequency at the center of the magnet, which is 123 MHz here - and we then consider frequencies slightly faster or slower than the reference. These frequency differences manifest as phase shifts for the spins in the rotating frame picture. (See Note 5.)

Below left, I've represented the resonance frequencies for three chunks of magnetization, one located at isocenter and one each for positions x1 and x2, relative to the Gx gradient illustrated above. In this vector diagram, equilibrium magnetization (before RF excitation) would be indicated by a vertical arrow pointed along the +z axis; all three chunks of magnetization would reside there. But after RF excitation, and following the imposition of the gradient Gx, the magnetization will be spread out and there will be local differences in resonance frequency, depending on position x:

LEFT: Magnetization in the rotating reference frame for three positions along x after imposition of a gradient Gx. RIGHT: The resonance frequencies in the lab frame (white) and in the rotating frame (colors) for the three positions.

In general then, the spatially-dependent frequency is given by:

Acquiring an MR signal in the presence of a gradient

Achieving the situation you've just encountered can be represented on a timing diagram, such as this one:

Starting with equilibrium magnetization (which isn't indicated, just assumed), a 90 degree RF excitation pulse rotates all proton spins in the sample from the z axis into the transverse plane of the rotating reference frame. The gradient, Gx is then activated, differentiating resonance frequencies with an x spatial dependence as we just saw. (See Note 6.) There's also T2 relaxation going on, but we will be able to ignore this for the purposes of spatial encoding; I include it below purely to show how the effect of the gradient and the innate relaxation have similar mathematical representations. Provided the decaying magnetization persists for enough time to record it - with our analog-to-digital converter, ADC - we're in good shape. (See Note 7.)

In the absence of a gradient, the transverse magnetization decays with T2 relaxation only:

where Mxy(t=0) is the magnetization immediately after the 90 degree excitation pulse.

If we turn on a gradient after the excitation pulse then there is an additional dephasing, resulting in further attenuation of net magnetization as well as providing the desired spatial encoding:

Spatial information is encoded because, as we've seen above, frequency (w) is a function of position, x when Gx is turned on. The term in yellow, arising from the gradient we've enabled, is just another dephasing term, albeit one that now contains some useful spatially-dependent information. (See Note 8.)

Now all we have to do is "decode" the spatial information resident in the yellow exponent. This can be accomplished quite simply by Fourier transformation of Mxy because, as we saw in Part Five, an FT will produce a one-dimensional plot of the frequency content of a time domain signal. If we ignore the T2 relaxation for simplicity (by assuming that it's very small during the gradient encoding and data acquisition period) then the signal intensity at each position along the one-dimensional frequency plot will be proportional to the spin density at that position, as represented by Mxy(t=0) in the equation above. Let's see this process in action.

One-dimensional imaging in pictures

We'll start with a cylinder filled with water so that the magnetization is uniform throughout, and zero outside. The x gradient will be applied perpendicular to the cylinder's axis:

The gradient creates a range of frequencies for each position, x across the cylinder's cross-section:

Next, we need to determine how much signal we will receive from each position along x, which we do with the Fourier transform. From the figure above let's suppose there are n=64 positions along x to consider. After FT the amount of signal in the frequency domain will be directly proportional to the spin density in each of the 64 planes orthogonal to the x direction. Clearly, at positions (such as x1 and x2 above) outside the sample there will be zero signal. Within the sample the spin density is uniform, but the shape of the cylinder means that each plane will 'sample' different volumes of water. The plane that samples the central section of the cylinder along x will obviously contain the highest number of protons (hence largest magnetization), thus will be the tallest portion of the frequency plot.

If we do this process 64 times and then join up the signal amplitudes in the frequency domain we end up with a one-dimensional image of the cross-section of the cylinder, as shown by the green ellipse below. I've indicated one plane, in red, through the cylinder to show the origin of the profile explicitly:

That was easy, right? Good. Then try two new samples on for size. Not only will we add a second water-filled sample having a trapezoidal cross-section, let's also measure a profile of the two objects along the y axis as well as x. For the y profile we'd do a separate experiment in which we apply a gradient, Gy instead of Gx. (See Note 9.)

A final point for today. In the last example I didn't mention the length of the two objects being imaged. What would be the relative heights of the cylinder and the trapezoid? Well, if you assume the light blue color indicates equivalent spin density for the water in each shape then, given that the peak intensities in the x profiles are the same, we would have to conclude that the lengths are equivalent. (The total number of spins in a plane through the center of each sample is the same.) Had I drawn the cylinder's profile along x with, say, twice the maximum amplitude of the trapezoid's then you could have concluded that the cylinder were twice the length of the trapezoid (again assuming equal spin density). If the length of the cylinder had been twice the length of the trapezoid, what would the y profile have looked like...?

Next post: gradient-recalled echoes and slice selection. After that, two-dimensional imaging!

----------------------

Notes:

1. We abbreviate magnetic field gradient to gradient pervasively in fMRI. There are other gradients in MRI, depending on what you're doing (e.g. temperature gradients), but 99.9% of the time when an MRI person uses the term gradient a magnetic field gradient is assumed.

2. If you're curious what MRI gradient coils look like and how they work, there's a simple introduction in the book, "Functional Magnetic Resonance Imaging" by Huettel, Song and McCarthy. In the second edition it's on pages 38-40. Modern gradient coil designs have come a long way from the basic conceptual designs that you'll mostly find in books, however. Vendors treat their designs like state secrets, and with good reason. The majority of MRI physicists would probably agree that stable, high fidelity, high powered gradients are the most critical aspect of overall MRI performance. Perhaps I'll do a separate post on gradients at a later date. Unless you're responsible for specifying and buying a new MRI scanner, or you're charged with doing your center's QA, you really don't need to understand the engineering of the gradients in order to do fMRI experiments.

3. On a typical 3 T scanner the maximum magnetic field gradient amplitude is around 40 mT/m. So the gradients are three orders of magnitude smaller than the main magnetic field. You don't need a whole lot of magnetic field strength to encode spatial information!

4. In a nutshell, this is what Paul Lauterbur did in 1973 to get himself a Nobel Prize (awarded in 2003, jointly with Sir Peter Mansfield). Okay, so Lauterbur had to do a little more work than this – he imposed a series of gradients at different orientations relative to two tubes of water, then used a back-projection algorithm to compute a 2D image of the tubes – but a linear modification of the Larmor equation was, in essence, his crucial insight. And it’s pretty darn funny in hindsight - at least it is to a reformed spectroscopist like myself - because for the first thirty-odd years of NMR, from the first successful measurements by Bloch and Purcell in 1946, everybody had been working really hard to get the flattest, most homogeneous magnetic field possible. Then along comes Lauterbur, tilts the magnetic field so that it has a known (linear) spatial dependence, invents MRI and asks the receptionist if he might now have his Nobel Prize, please. No wonder Damadian was so ticked off.

Lest you should get despondent when your latest fMRI paper gets rejected, apparently Lauterbur ran into obstacles trying to patent his idea, and then had his first attempt to publish in Nature rejected. Ah, peer review. You have to laugh, don't you? Anyway, here's Lauterbur's image, the very first MRI:

From: P.C. Lauterbur, Nature 242, 190-1 (1973).

5. The rotating reference frame, denoted by primes, rotates about the lab frame such that +z and +z' are parallel and aligned (by convention) along the main magnetic field direction, B0. In this representation the x and y lab axes are whirling around z at 123 MHz, rendering the rotation of the magnetization, and the rotating axes x' and y', static in this view. Then any precessional frequency changes of magnetization will appear as phase shifts relative to the (now static) x' and y' axes.

6. On my magnet at least, the magnetic field gradients form a right-handed set. That is, if you hold up the index finger, middle finger and thumb of your right hand so that they are orthogonal, your index finger indicates +Z, the magnetic field direction (towards the front of the magnet), your middle finger indicates the +X direction, which is the subject's left (assuming the subject is in head first, supine), and your thumb indicates +Y, towards the top of the magnet.

7. If you are already familiar with pulse sequences then you will immediately spot the lack of a gradient-recalled echo, which is how we usually acquire MRI data for several experimental reasons, and it's T2* rather than T2 relaxation that's happening. I'm keeping it simple for now. Worry not, the GRE will be along in the next post. Likewise, I'm not using slice selection here because it's not necessary yet either. I want to get through the in-plane k-space representation first. We're going to go in the order: frequency encoding, phase encoding, slice selection.

8. I want to reinforce the point that the term in yellow, arising from the gradient we've enabled, is just another dephasing term. If we weren't interested in its spatial encoding properties (via the frequency distribution along x) then there is another "use" for this gradient pulse: it enhances signal attenuation. Furthermore, in MRI a gradient is a gradient is a gradient! Whether you are able to control it from the scanner console, as with Gx, or it arises intrinsically because of the properties of your sample (e.g. your subject's skull, or his blood vessels), the effect is the same as far as the spins are concerned. Really, the only differences between our linear gradients and the gradients that arise intrinsically is this aspect of control, the known (linear) spatial dependence of Gx,y,z, and the length over which the gradients act. Gx,y,z are imposed across the entire sample, whereas the subject's skull probably has quite local effects. Beyond these differences, however, the spins don't care. If they experience a gradient arising from the skull or from Gx they will dephase all the same. This is good and bad news; bad news because it will cause difficulties for image localization, good news because it creates this thing called BOLD and enables one version of fMRI. More on all of this stuff in many later posts.

9. Why can't we turn on Gx and Gy at the same time, and why must we acquire each profile in separate experiments? If Gx and Gy are on simultaneously we would then have one gradient that is the vector sum of Gx and Gy, not separate x and y gradients! The combined gradient would create another profile in a direction 45 degrees between x and y.

MRI Claymation

2011-07-02T20:07:00.000-07:00

It's a long weekend here in the US of A, it's hotter than Hades here in northern California (they said there would be a fog-cooled sea breeze! I want my money back!), and I am only halfway through the next post in the background physics series on account of having spent a very pleasant week in Quebec at the Human Brain Mapping conference. So, in lieu of anything more useful at short notice, I thought I'd share a truly awesome video I just found online, courtesy of Andre van der Kouwe and colleagues at MGH. The first two minutes demonstrate the method - surface renderings from MRIs of clay figures - and then it gets really fun: MRI making an image of itself.

427 views in two years simply doesn't do this work justice. Let's fix that!