Similarity functions include correlation coefficient, sum
of absolute differences, and [normalized] mutual information.
Transformations
include various types of affine (non-warping) operations. The
Optimization algorithms include Powell's algorithm, a parabolic fit to
random points, simulated annealing, or no operation (useful for
plotting). Interpolation methods include nearest neighbor,
tri-linear and a random jitter.
Similarity functions: The
correlation coefficient is defined as E{(f-m
f)*(g-m
g)}/(E{(f-m
f)
2}*E{(g-m
g)
2})
1/2,
where f and g are the pixel values two images (F and G), E{} is the
average over the pixels in the image, and m
f and m
g
are the respective mean values. The correlation coefficient is
just a cross correlation normalized for the power in the images.
The sum of absolute differences
is defined as 1-sum{|(f-m
f)-(g-m
g)|}/sum{|f-m
f|+|g-m
g|},
where sum{} is the sum over the pixels in the image.
Both the
correlation coefficient and the sum of absolute differences assume that
the two images are on the same scale. White objects in F
correspond to white objects in G, and black objects in F correspond to
black objects in G. These two similarity functions are most
appropriate for aligning images from the same modality, e.g. CT to
CT. The following example shows a CT and a low resolution
transmission map.
[Normalized] mutual information is defined as
I(F,G)/H(F,G), where I(F,G) is the mutual information of two images, F
and G, and H(F,G) is the joint entropy of the two images. The
mutual information is defined as
sum{pdf(f,g)*log[pdf(f,g)/pdf(f)*pdf(g)]}, where pdf(f,g) is a
two-dimensional histogram of pixel values in images F and G, and pdf(f)
and pdg(g) are the marginal one-dimensional histograms of the pixel
values in the respective images. The joint entropy is defined as
sum{pdf(f,g)*log[1/pdf(f,g)]}.
The advantage of mutual
information is that the two images do not need to be on the same
scale. White objects in F can correspond to back objects in G,
etc. The requirement is only that regions of a single intensity
in F correspond to regions which have a single intensity in G (or a
small number of intensities in G). The two intensity values do
not need to be the same, rather the region with a single intensity
value in F needs to correspond to a region in G which has a single
intensity value (or a small number of intensities). Mutual
information is particularly useful for multimodality comparison, e.g.
CT to MRI. Random jitter interpolation is most appropriate for
mutual information. The following example shows the first image
from mri-stack.tif and the same image which was color coded, converted
to RGB, converted to 8-bit color, inverted, and redisplayed in black
and white. Although the intensities have been altered, the same
regions are recognizable; this feature is the key feature for the
mutual information approach.
Transformation: Various
affine (non-warping) transformations are allowed. Limiting the
number of parameters which are optimized increases the probability of
success. I have experimented most with "Translate", "Scale",
and "Rotate". I have experimented much less (and much less
successfully) with combinations of operations. "Shear" has not
been adequately debugged.
Optimization algorithm:
The optimization algorithm is a procedure for changing the
transformation. Based on the value of the similarity function, a
sequence of transformations is tested in an attempt to find the
transformation which results in maximizing the similarity
function.
Powell's algorithm and simulated annealing are taken from
Numerical
Recipes.
Powell's algorithm does a sequence of 1-dimensional searches for a
maximum in the similarity function.
Simulated annealing uses a downhill simplex algorithm which has been
modified so that it will occasionally go up hill. Up hill moves
occur with a probability determined by the simulated annealing
algorithm. The idea of allowing up hill moves is to allow the
algorithm to get out of local minima. As the optimization
progresses, a "temperature" decreases, and this "temperature"
determines how likely it is for the algorithm to go up a hill based on
the ratio of the hill size to the "temperature". It jumps over
big hills initially looking for deep valleys; subsequently, it jumps
over lower hills finding the deepest point within the valley it found
during the earlier more vigorous jumping. The downhill simplex
algorithm requires optimization of at least two parameters, so this
algorithm cannot be used if only a single parameter is being optimized.
The parabolic fit algorithm is just a parabolic fit to a set of
samples at random points over the rage specified. The parabolic
fit is performed for each parameter separately. It often
provides the most accurate sub-pixel registration. It tends to be
least affected by the method used for interpolation.
"No
operation" may be used to plot without first running an optimization.
Interpolation: After
transformation, the new pixel locations will not in general line up
with the old pixel locations. Nearest neighbor interpolation
uses the value from the pixel which is nearest to the position of the
new pixel. Tri-linear interpolations uses an average of the
surround pixel values weighted by the distance of the new location to
the old locations. Random jitter uses the value of one of the
surrounding pixel values chosen at random depending upon the distance
to the original pixel locations. The random jitter is uniform
over one pixel spacing.
The similarity
functions tend to have local maxima where the matrix points are
aligned, especially tri-linear interpolation. Interpolation
using random jitter tends to de-emphasize these
maxima. Random jitter is particularly appropriate for the mutual
information similarity function. Slice 2 is redisplayed using the
selected interpolation method in order to show the effects of the
interpolation. The redisplayed image is most useful for showing
the effects of random jitter. (This interpolation is only used
for the first redisplay.)
Parameters: For different
options, various parameters need to be entered. Only certain of
the parameters are used depending upon the options selected. The
initials connect the options with the parameters. It's a bit
quicker, if somewhat more confusing to have all the options on the
initial dialog rather than having a sequence of dialogs.
If "
CC: Std for smoothing" is
non-zero when using the correlation coefficient similarity function,
then the second stacks is smoothed with a Gaussian defined by this
standard deviation. Since this is a 3-dimensional operation, it
can take time (and memory). Smoothing has been described as a
method of overcoming local maxima caused by interpolation; however, in
my hands it seems not to get rid of the local maxima, merely to make
the local maxima smaller. Making the maxima smaller doesn't
really help with Powell's algorithm, although presumably it could
help with simulated annealing.
"
PW, SA: stopping rule (pixels)"
gives a stopping rule for Powell's algorithm and simulated
annealing. When changes are less than the number of pixels given
by this parameter, the algorithm stops. In
Numerical Recipes,
stopping is defined by machine precision. Since interpolation
limits registration accuracy to something on the order of a couple of
tenths of a pixel, I have added a stopping rule in terms of pixels.
"
PW, SA: minimum overlap (%)"
gives the minimum amount of overlap between images. This minimum
and the change penalty affect the similarity function not the
optimization algorithms, so it is still possible to get to some run
away situations.
"
PF: range" and "
PF: samples" give the range and the
number of samples used in the parabolic fit optimization
algorithm. A number of samples is chosen at random over a
particular range in pixels. The values of the similarity function
are fit to a parabola and the maximum is determined.
"
SA: temperature" is the
starting temperature for the simulated annealing algorithm. This
temperature should be similar to the changes in the value of the
similarity function.
"
SA: decrease per step (0,1)"
is fraction by which the temperature is decreased during each step of
the simulated annealing algorithm. This parameter and the next
determine how fast the temperature will decrease.
"
SA: moves per step" is the
number of steps the simulated annealing algorithm makes between changes
in the temperature. It might be appropriate to increase this
parameter or decrease the former parameter when optimizing for a larger
number of transformation parameters.
Show moves: Show moves is
for debugging.
Change penalty: The
similarity function can be penalized for large changes in order to keep
the optimization form going to unreasonable solutions. The
penalty is defined as halfPenalty
2/(distance
2+halfPenalty
2),
where "halfPenalty" is the distance where the penalty is 1/2, and
"distance" is the distance in pixels of the test movement. The
penalty, which
varies from 1 to 0 as the distance becomes larger, is multiplied by the
similarity function. (The optimization parameters are scaled so
that a value of 1 corresponds to about a change in the image of 1
pixel.) The change penalty is most effective (changes the most)
near the halfPenalty point. If the algorithm manages to get to a
distance which is large with respect to halfPenalty, run away can occur.
Plot: Plotting will show
a plot of the similarity function along each of the parameters which is
being optimized. It can sometimes be helpful in identifying a
local maximum.
6. External
If Align3_TP is installed with argument "External" options appear which
allow storing and retrieving data. The purpose of this option is
to allow an external registration program to be used. For
example, the images could be roughly aligned with this plugin,
precisely aligned with an external progam, and then the alignment could
be viewed with this plugin. This option is not
further described in this manual.
7. Automatic
Align3_TP is intended to be run in an interactive mode by an imaging
expert. However,
if Align3_TP is intalled with argument "Automatic", it will run in an
automatic mode. In automatic mode, it will follow a script
contained in "plugins/Align Stacks/Align3_TP.script". Command
templates are found in "plugins/Aling
Stacks/Align3_TP_Template.script". The template file tersely
describes the script syntax, and which commands are available for
scripting. The script files need to be fairly accurate, and there
is little error checking. This option has not been extensively
debugged.
Show 3D region & Set 3D region
"Show 3D region" and "Set 3D region" set and show two dimensions of a
3D rectangular region of interest using an ImageJ ROI. They work
only when the view is axial, coronal, or sagital, e.g. lined up with
the data matrix. "Show 3D region" and "Set 3D region" are used
with volume registration and with the external registration option.