Estimation and removal of movement-by-susceptibility induced variance in fMRI time series.

- Rationale
- Previous Work
- Assumptions behind Unwarp
- The "Realign & Unwarp" function
- References
- Acknowledgments
- Downloads
- Installation instructions

Even after realignment there is considerable variance in fMRI time series that covary with, and is most probably caused by, subject movements. It is also the case that this variance is typically large compared to experimentally induced variance. Anyone interested can include the estimated movement parameters as covariates in the design matrix, and take a look at an F-contrast encompassing those columns. It can be quite dramatic. The result is loss of sensitivity, and if movements are correlated to tas specificity. I.e. we may mistake movement induced variance for true activations. Because the movement induced variance is often very large compared to "true" activations false positives may ensue even from a relatively modest correlation between task and movement.

The problem is well known, and several solutions have been suggested. A quite pragmatic (and conservative) solution is to include the estimated movement parameters (and possibly squared) as covariates in the design matrix. Since we typically have loads of degrees of freedom in fMRI we can usually afford this. The problem occurs when movements are correlated with the task, since the strategy above will discard "good" and "bad" variance alike.

The "covariate" strategy described above was predicated on a model where variance was assumed to be caused by "spin history" effects, but will work pretty much equally good/bad regardless of what the true underlying cause is.

Others have assumed that the residual variance is caused mainly by errors caused by the interpolation kernel in the resampling step of the realignment. One has tried to solve this through higher order resampling (huge Sinc kernels, or k-space resampling).

The "Unwarp" toolbox is based on a different hypothesis regarding the residual variance. EPI images are not particularly faithful reproductions of the object, and in particular there are severe geometric distortions in regions where there is an air-tissue interface (e.g. orbitofronal cortex and the anterior medial temporal lobes). In these areas in particular the observed image is a severely warped version of reality, much like a funny mirror at a fair ground. When one moves in front of such a mirror ones image will distort in different ways and ones head may change from very elongated to seriously flattened. If we were to take digital snapshots of the reflection at these different positions it is rather obvious that realignment will not suffice to bring them into a common space.

The situation is similar with EPI images, and an image collected for a given subject position will not be identical to that collected at another. Hence, even after a "successful" realignment there will be residual variance caused by the object having different shape at different time points. We call this effect susceptibility-by-movement interaction. "Unwarp" is predicated on the assumption that the susceptibility-by-movement interaction is responsible for a sizeable part of residual movem nt related variance.

Assume that we know how the deformations change when the subject changes position (i.e. we know the derivatives of the deformations with respect to subject position). That means that for a given time series and a given set of subject movements we should be able to predict the "shape changes" in the object and the ensuing variance in the time series. It also means that, in principle, we should be able to formulate the inverse problem, i.e. given the observed variance (after realignment) and known (esti ated) movements we should be able to estimate how deformations change with subject movement.

We have made an attempt at formulating such an inverse model, and at solving for the "derivative fields". A deformation field can be thought of as little vectors at each position in space showing how that particular location has been deflected. A "derivative field" is then the rate of change of those vectors with respect to subject movement. Given these "derivative fields" we should be able to remove the variance caused by the susceptibility-by-movement interaction. Since the underlying model is so re tricted we would also expect experimentally induced variance to be preserved. Our experiments have also shown this to be true. Indeed one particular experiment even indicated that in some cases the method will reintroduce experimental variance that had been obliterated by movement related variance.

In theory it should be possible to estimate also the "static" deformation field, yielding an unwarped (to some true geometry) version of the time series. In practice that doesn't really seem to work, hence the method deals only with residual movement related variance induced by the susceptibility-by-movement interaction.

This means that the time-series will be undistorted to some "average distortion" state rather than to the true geometry. If one wants additionally to address the issue of anatomical fidelity one should combine Unwarp with a measured field-map. A field-map in the format that Unwarp expects can be created using the FieldMap toolbox.

The description above can be thought of in terms of a Taylor expansion of the field as a function of subject movement. Unwarp alone will estimate the first (and optionally second, see below) order terms of this expansion. It cannot estimate the zeroth order term (the distortions common to all scans in the time series) since that doesn't introduce (almost) any variance in the time series. The measured fieldmap takes the role of the zeroth order term. Refer to the FieldMap toolbox and the documents FieldMap.man and FieldMap_principles.man for a description of how to obtain fieldmaps in the format expected by Unwarp.

If we think of the field as a function of subject movement it should in principle be a function of six variables since rigid body movement has six degrees of freedom. However, the physics of the problem tells us that the field should not depend on translations nor on rotation in a plane perpendicular to the magnetic flux. Hence it should in principle be sufficient to model the field as a function of out-of-plane rotations (i.e. pitch and roll). One can object to this in terms of the effects of shimming (object no longer immersed in a homogenous field) that introduces a dependence on all movement parameters. In addition SPM/Unwarp cannot really tell if the transversal slices it is being passed are really perpendicular to the flux or not. In practice it turns out thought that it is never (at least we haven't seen any case) necessarry to include more than Pitch and Roll. This is probably because the individual movement parameters are typically highly correlated anyway, which in turn is probably because most heads that we scan are attached to a neck around which rotations occurr.

On the subject of Taylor expansion we should mention that there is the option to use a second-order expansion (through the defaults) interface. This implies estimating also the rate-of-change w.r.t. to some movement parameter of the rate-of-change of the field w.r.t. some movement parameter (colloquially known as a second derivative). It can be quite intresting to watch (and it is amazing that it is possible) but rarely helpful/necessarry.

In addition to inducing residual (after realignment) movement-related variance, movement-by-susceptibility distortion changes may bias the estimation of the movement. For each iteration Unwarp gets a better idea of the true shape of each scan, and can potentially get a better estimate of movement. In Unwarp there is an option (default) to re-estimate (do a new realign) the movements between each iteration of estimating the fields. Our testing indicates that this is a good idea.

In the defaults there is also an option to include Jacobian intensity
modulation when estimating the fields. "Jacobian intensity
modulation" refers to the dilution/concentration of intensity that ensue
as a consequence of the distortions. Think of a semi-transparent coloured
rubber sheet that you hold against a white background. If you stretch a part of
the sheet (induce distortions) you will see the colour fading in that
particular area (you can also think of future appearance of a tattoo obtained
when young and slim). In theory it is a brilliant idea to include also these
effects when estimating the field (see e.g. Andersson et al, NeuroImage
20:870-888). In practice for this specific problem it is ** NOT** a good idea.

It should be noted that this is a method intended to correct data afflicted by a particular problem. If there is little movement in your data to begin with this method will do you no good. If on the other hand there is appreciable movement in your data (>1mm or >1deg) it will remove some of that unwanted variance. If, in addition, movements are task related it will do so without removing all your "true" activations.

The method attempts to minimise total (across the image volume) variance in the data set. It should be realised that while (for small movements) a rather limited portion of the total variance is removed, the susceptibility-by-movement interaction effects are quite localised to "problem" areas. Hence, for a subset of voxels in e.g. frontal-medial and orbitofronal cortices and parts of the temporal lobes the reduction can be quite dramatic (>90%).

It should also be noted that there are several reasons for residual movement related variance. Notably:

- Susceptibility-distortion-by-movement interaction: Is what we deal with here.
- Susceptibility-dropout-by-movement interaction: The susceptibility induced field inhomogenity will, in addition to distortions, cause signal loss due to through-plane dephasing (which will not be rephased by encoding gradients that are all in-plane). This signal modulation will also change with subject position.
- Spin-history effects: The signal will depend on how much of
longitudinal magnetisation has recovered (through T
_{1}relaxation) since it was last excited (short T_{R}→low signal). Assume we have 42 slices, a T_{R}of 4.2seconds and that there is a subject z-translation in the direction of increasing slice # between one excitation and the next. This means that for that one scan there will be an effective T_{R}of 4.3seconds, which means that intensity will increase. - Slice-to-vol effects: The rigid-body model that is used by most motion-correction (e.g. SPM) methods assume that the subject remains perfectly still for the duration of one scan (a few seconds) and that any movement will occurr in the few μs/ms while the scanner is preparing for next volume. Needless to say that is not true, and will lead to further apparent shape changes.

Unwarp will deal only with the first of these, which means that even after that correction the other components will remain. It is difficult to say which of these effetcs dominate. It will depend on ones scanner, and even on the specific data set.

Friston KJ, Williams SR, Howard R, Frackowiak RSJ and Turner R (1995)
Movement-related effect in fMRI time-series. *Magn Reson Med* **35**:346-355

Jezzard P and Balaban RS (1995) Correction for geometric distortions in
echoplanar images from B0 field variations. *Magn Reson Med* **34**:65-73

Wu DH, Lewin JS and Duerk JL (1997) Inadequacy of motion correction
algorithms in functional MRI: Role of susceptibility-induced artefacts.
*J Magn Reson Imag* **7**:365-370

The software was developed in part while JLRA was supported by STINT (Stiftelsen for internationalisering av högre utbildning och forskning.)

This software has been tested by the authors on data collected on two different make MR scanners. Undoubtedly there are numerous bugs that will surface when more people start using it. We are grateful for bug reports, and other input, and the more we receive the faster bugs will be weeded out.

Unwarp Toolbox for SPM5 can be downloaded with the latest version of SPM5

Since Unwarp is a part of main SPM2 this distribution will (eventually, when bugs are fewer) be part of the general SPM updates. Meanwhile we recommend the following method for installing the update.

- Create a local directory, e.g. /usr/local/
*'whatever'*/Unwarp_updates - Download the Unwarp_ddmmyy.tar.gz file into that directory
- Unzip and untar, on a Linux machine

*..]$ gunzip Unwarp_ddmmyy.tar.gz*

*..]$ tar xvf Unwarp_ddmmyy.tar.gz* - If you are on a Macintosh, or some other exotic brand, you will
need to compile spm_dilate_erode.c. It should work to simply
*cd*into your Unwarp_updates directory, start Matlab and at the Matlab prompt type

>>*mex spm_dilate_erode.c* - Make sure that /usr/local/
*'whatever'*/Unwarp_updates is included in your Matlab path ahead of your main SPM distribution. - If (very optional) you want the help in the SPM graphical help (for the "Realign & Unwarp" button) to reflect Realign and Unwarp you should rename mod_spm.m spm.m and put it in your main SPM distribution directory. DONT rename it and leave it in the Unwarp_updates directory!

Chloe Hutton (chutton at fil.ion.ucl.ac.uk)

John Ashburner (john at fil.ion.ucl.ac.uk)

Robert Turner (r.turner at fil.ion.ucl.ac.uk)

Karl Friston (karl at fil.ion.ucl.ac.uk)