Transcript L16_SPM_Chap15 - Research Imaging Institute

Statistical Parametric Mapping
Lecture 16 - Chapter 15
Registration, brain atlases and cortical
flattening
Textbook: Functional MRI an introduction to methods, Peter Jezzard, Paul
Matthews, and Stephen Smith
Many thanks to those that share their MRI slides online
Structural Analysis
• Physical meaning of functional findings
• Registration
– image-to-image
– image-to-atlas
– combine functional and anatomical images
• 3-D Affine, warping
• Flattening of cortex
Figure 15.2 An example of the unwarping transformation to remove geometric distortion.
Acquire field map image (Jezzard and Balaban 1995). Calculate local geometric
distortion due to field changes and warp to correct. Use Jacobian (local ratio of volume
elements) to correct for signal intensity when stretching or compressing.
Figure 15.3 An illustration of automatic removal of non-brain structures using BET
(Smith 2000).
non-brain structures vary by individual
and include additional tissue types (bone,
muscle, fat)
Start with small sphere within brain
and inflate until certain boundary
conditions are met.
Figure 15.4 An example of bias field removal; original image, estimated bias field
and restored image (Zhang et al. 2001).
At higher field strengths RF penetration issues lead to RF inhomogeneity (depend on
shape of transmit and receive coils). The bias field image is calculated and multiplied
by the uncorrected image for correction.
• iterate segmentation and bias field correction until acceptable bias field
obtained (Zhang et al. 2001).
• search for bias field that maximizes some measure of ‘sharpness’ or entropy in
the image (Sled et al. 1998). Uses joint histogram (see Fig 15.13).
Figure 15.5 An illustrative example showing a translation of 20 units to the right.
Registration includes translation, rotation, scaling for 9-parameter affine transform
and adds skewing for 12-parameter transform.
Figure 15.6 Various examples of linear transformations of an original image on the left.
The order of rotations, translations, scaling and skew are important. However, two
different orderings that lead to the same 4x4 transform matrix produce the same effect.
• Rigid-body (6 degrees of freedom) - translation, rotation only
• Similarity (7 DOF) - translation, rotation, single global scaling
• Affine (9 or 12 DOF) - translation, rotation, scale, (skew)
9-Parameter Affine Transform
x' 1
  
y' 0
z' 0
  
1  0
0
1
0
0
0 Txx
 
0 Tyy
1 Tzz 
 
0 1 1
x' 1
0
0
  
y' 0 cos( x ) sin(  x )
z' 0 sin(  x ) cos( x )
  
0
0
1  0
4x4 translation
matrix T (3
parameters)
0x
 
0y
0z 
 
11
rotation matrix Rx for
rotation about x-axis
(1 parameter)
x' Sx
  
y' 0
z' 0
  
1  0
x' cos( y )0 sin(  y )
  
1
0
y'  0
z' sin(  y ) 0 cos( y )
  
0
0
1   0
0
Sy
0
0
0x
 
0y
0z 
 
11
rotation matrix Ry for
rotation about y-axis
(1 parameter)

0
0
Sz
0
0x
 
0y
0z 
 
11
4x4 scale matrix S
(3 parameters)
x' cos( z ) sin(  z )
  
y' sin(  z ) cos( z )
z'  0
0
  
0
1   0
0
0
1
0
0x
 
0y
0z 
 
11
rotation matrix Rz for
rotation about z-axis
(1 parameter)
M= [S][Rx][Ry][Rz][T]
Above ordering does translations first to match origins, then rotations about new
origin, and finally scaling of the aligned image. Most software used 9 or 12
parameters (sometimes called degrees of freedom) for registration.
12-Parameter affine transform matrix used
with 3-D images
 x'  m11 m12 m13 m14   x 
 y ' m21 m 22 m23 m24  y 
 
 
 z '   m31 m32 m33 m34   z 
  
 
0
0
1  1 
1  0
3 each rotations, scales, shears
m14, m24, and m34 are x, y, and
z translations
High DOF warping of brain image volumes
Topology preserving
means that joined
structures remain
joined after warp.
Violation leads to
folding over, holes, and
disjoint tissues.
Deformation vector
matrix matched to
image used. Local
smoothness of Jacobian
can be used to preserve
topology.
Figure 15.7 Examples of topology preserving (left) and non-preserving (right) warps applied
to the original image shown if figure 15.6. Note hole near ventricle (upper right) and splitting
of the lateral ventricles (lower right).
Higher order
interpolation may be a
2nd or 3rd order
polynomial fit through
the initial data points.
Sinc function
interpolation is
considered the most
accurate, but problem
with extent that must be
truncated.
Initial data points
linear interpolation
resampled data points
What to do at the ends?
Figure 15.8 A 1-D example of linear interpolation. Solid lines show alignment of initial
data points. Lines between closed circles are linear trend between adjacent data points.
Dotted lines show the the alignment of desired data points.
Nearest Neighbor (NN)
Trilinear
Figure 15.9 Example showing block-like artificial boundaries when NN interpolation is
used with low resolution image. Zooming of original 8x8x8 mm spacing matrix to a 1x1x1
mm display emphasizes this. Trilinear interpolation leads to a smoother and perhaps more
realistic image for this low resolution image.
Nearest Neighbor
Trilinear
Figure 15.10 Example of ‘artifacts’ (diagonal discontinuities) when NN interpolation is used
with rotation. Trilinear interpolation produces a smoother result.
• For many images trilinear interpolation is often adequate.
• NN used to conserve histogram only.
• sinc interpolation used when highest quality is critical.
Manual or Semi-Automated Registration
• Identifying landmarks
• point-like (AC and PC)
• centroid of highly penetrant features
• mid-sagittal plane
• fiducials (good for MRI, CT, PET, SPECT registration)
• Mango image processing application
• landmark based (load landmark coordinates from file)
• point matching (pick corresponding landmarks in two
images)
• atlas based (Talairach atlas method called SN which
uses key landmarks of the atlas brain)
Automated Registration
• Remove interoperator variability
• Automated landmark matching (fiducials)
• Iteratively adjust some coordinate transform using some
measure of similarity between source and target image
• similarity function/cost function as terminology
• Measuring Similarity
• intra-modality
• inter-modality
Intra-modal Registration
image pairs should map to the same or similar intensity and
have similar contrast mechanism
MeanAbsoluteDifference 
Cost Functions
(minimize)

Similarity
Function
(maximize)

MeanSquareDifference 
NormalizedCorrelation 
1
 IAi  IBi
N i
2
1
IAi  IBi 

N i


(1/N) IAi  I Ai
i


1
 IAi  I Ai IBi  I Bi
N i

2


(1/N) IBi  I Bi
i
IAi is voxel value ‘i’ in image A and IBi is the voxel value in image B, and we are
matching image A to image B.

2
Image 1
Image 2
Image 2 – Image 1
white --> +
grey --> 0
black --> -
Figure 15.11 Example showing difference image formed for large rotation (top),
small rotation (middle), and different signal contrast and levels.
Inter-modal Registration
Image pairs do not map to the same or similar intensity and/or have
different contrast mechanism. Assumes similar image regions have
similar intensities.
WoodsFunction  
Cost Function
(minimize)
k
nk  k
N k
where k2 is the variance, k the mean and nk is the number of voxels in region k

n k  k2
CorrelationRatio  1 
2
k N 
Similarity
Functions
(maximize)

MutualInformation  H(IA )  H(IB )  H(IA ,IB )
n 
n
whereH(IA ,IB )   ij log  ij 
N 
ij N
and 2 is the total variance across the image. H(IA) and H(IB) are individual image
histograms, H(IA,IB) thejoint entropy and nij is the number of voxels in each
histogram bin (i,j). Recall that entropy (H) is a positive measure of disorder.
Region 1
Region 2
Regions (white)
are used to
evaluate Woods
cost function in
image to be
transformed.
The assumption is
that if the regions
which represent
similar tissues in
the target brain
have lower
variance in the
source brain then
match is
improving.
Region 3
Region 4
Figure 15.12 Illustration of region formed to support fitting using the Woods cost function.
Here k would be four.
T2 weighted
~ perfect
alignment
T1 weighted
Joint Histogram
T2
T1
Intermediate
alignment
T2
T1
poor
alignment
T2
T1
Figure 15.13 Illustration of typical joint histograms for T1 and T2 weighted images with
differing spatial mismatch due to rotation. Entropy calculated from joint histogram.
Optimization of fitting
Figure 15.14 Illustration, using similarity measures calculated from a real image pair,
showing both local maxima and global maximum.
Optimization of fitting
1 mm
2 mm
4 mm
8 mm
Figure 15.15 Example showing image with 1x1x1 mm voxels and three subsamplings to
2 mm, 4 mm, and 8 mm sample spacings.
Multi-scale approach analysis starting at 8 mm and progressing stepwise to 1 mm helps
to avoid local maxima or minima.
-faster processing at lower resolution
-each successive step seeded by better estimate from lower resolution step
Tested Software for Registration
Automated Image Registration (AIR) - Woods et al. 1993 (Woods function)
FMRIB Linear Image Registration Tool (FLIRT) - Jenkinson and Smith 2001
(various options but default is correlation ratio)
MRITOTAL - Collins et al. 1994 (Mutual Information)
SPM - Friston, Ashburner, et al. 1995 (MSD and modality specific templates)
UMDS - Studholme et al. 1996 (Mutual Information)
Many other software solutions are available for registration so look for the one that
best suits your needs.
Fox, et al
Wash. U.
•
•
•
•
•
•
Lateral Skull X-Ray
Horizontal Grid Lines
A-P Dimension
S-I Dimension
Transform PET
Talairach Atlas 1967
– AC-PC line
– Origin
Source
Transformed
Transform
xFtransformed
DF
Compare
Ftarget
Spatial Normalization Algorithm
1.
2.
3.
4.
5.
y
Target
Set Transformed = Source.
Extract features (F's).
Compare features (DF).
Transform source to match features.
Repeat 2-4 until done.
The schematic plan for manual spatial normalization of brain images.
+y
Anterior-Right
+y
Anterior-Right
+x
Axial views of the brain before (left) and after (right) alignment to the mid-sagittal
plane (red line). The image was rotated clockwise and translated to the right in
this example. Note the marker used to identify the right side of the patient.
+x
+z
Superior-Right
+z
Superior-Right
+x
Coronal views of brain images before (left) and after (right) alignment to the
mid-sagittal plane (red line). The brain was rotated counter clockwise and
translated to the right.
+x
Mid-sagittal section views of the brain before (left) and after (right) AC-PC alignment
using a four-point fitting method. The four landmarks are the anterior-inferior margin of
the corpus callosum (blue), inferior margin of the thalamus nucleus (yellow), the
superior colliculus (green), and the apex of the cerebellum.
Axial View
+y
-x
Sagittal View
+z
+y
Axial and sagittal views of the bounding box after manual adjustment to match the bounding
limits of the cerebrum. The anterior, posterior, left, right, superior,and inferior bounds are
illustrated. Bounds do not generally fall within any one view of the brain.
Sagittal View
Axial View
AC
-x
PC
AC
pb
+y
sc
PC
In the axial view the anterior commissure (AC) and posterior commissure (PC) appear as
thin white lines connecting white matter between hemispheres. In the sagittal view the
AC is a conspicuous white, slightly elliptical structure, and the PC is at the elbow
between the pineal body (pb) and superior colliculus (sc).
Convex Hull Surface-Based Registration
15O-water
18F-FDG
Before
Talairach Convex Hull Template
After
Lancaster et al. 1999
MRI
Average Brain Templates Used for Registration
A
B
High DOF Warp
Octree Spatial Normalization
(OSN)
A. Conv. Hull global SN)
B. OSN regional SN brain
surface only
C
C. OSN regional SN
following GM, WM, and
CSF
D. Target image
D
Atlases
Figure 15.16 Example of axial, sagittal, and
coronal section MR images with points used to
specify Talairach space.
•
•
•
•
x-y-z origin at the anterior commissure (AC)
mid-sagittal plane is the y-z plane.
y-axis defined by AC-PC line
bounding box of aligned cerebrum to match that
of the 1988 Talairach atlas.
Talariach's Coordinate System
Z = +1 mm
• AC-PC line
• AC as origin
• Bounding Box
– 136 x 172 x 118 mm
• Right-handed system
Origin
(AC)
Montreal Neurological Institute (MNI) Coordinates
Using FLIRT to fit a T1W MR image to the MNI305 3-D average brain template
(template brain feature outline indicated by red lines).
Note the large rotation about the y-axis indicated in the left image.
Automated Atlas Labels Using Talairach Coordiantes
Z = +1
Level 1- Hemisphere
Level 2 - Lobar
Level 3 - Gyral
Level 4 - Tissue Type
Level 5 - Cell Type
Probabilistic Atlas
(Toga et al. 2006).
LONI Probabilistic Brain Atlas
Maximum probability labels for
various brain areas following
fitting to the ICBM452 brain
template using AIR with nonlinear warping.
Data based on brain images
from 40 healthy subjects.
Segmentation of each region
done by experts with review by
neuroanatomists.
Figure 15.17 Example images showing
the different stages of flattening in one
particular approach (Dale et al. 1999;
Fischl et al. 1999), courtesy of R.
Tootell and N. Hadjikhani. Inflation
removes the main folds of the sulci and
gyri, and flattening produces a planar
surface on which different functional
areas are shown.
Cuts are necessary to make flat maps
of the cortical surface.
Figure 15.18 Example of ‘flattened
activation’ courtesy of R. Tootell and N.
Hadjikhani. The different parts of the visual
cortex are identified using phase-encoded
simulation.