Spatial Preprocessing

1. Spatial Preprocessing Ged Ridgway University College London Thanks to: John Ashburner, Klaas Enno Stephan, and the FIL Methods Group.Meike Grol, Chloe Hutton, Jesper Andersson, Floris de Lange,Miranda van Turennout, Lennart Verhagen, …

2. Terminology of fMRI Condition B Condition A t Structural (T1) images: - high resolution - to distinguish different types of tissue • Functional (T2*) images: • - lower spatial resolution • to relate changes in BOLD signal • to an experimental manipulation  Time series: A large number of images that are acquired in temporal order at a specific rate

3. Terminology of fMRI subjects sessions runs single run volume slices voxel TR = repetition time time required to scan one volume

4. Terminology of fMRI Scan Volume: Field of View (FOV), e.g. 192 mm Slice thickness e.g., 3 mm Axial slices Matrix Size e.g., 64 x 64 3 mm In-plane resolution 192 mm / 64 = 3 mm 3 mm Voxel Size (volumetric pixel) 3 mm

5. fMRI Time-series Movie

6. Overview of SPM Analysis Statistical Parametric Map Design matrix General Linear Model Parameter Estimates fMRI time-series MotionCorrection Smoothing SpatialNormalisation Anatomical Reference

7. Preprocessing in SPM • Realignment • With non-linear unwarping for EPI fMRI • Slice-time correction • Coregistration • Normalisation • Segmentation • Smoothing

8. The many guises of affine registration • Manual reorientation • Rigid intra-modal realignment • Motion correction of fMRI time-series • Within-subject longitudinal registration of serial sMRI • Rigid inter-modal coregistration • Aligning structural and (mean) functional images • Multimodal structural registration, e.g. T1-T2 • Affine inter-subject registration • First stage of non-linear spatial normalisation • Approximate alignment of tissue probability maps

9. Other types of registration in SPM • Non-linear spatial normalisation • Registering different subjects to a standard template • Unified segmentation and normalisation • Warping standard-space tissue probability maps to a particular subject (can normalise using the inverse) • High-dimensional warping • Modelling small longitudinal deformations (e.g. AD) • DARTEL • Smooth large-deformation warps using flows • Normalisation to group’s average shape template

11. Manual reorientation Image “headers” contain information that lets us map from voxel indices to “world” coordinates in mm Modifying this mapping lets us reorient (and realign or coregister) the image(s)

12. Manual reorientation

13. Manual reorientationInterpolation Nearest Neighbour (Bi)linear (Bi)linear Sinc

14. Interpolation • Applying the transformation parameters, and re-sampling the data onto the same grid of voxels as the target image (usually) • AKA reslicing, regridding, transformation, and writing (as in normalise -write) • Nearest neighbour gives the new voxel the value of the closest corresponding voxel in the source • Linear interpolation uses information from all immediate neighbours (2 in 1D, 4 in 2D, 8 in 3D) • NN and linear interp. correspond to zeroth and first order B-spline interpolation, higher orders use more information in the hope of improving results • (Sinc interpolation is an alternative to B-spline)

15. Manual reorientation – Reslicing Reoriented(1x1x3 mm voxel size) Resliced (to 2 mm cubic)

16. Quantifying image alignment • Registration intuitively relies on the concept of aligning images to increase their similarity • This needs to be mathematically formalised • We need practical way(s) of measuring similarity • Using interpolation we can find the intensity at equivalent voxels • (equivalent according to the current transformation parameter estimates)

17. Voxel similarity measures Pairs of voxel intensities • Mean-squared difference • Correlation coefficient • Joint histogram measures

18. Intra-modal similarity measures • Mean squared error (minimise) • AKA sum-squared error, RMS error, etc. (monotonically related). AKA “least squares” • Assumes simple relationship between intensities • Optimal (only) if differences are i.i.d. Gaussian • Okay for fMRI realignment or sMRI-sMRI coreg • Correlation-coefficient (maximise) • Slightly more general, e.g. T1-T1 inter-scanner • Invariant under affine transformation of intensities • AKA Normalised Cross-Correlation, Zero-NCC

19. Automatic image registration • Quantifying the quality of the alignment with a measure of image similarity allows computational estimation of transformation parameters • This is the basis of both realignment and coregistration in SPM • Allowing more complex geometric transformations or warps leads to more flexible spatial normalisation • Let’s look at an example of realignment next • We will return to coregistration later...

20. Option of second pass, re-registering to mean from first pass, improves results if first image is not perfectly representative of all others

22. How much motion is significant? 0.05 pixel shift (187 mm) 0.25 pixel shift (~1 mm) Edge: DS(x) ~ 70-90% 3 - 5% DS Inside: DS(x) ~ 10-20% 3 - 5% DS Slide from Kent Kiehl’s 2005 (USA) Course Workshop

23. Motion in fMRI • Even small head movements can be a major problem: • Movement artefacts add up to the residual variance and reduce sensitivity. • Data may get completely lost if sudden movements occur during a single volume • Movements may be correlated with the task performed • Minimising movements is one of the most • important factors for ensuring good data quality

24. fMRI motion prevention • Constrain the volunteer’s head • Instruct him/her explicitly to remain as calm as possible, not to talk between sessions, and swallow as little as possible • Do not scan for too long – everyone will move after while!

25. fMRI motion correction • Each volume in the time-series is realigned • Using rigid-body registration • Assumes that brain shape doesn’t change • Least squares cost function • Realigned images must be resliced for analysis • Not necessary if they will be normalised anyway • An example of severe motion…

28. Residual errors in realigned fMRI Even after realignment, as much as 90% of the variance can be accounted for by effects of movement • This can be caused by e.g.: • Movement between and within slice acquisition • Interpolation artefacts due to resampling • Non-linear distortions and drop-out due to inhomogeneity of the magnetic field Incorporate movement parametersas confounds in the statistical model

29. Movement-by-distortion interaction • Subject disrupts B0 field, rendering it inhomogeneous • distortions occur along phase-encode direction • Subject moves during EPI time series • Distortions vary with subject orientation • Shape varies

30. Co-registration example • Rigid intra-modality inter-subject registration • Purely to illustrate concepts • Coreg more commonly used to register a structural image to a mean functional one • Allows more accurate anatomical localisation • Assists spatial normalisation (see later) • Examples in SPM8 Manual, Chapters 28, 29

31. Reference AKA templateAKA fixed image (single_subj_T1) Source AKA template (!)AKA moving image(auditory/sM00223)

32. spm_defaults; global defaults opts = defaults.coreg.estimate opts = cost_fun: 'nmi' sep: [4 2] tol: [1x12 double] fwhm: [7 7] opts.cost_fun = 'ncc'; % intra-modal x = spm_coreg(ref, src, opts) % … iteration 1... 3.72 0 0 0 0 0 | -0.50731 3.72 15.95 0 0 0 0 | -0.55893 3.72 15.95 19.61 0 0 0 | -0.62815 3.72 15.95 19.61 -0.1837 0 0 | -0.6434 3.72 15.95 19.61 -0.1837 -0.00151 0 | -0.64341 3.72 15.95 19.61 -0.1837 -0.00151 -0.00063 | -0.64341 iteration 2... ... iteration 8... 2.447 20.97 28.73 -0.1905 -0.003633 -0.01026 | -0.71486 2.447 20.97 28.73 -0.1905 -0.003633 -0.01026 | -0.71486 2.447 20.97 28.74 -0.1905 -0.003633 -0.01026 | -0.71487 2.447 20.97 28.74 -0.1905 -0.003633 -0.01026 | -0.71487 2.447 20.97 28.74 -0.1905 -0.003633 -0.01026 | -0.71487 2.447 20.97 28.74 -0.1905 -0.003633 -0.01031 | -0.71487 x = 2.52 20.86 28.79 -0.20 -0.02 -0.01 M = spm_matrix(x) M = 0.999 -0.011 -0.016 2.522 0.008 0.981 -0.195 20.857 0.018 0.195 0.981 28.790 0 0 0 1.000 det(M) ans = 1.000 Svx2Rvx = inv(ref.mat)*inv(M)*src.mat

33. Inter-modal similarity measures • Often derived from joint and marginal entropies • Entropies via probabilities from histograms • Mutual Information (maximise) • MI= ab P(a,b) log2 [P(a,b)/( P(a) P(b) )] • Related to entropy: MI = H(a) + H(b) - H(a,b) • H(a) = -a P(a) log2P(a) • H(a,b) = -a P(a,b) log2P(a,b) • Normalised MI = (H(a) + H(b)) / H(a,b)

34. Joint and marginal histograms

35. Joint histogram based registration After deliberate misregistration(10mm relative x-translation) Initially registered T1 and T2 templates Joint histogram sharpness correlates with image alignmentMutual information and related measures attempt to quantify this

37. Overview of SPM Analysis Statistical Parametric Map Design matrix General Linear Model Parameter Estimates fMRI time-series MotionCorrection Smoothing SpatialNormalisation Anatomical Reference

38. Spatial Normalisation

39. Spatial Normalisation - Reasons • Inter-subject averaging • Increase sensitivity with more subjects • Fixed-effects analysis • Extrapolate findings to the population as a whole • Mixed-effects analysis • Make results from different studies comparable by aligning them to standard space • e.g. The T&T convention, using the MNI template

40. Standard spaces The Talairach Atlas The MNI/ICBM AVG152 Template The MNI template follows the convention of T&T, but doesn’t match the particular brainRecommended reading: http://imaging.mrc-cbu.cam.ac.uk/imaging/MniTalairach

41. Coordinate system sense z z y x y y z x x • Analyze™ files are stored in a left-handed system • Talairach space has the opposite (right-handed) sense • Mapping between them requires a reflection or “flip” • Affine transform with a negative determinant Rotated example

42. Spatial normalisation – Issues We want to warp the images to match functionally homologous regions from different subjects • However… • No exact match between structure and function • Different brains – different structural topology • Computational problems (local minima, etc.) Compromise by correcting gross differences followed by smoothing of normalised images

43. Spatial Normalisation – Procedure • Start with a 12 DF affineregistration • 3 translations, 3 rotations3 zooms, 3 shears • Fits overall shape and size • Refine the registration withnon-linear deformations • Algorithm simultaneously minimises • Mean-squared difference (Gaussian likelihood) • Squared distance between parameters and their expected values (regularisation with Gaussian prior)

44. Spatial Normalisation – Warping Deformations are modelled with a linear combination of non-linear basis functions

45. Spatial Normalisation – DCT basis The lowest frequencies of a 3D discrete cosine transform (DCT) provide a smooth basis spm_dctmtx(5,5) ans = 0.447 0.602 0.512 0.372 0.195 0.447 0.372 -0.195 -0.602 -0.512 0.447 0.000 -0.633 -0.000 0.633 0.447 -0.372 -0.195 0.602 -0.512 0.447 -0.601 0.512 -0.372 0.195 % Note, pinv(x)=x’, projection P=x*x’ P{n} = x(:,1:n)*x(:,1:n)’ P{N} == eye(N) P{n<N} projects to smoother approx. plot(spm_dctmtx(50, 5))

46. Spatial Normalisation– Templates and masks PET T1 Transm T2 T1 PD PD PET EPI A wider range of contrasts can be registered to a linear combination of template images. Spatial normalisation can be weighted so that non-brain voxels do not influence the result. More specific weighting masks can be used to improve normalisation of lesioned brains.

47. Spatial Normalisation – Results Affine registration Non-linear registration

48. Optimisation • Find the “best” parameters according to an “objective function” (minimised or maximised) • Objective functions can often be related to a probabilistic model, e.g. Maximum Likelihood Global optimum(most probable) Objective function Local optimum Local optimum Value of parameter

49. Optimisation – poor local optima

50. Optimisation – regularisation • The “best” parameters according to the objective function may not be realistic • In addition to similarity, regularisation terms or constraints are often needed to ensure a reasonable solution is found • Also helps avoid poor local optima • These can be considered as priors in a Bayesian framework, e.g. converting ML to MAP: • log(posterior) = log(likelihood) + log(prior) + c