ARTICLE LINKS:
Abstract | References=20 (22) | View=20 full size inline images

Journal of = Computer Assisted=20 Tomography: Volume=20 22(1) January/February 1998 = pp 139-152

Automated Image Registration: I. General = Methods and=20 Intrasubject, Intramodality Validation

Woods, Roger P.; Grafton, Scott T.; Holmes, = Colin J.;=20 Cherry, Simon R.; Mazziotta, John C.

From the Division of Brain Mapping (R. P. = Woods, C. J.=20 Holmes, and J. C. Mazziotta) and Departments of Neurology (R. P. Woods, = C. J.=20 Holmes, and J. C. Mazziotta), Pharmacology (S. R. Cherry and J. C. = Mazziotta),=20 and Radiology (J. C. Mazziotta), UCLA School of Medicine, Los Angeles, = CA, and=20 Departments of Neurology and Nuclear Medicine, Emory University School = of=20 Medicine, (S. T. Grafton), Atlanta, GA, U.S.A.

Address correspondence and reprint requests = to Dr. R.=20 P. Woods at Department of Neurology, 710 Westwood Plaza, Rm. 3-145, Los = Angeles,=20 CA 90095, U.S.A.

Article Outline=20

Figures/Tables=20

Abstract TOP=20

Purpose: We sought to describe and = validate an=20 automated image registration method(AIR 3.0) based on matching of voxel=20 intensities.

Method: Different cost functions, = different=20 minimization methods, and various sampling, smoothing, and editing = strategies=20 were compared. Internal consistency measures were used to place limits = on=20 registration accuracy for MRI data, and absolute accuracy was measured = using a=20 brain phantom for PET data.

Results:=20 All strategies were consistent with subvoxel accuracy for intrasubject,=20 intramodality registration. Estimated accuracy of registration of = structural MRI=20 images was in the 75 to 150 =CE=BCm range. Sparse data sampling = strategies reduced=20 registration times to minutes with only modest loss of accuracy.

Conclusion: The registration algorithm = described is a=20 robust and flexible tool that can be used to address a variety of image=20 registration problems. Registration strategies can be tailored to meet = different=20 needs by optimizing tradeoffs between speed and accuracy.

In 1992 Woods et al. (1)=20 described a method for aligning PET images using a calculus-based = minimization=20 procedure and voxel intensities. This technique also proved useful for=20 intramodality registration of MR images (2)=20 and was extended to allow cross-modality registration of PET and MR = images (3).=20 A cost function based on the uniformity of the ratio of one image to the = other=20 served to guide registration through iterative univariate calculus-based = minimization. The method has compared favorably with other registration=20 techniques(4,5)=20 and has been distributed to many laboratories as part of the Automated = Image=20 Registration (AIR) package.

The AIR package has subsequently been = revised to=20 increase its speed and accuracy and to expand its ability to address a = broader=20 range of registration problems. The univariate minimization algorithm = has been=20 replaced by more robust multivariate methods; the cost function and = rigid body=20 spatial transformation models have been supplemented with alternative=20 approaches; and interpolation routines more appropriate for MR images = have been=20 added. Consequently, earlier references (1,3)=20 no longer accurately characterize the mathematics or the performance of = the=20 current version of the AIR package.

The purposes of this article are (a) to = describe the=20 mathematical basis of the registration strategy used by AIR 3.0 and its=20 relationship to several similar techniques; (b) to systematically = compare=20 different minimization procedures, smoothing strategies, editing = strategies,=20 sampling strategies, and interpolation strategies for intrasubject,=20 intramodality registration; (c) to describe a method for combining all = possible=20 pairwise registrations of a set of images to generate more accurate = registration=20 results; and (d) to describe the use of internal inconsistencies among = redundant=20 pairwise registrations to place limits on true registration accuracy in = the=20 absence of known gold standards. Intermodality registration will not be=20 addressed, and intersubject registration is described and validated = separately(6).

REGISTRATION ALGORITHM TOP=20

Figure=20 1 provides a schematic overview of the registration strategy used by = AIR=20 3.0. After optional smoothing or interpolation to cubic voxels, one = image, which=20 will be referred to as the reslice image, is resampled to match the = other image,=20 referred to as the reference image. Resampling is based on the current=20 parameters of the spatial transformation model and also requires an=20 interpolation model to compute voxel intensities. After thresholding to = exclude=20 voxels outside the head and optional editing to exclude voxels outside = the brain=20 in the reference image, a cost function reflecting the similarity of the = two=20 images is computed. For linear spatial transformation models, biases are = avoided=20 by reversing the roles of the reslice and reference image and inverting = the=20 spatial transformation to compute a second estimate of the cost = function, which=20 is then averaged with the first. This bias elimination procedure can be=20 optionally omitted in AIR but was always used here. To improve speed, = the cost=20 function is initially computed for only a limited sampling of the voxels = (the=20 default is every 81st voxel) and sampling is increased with subsequent=20 iterations (the default is by factors of three to reach a final sampling = of=20 every voxel). The derivatives of the cost function with respect to the=20 parameters of the spatial transformation model are computed and are used = to=20 compute new parameters and iteratively minimize the cost function. = Termination=20 criteria are tested with each iteration to decide whether to continue = iterating,=20 to increase sampling, or to stop. The spatial transformation that = produced the=20 lowest value of the cost function is stored and can be used to produce=20 registered images. If desired, the parameters that produce the optimal=20 transformation can also be stored independently. Substantial = modifications as=20 compared with the original AIR algorithm (1)=20 are described in the following sections.

FIG. = 1.=20 Schematic diagram of the registration algorithm. Boxes shown in = dashed=20 lines represent optional procedures. One of the images is = initially=20 designated as the reference study and the other as the reslice = study.=20 Although not shown in the diagram, the cost function is computed a = second=20 time with the roles of reference and reslice study interchanged, = and it is=20 the average of these two cost function estimates that is = minimized.=20 Regions of the schematic images shown in solid black are excluded = from=20 analysis (a) because they are below the specified threshold, (b) = because=20 they are excluded by the optional editing mask, (c) because they = are not=20 part of the current sampling set, or (d) because they are outside = the=20 field of view of the reslice study. These exclusions actually = occur before=20 or during interpolation of the reslice study, but this is not = indicated in=20 the diagram to preserve conceptual clarity of the other steps. = When the=20 sampling density is increased, the spatial transformation = parameters that=20 gave the optimal value for the cost function at the prior sampling = density=20 serve as the new starting point for=20 minimization.

Optional Image Smoothing TOP=20

Smoothing noisy PET images before = registration=20 improves accuracy(1).=20 Whether to smooth MRI data will be addressed in Validation Studies. An = optional=20 fast Fourier, Gaussian convolution routine replaces the box smoothing = strategy=20 used previously. The width of the smoothing kernel is specified = independently=20 for each image axis. If applied, smoothing is performed before = resampling,=20 thresholding, or application of any editing masks.

Optional Interpolation of Reference = Study to=20 Cubic Voxels TOP=20

In the original AIR algorithm, the = reference volume=20 always consisted of voxels that had been linearly interpolated to make = them=20 cubic. This is now an option, but not the default. Instead, anisotropic = voxel=20 sizes are taken into account mathematically while still ensuring the = integrity=20 of the selected spatial transformation model in real world coordinates. = The=20 tradeoffs between speed and accuracy as a function of whether the = reference=20 volume is interpolated will be addressed by the PET validation=20 studies.

Interpolation Model TOP=20

To compare the images being registered, = one image must=20 be resampled according to the parameters of the spatial transformation = model.=20 This requires interpolation of intensities at locations between the = voxel=20 locations represented in the original image. After registration is = complete,=20 final images must be created, which again requires interpolation to = compute=20 resampled voxel intensities. It is possible, but not necessary, to use = the same=20 interpolation model in both contexts. In the distribution version of AIR = 3.0,=20 trilinear interpolation remains the only model used to compute the cost = function=20 during minimization, whereas any of several models including nearest = neighbor,=20 trilinear, sinc, and chirp-z interpolation(7)=20 can be used to create the final images. Windowing(8)=20 and scan line decomposition (9)[with=20 precautions to avoid aliasing (10)]=20 are also available in the final resampling algorithm to improve speed. = Hybrid=20 models are also included for multislice data sets that are not = bandlimited along=20 the axes between planes (8).=20 A special implementation of AIR that uses windowed sinc interpolation to = compute=20 the cost function is described in Validation Studies.

Thresholding, Optional Editing, and = Bias=20 Elimination TOP=20

Simple thresholding can be used to exclude = voxels from=20 outside the body that provide no useful spatial information, but it may = also be=20 advantageous to exclude voxels from the scalp, skull, and dura when = computing=20 the cost function for MRI data (11).=20 Exclusion of these structures could be accomplished by simply = registering edited=20 images, but tendencies to align the artificial edges created by the = editing=20 process could be problematic. This can be avoided by allowing the = algorithm to=20 apply user-generated mask files to the images. When an image is serving = as the=20 reference image, only voxels that are non-zero in the associated mask = file=20 contribute to the cost function, and the mask associated with the other = image is=20 ignored. The second mask file is used when the roles of the reference = and=20 reslice files are exchanged to compute an unbiased cost function. Since = edited=20 versions of both images are never compared directly with one another, = alignment=20 of artificial edges created by editing cannot lower the cost function. = The=20 effect of editing nonbrain structures on registration accuracy will be = addressed=20 in Validation Studies.

Cost Functions TOP=20

The cost function provides the algorithm = with a=20 quantitative measure of how well the images are registered. AIR 3.0 = allows a=20 choice of three different cost functions. The first cost function, = referred to=20 here as the ratio image uniformity (RIU) cost function, is identical to = the one=20 described previously(1).=20 A resampled image is divided by the image to which it is being = registered on a=20 voxel-by-voxel basis to create a ratio image, and the uniformity of this = ratio=20 image is measured by computing its standard deviation. The standard = deviation is=20 then divided by the mean ratio to provide a normalized cost function = value.=20 Minimization of the cost function increases the uniformity of the ratio = image,=20 which is independent of global intensity scaling of the original images, = and=20 improves registration.

The second cost function assumes that the = images being=20 registered have already been properly adjusted for global intensity = differences=20 and uses a least-squares approach similar to that described by Hajnal et = al.(8),=20 and subsequently adopted by Friston et al.(12).=20 This cost function will be referred to as the least-squared difference = image=20 (LS) cost function. If no intensity rescaling is needed, the spatially = resampled=20 reslice image should be almost identical to the reference image when the = images=20 are well registered. The difference between the resampled reslice image = and the=20 reference image is computed at each voxel, and the square of this = difference is=20 averaged across voxels to generate this cost function.

The third cost function is similar to the = LS cost=20 function but adds an extra parameter to the minimization that allows = global=20 intensity rescaling of the images relative to one another. This approach = has=20 been advocated by Alpert et al. (13)=20 and by Snyder (14)=20 and will be referred to as the scaled least-squared difference image = (SLS) cost=20 function.

Minimization Procedure and Spatial=20 Transformation Model TOP=20

The objective of the minimization = procedure is to=20 conduct an efficient search of the parameter space defined by the = mathematical=20 spatial transformation model to minimize the cost function. By default, = AIR=20 starts this search with parameters that will align the exact centers of = the two=20 image sets without additional real world rotation, translation, or = scaling.=20 Alternatively, any set of starting parameters can be explicitly = declared. All of=20 the work described herein uses a rigid body spatial transformation = model, which=20 has six parameters, so the minimization procedure must search a 6D = parameter=20 space to find the optimal rigid body registration. The original AIR = algorithm=20 conducted this multidimensional search by sequentially performing = undimensional=20 minimizations, each of which corresponded to only one of the six = parameters. The=20 AIR 3.0 algorithm replaces this strategy with two variants of a = multivariate=20 calculus-based minimization procedure. The first variant, which will be = referred=20 to as full Newton-type minimization, is similar to one briefly described = for 12=20 parameter intersubject registration(15)=20 but has been generalized so that it is applicable to any spatial = transformation=20 model. The method assumes that the cost function(as a function of the = spatial=20 transformation parameters) can be approximated near its minimum by a = parabolic=20 surface. This parabolic surface can be fully characterized by a vector = b, consisting of the first partial = derivatives of the=20 cost function with respect to each parameter at a given point in = parameter=20 space, and a Hessian matrix A, = consisting of=20 the second partial derivatives of the cost function with respect to each = pair of=20 parameters at the same point in parameter space. If x is a vector representing the adjustment of = each=20 spatial transformation parameter needed to reach the minimum of the = parabolic=20 surface, x can be derived fromb and A by = solving the=20 simple matrix equation(16)=20

A*x =3D -b =

The first and second derivatives of the = cost function=20 with respect to the transformation parameters are all calculated = analytically=20 and are a function of the spatial transformation model, the = interpolation model,=20 and the particular cost function. New estimates of the minimum are = calculated=20 iteratively until the termination criteria described in the next section = are=20 met. In the event that the Hessian matrix A is=20 not positive definite(indicating a saddle point or tendency toward a = maximum=20 rather than a minimum), the algorithm either increases the sampling = density or,=20 if sampling is already maximal, generates a warning message and = terminates.

The second minimization procedure is an = approximation=20 to the first one. It assumes that the second derivatives of the = interpolated=20 voxel values with respect to the sampling coordinate locations are all = zero.=20 This decreases the computation time per iteration and reduces the = likelihood of=20 encountering problematic local maxima or saddle points. Such = approximations are=20 commonplace in minimization algorithms. For example, the = Levenberg-Marquardt=20 algorithm(16)=20 makes this same assumption and also assumes that the second derivatives = of the=20 spatial transformation model with respect to its parameters are all = zero. This=20 second minimization method will be distinguished from full Newton-type=20 minimization by appending LM after the cost function(i.e., RIU-LM, = LS-LM, and=20 SLS-LM) and for the sake of brevity will be referred to as = Marquardt-like=20 minimization. Full Newton-type minimization should be assumed as the=20 default.

Termination Criteria TOP=20

If the associated parabolic approximations = were exact,=20 the full Newton-type minimization method would proceed to the minimum of = the=20 cost function in a single iteration. However, they are not exact, and = multiple=20 iterations are required. Criteria are necessary to decide when iteration = should=20 terminate. For the original AIR algorithm (1),=20 the criteria were based on the fact that the first derivative of the = cost=20 function with respect to each parameter should be very close to zero = near the=20 minimum. The primary termination criterion of the new algorithm takes = advantage=20 of the fact that the derivatives used to compute the location of the = minimum of=20 a multidimensional parabolic surface can also be used to predict the = change in=20 the cost function associated with moving from the current point in = parameter=20 space to the predicted minimum. The predicted change (=CE=94c) is given by the matrix equation. Equation=20

Equation = 4C

where b, = x, and=20 A are defined as before(16,=20 p 414). So long as A is positive = definite, the=20 predicted change will be negative and should get progressively closer to = zero as=20 the true minimum of an arbitrary differentiable surface is approached. = Thus, the=20 predicted cost function change provides a single numerical value that=20 simultaneously incorporates information about all of the parameters in = units=20 that are independent of the particular spatial transformation model. = Termination=20 occurs when the predicted cost function change drops below a = prespecified small=20 value.

To provide finer control, the new = algorithm also has=20 secondary termination criteria. The two secondary termination criteria = are (a)=20 the total number of iterations performed at a given sampling interval = and (b)=20 the number of iterations performed without improvement in the actual = cost=20 function at a given sampling interval. In another modification of the = original=20 methodology, the new algorithm always retains a copy of the parameters = that=20 resulted in the lowest actual value of the cost function at the current = sampling=20 density and uses this best value, rather than the most recent value, = when=20 initializing the next sampling density or saving the final results.

Selection of appropriate primary = termination criteria=20 will be addressed in Validation Studies. For secondary termination = criteria, a=20 default value of 25 total iterations or 5 iterations without any = improvement in=20 the actual cost function should be assumed unless otherwise = specified.

Implementation and Distribution TOP=20

AIR 3.0 is written entirely in C and = requires no=20 proprietary code or third party packages. The registration algorithm can = be=20 compiled to use either 8 or 16 bit internal image representation. The=20 benchmarking work described here for PET images was done on a SunSPARC = 10=20 workstation, compiled for 8 bit internal representation using the = standard Sun C=20 compiler. The MRI validation work was performed on a Power Macintosh = 8500/120=20 running the MachTen UNIX environment, compiled for 16 bit images. Source = code=20 for the algorithm and supporting software are available to the research=20 community for research purposes free of charge. Documentation and = information=20 about downloading the software are available on the World Wide Web at = the=20 Universal Resource Locator=20 (URL):http://bishopw.loni.ucla.edu/AIR3/index.html

VALIDATION STUDIES TOP=20

Two separate sets of validation studies = will be=20 described. The first uses a previously described PET phantom data set (1)=20 that allows registration accuracy to be established with independent = gold=20 standards. The second uses high resolution MR scans of a single subject = to=20 evaluate intrasubject MR registration. As often happens with real data = sets,=20 gold standards are not available for the MRI validation study, and = special=20 attention is given in both validation studies to the use of internal=20 inconsistency measures as an alternative approach to validation in this = setting.=20 Before describing the individual validation studies, some terms and = methods=20 common to both studies are described.

Discrepancies, Errors, and Internal = Inconsistencies TOP=20

Several terms will be used repeatedly and = warrant=20 explicit definition. The local discrepancy between two different = transformations=20 for registering a pair of images will be defined as the 3D distance = between the=20 locations to which a given voxel is mapped by the two transformations. = The mean=20 discrepancy between two different transformations for registering a pair = of=20 images will be defined as the average of the local discrepancies over = all voxels=20 that constitute those parts of the brain that are present in both = images. The=20 corresponding plural term mean discrepancies will refer collectively to = a set=20 that specifically includes the mean discrepancy of every possible unique = pairwise registration of a group of images of the same object. The term = global=20 mean discrepancy will refer to the average mean discrepancy across such = a=20 collective set. The terms maximum discrepancy, maximum discrepancies, = and global=20 maximum discrepancy will be defined analogously using the appropriate = maximum,=20 rather than the mean, across voxels or registrations. The unqualified = term=20 discrepancies will refer collectively to mean discrepancies and maximum=20 discrepancies, and the term global discrepancies collectively to the = global mean=20 discrepancy and the global maximum discrepancy.

If, and only if, one of the two = transformations (or=20 sets of transformations) being compared is based on independent gold = standard=20 measurements, the terms error or errors will analogously replace the = words=20 discrepancy or discrepancies in all the above definitions. If, and only = if, a=20 set of transformations is being compared with reconciled mean=20 transformations(defined in the next section) derived from that same set, = the=20 terms internal inconsistency or internal inconsistencies will = analogously=20 replace the word discrepancy or discrepancies in all the above = definitions.=20 Discrepancies will always be characterized by indicating the two = distinct=20 registration strategies to which the discrepancy applies, whereas errors = and=20 internal inconsistencies are properties of a single registration=20 strategy.

Derivation of Reconciled Mean = Transformations=20 from Pairwise Registrations TOP=20

If a perfect registration method were used = to perform=20 all possible pairwise registrations of a rigid body, the resulting=20 transformations would be completely internally consistent. For example, = the=20 direct pairwise registration of image A to image C would be identical to = the=20 combined results for registering image A to image B and for registering = image B=20 to image C. In this sense, a perfect set of all pairwise registrations = should=20 contain redundant information about the interrelationships between the = images.=20 With imperfect registration, direct pairwise registration of image A to = image C=20 will not be identical to the combined results for registering image A to = B and=20 image B to C due to errors. Conceptually, each imperfect pairwise result = can be=20 viewed as the combination of the correct result and an error associated = with=20 that particular pair. Better estimates of the correct results should be=20 attainable by defining the minimal set of nonredundant parameters = sufficient to=20 derive all pairwise registrations and then adjusting this set to somehow = minimize local discrepancies between the predicted pairwise results and = those=20 that are actually observed. For N images, N - 1 six parameter rigid body transformations = suffice to=20 derive all possible pairwise registrations. The completely internally = consistent=20 set of transformations defined by these six (N - 1)=20 parameters will be referred to here as the reconciled mean = transformations.

AIR 3.0 includes an algorithm that = computes reconciled=20 mean transformations from an existing set of all possible pairwise=20 registrations. The algorithm initializes the six (N -=20 1) parameters by assuming that the images are all already perfectly = registered=20 and refines these initial estimates iteratively. From a computational=20 stand-point, it is advantageous to minimize the summed squares of the = distances=20 defined by each local discrepancy between the observed and predicted=20 transformations rather than the sum of the distances themselves. Squared = local=20 discrepancies are computed for all voxels that correspond to parts of = the brain=20 that are represented in all of the images, and these squared values are = summed=20 across voxels and across all possible pairwise registrations to generate = a=20 global summed squared discrepancy. Appropriate brain voxels need only be = explicitly identified once in a single image since their locations can = be=20 remapped to other images using the current estimates of the reconciled = mean=20 transformations. Full Newton-type minimization is used to iteratively = minimize=20 the global summed squared discrepancy by adjusting the six (N - 1) parameters until the predicted change for = the next=20 iteration is <10^-10. The algorithm = reports the=20 final optimized global discrepancy, normalized for the number of image = pairs and=20 the number of landmark locations, and writes out the reconciled mean=20 transformations.

Statistical Comparisons TOP=20

Distributions of mean (or maximum) errors = or internal=20 inconsistencies among all possible pairwise registrations generated by = two=20 different registration strategies will be compared with one another = using the=20 two sided two sample Kolmogorov-Smirnov test implemented in the = statistical=20 package, S-Plus(MathSoft, Seattle, WA, U.S.A.). This test evaluates = whether two=20 samples are drawn from the same distribution by comparing their = empirical=20 cumulative distribution functions (16)=20 and does not presuppose any particular shape for the underlying = distribution.=20 All unqualified or implicit references to significance testing will = pertain=20 specifically to Kolmogorov-Smirnov tests applied to mean or maximum = errors or=20 internal inconsistencies. Results reported as significant will indicate = a p=20 value of<0.05. All results to be reported as significant have been = verified=20 graphically to involve a consistent shift in one of the empiric = cumulative=20 distribution functions rather than an isolated change in shape of one of = the=20 empiric cumulative distribution functions.

PET Phantom Validation TOP=20

AIR 3.0 was validated for intramodality, = intrasubject=20 registration of PET data using the same brain phantom data set used to = validate=20 the original AIR algorithm (1).=20 It consists of 31 PET scans of the Hoffman brain phantom (17)=20 acquired at various known scanner gantry and bed positions. Image pairs = with=20 rotational misalignment of>30=C2=B0 and translational misalignment up = to 10 mm are=20 included. The amount of radioisotope in the phantom was selected to = simulate the=20 poor counting statistics of 2D PET H₂ 15O studies. Inflatable balloons containing = higher levels=20 of activity were deflated partway through the data acquisition to = simulate focal=20 changes in activity. The data set was acquired in 2D mode using a = Siemens/CTI=20 831-08 tomograph (Siemens, Hoffman Estates, IL, U.S.A.). Images were=20 reconstructed with a Shepp reconstruction filter with a roll-off = frequency of=20 0.16 mm^-1 to generate images with a full = width at=20 half-maximum (FWHM) in-plane resolution of 6.1 mm. Voxels in the = reconstructed=20 images were 1.745 =C3=97 1.745 =C3=97 6.75 mm, and the image matrix = dimensions were 128 =C3=97=20 128 =C3=97 15 planes. All images were smoothed with a 2D isotropic = Gaussian filter to=20 an in-plane resolution of 10 mm. During smoothing, images were scaled to = map the=20 hottest voxel in an image to the highest representable 8 bit value.

The goals of the phantom validation = studies were (a)=20 to compare AIR 3.0 with the original version of AIR: (b) to compare the=20 different cost functions implemented in AIR 3.0 with one another; (c) to = determine whether interpolation of the reference image to cubic voxels = improves=20 registration accuracy in data sets with highly anisotropic voxel sizes; = (d) to=20 determine whether sparse sampling of the data can be used to increase=20 registration speed without adversely influencing registration accuracy; = (e) to=20 test the hypothesis that reconciling all possible pairwise registrations = improves registration accuracy; and (f) to evaluate internal = inconsistencies as=20 a metric for registration accuracy.

To address these issues, all 465 unique = pairs of the=20 31 phantom images were registered using several different strategies = (see Table=20 1). These strategies also included the original AIR algorithm. In = all cases,=20 an 8 bit threshold value of 55 served to segment the reference image = into brain=20 and nonbrain values. Gold standard registration parameters were computed = for=20 each pair of images based on the known scanner gantry and bed positions: = All 31=20 of the original data sets were resampled into a single common space = using the=20 gold standard transformation parameters and averaged(with weightings = necessary=20 to adjust for missing data from outside the field of view) to produce a = single=20 mean phantom image set. This image set was manually edited to remove all = peripheral nonbrain voxels. The resulting brain mask was projected back = onto the=20 original 31 data sets using the gold standard registration parameters = and used=20 to distinguish brain from nonbrain voxels when measuring errors and = internal=20 inconsistencies. For each registration strategy, mean and maximum errors = were=20 computed. The mean registration times were recorded for each strategy.=20 Reconciled mean transformations and mean internal inconsistencies were = derived,=20 and the mean errors of the reconciled mean transformations were computed = for=20 each strategy.

TABLE = 1. Average registration times, global errors, and = global=20 internal inconsistencies for PET phantom validation=20

To evaluate whether the focal activity = changes=20 simulated by inflating or deflating the balloons within the phantom = affected=20 registration accuracy, transformations derived using a given strategy = were=20 subdivided into those in which the balloons were in the same state in = both scans=20 (either both inflated or both deflated) and those in which the balloons = were in=20 different states. The error distributions were then compared across = these two=20 subgroups using Kolmogorov-Smirnov tests. This was repeated for each=20 registration strategy independently.

Intrasubject MRI Validation TOP=20

A single normal subject was scanned eight = times=20 consecutively on a Phillips 1.5 T MR scanner. Sagittal volumes of 140 1 = mm=20 slices were acquired with a field of view of 256 =C3=97 204 mm. A 3D = spoiled GRASS=20 sequence [TR/TE =3D 18/10 ms, flip angle 30=C2=B0, NSA (NEX) 1, flow = compensation] was=20 used with a total scan time per volume of 10 min 50 s. The subject (one = of the=20 authors) was highly motivated to prevent any head movements during the = course of=20 any given acquisition. Head packing was used to help prevent movements = during=20 acquisition, but no precautions were taken to prevent small head = movements=20 between acquisitions. The data were stored at 12 bit precision in 16 bit = format=20 with quantitative preservation of absolute scaling across studies. A set = of 160=20 planes covering the brain from the bottom of the cerebellum to the top = of the=20 brain was selected in the first data set. The same planes (with respect = to the=20 scanner, not with respect to the anatomy) were selected in the other = seven data=20 sets, so any apparent movement between scans accurately reflects true = movement=20 by the subject. The images were similarly reduced along the right-left = axis of=20 the brain so that only the ears were visible in the most extreme = sagittal planes=20 of the first data set. The resulting data sets consisted of isotropic 1 = mm=20 voxels with dimensions of 256 =C3=97 160 =C3=97 160 voxels. All 28 = possible unique=20 pairwise registrations of the data were performed using several = different=20 registration strategies. The goals were (a) to compare the speed, = internal=20 consistencies, and results of the three cost functions; (b) to = investigate the=20 effects of smoothing on registration; (c) to investigate the effects of = editing=20 nonbrain structures on registration; (d) to determine the tradeoffs = associated=20 with different convergence criteria, minimization procedures, and data = sampling=20 strategies; and (e) to use internal inconsistencies to establish limits = for true=20 registration accuracy. The effects of cost function, smoothing, and = editing of=20 the images were investigated first. Isotropic, 3D smoothing filters with = a FWHM=20 of 0, 2.0, 4.0, and 8.0 mm were applied. The two images being registered = were=20 always filtered identically. To ensure consistency, manual editing to = remove=20 nonbrain structures was performed only once on an average of the eight = unedited,=20 unsmoothed images after co-registration with the LS cost function. The = results=20 of manual editing were then transformed back into the eight native image = spaces=20 (using the inverse of the transformations used for co-registration) in = the form=20 of masks. The registration algorithm masked each unedited image = internally to=20 avoid any tendency to align artificial edges created by the editing = process. The=20 default sampling intervals, convergence criteria, and minimization = procedure=20 were used. A 16 bit voxel value of 480 served to separate the head = (brain and=20 surrounding tissues) from background. For each registration strategy, = the mean=20 and maximum internal inconsistencies were measured after computing the=20 corresponding set of reconciled mean transformations. Mean and maximum=20 discrepancies between selected pairs of strategies were also = computed.

The effects of different primary = convergence=20 thresholds, sampling densities, and the two different minimization = procedures=20 for each of the cost functions were investigated using unedited and = unsmoothed=20 data. A modified version AIR was also used to investigate whether = computation of=20 the LS cost function using windowed sinc interpolation would improve the = internal inconsistency of the results. A cosine half-bell windowing = technique=20 identical to that described by Hajnal et al. (8)=20 was used with full 3D sinc interpolation. Because sinc interpolation is=20 extremely slow, the sinc version of the algorithm was always initialized = with=20 optimum parameters obtained using trilinear = interpolation.

RESULTS TOP=20

PET Phantom TOP=20

Table=20 1 shows the average registration times, global mean errors, and = global=20 maximum errors for the various PET registration strategies. It also = shows the=20 global internal inconsistencies associated with each registration = strategy and=20 the global mean errors of the reconciled mean transformations. All = methods=20 resulted in subvoxel accuracy with global maximum errors of <2.0 = mm.

As compared with using an uninterpolated = reference=20 image, interpolation of the reference image to cubic voxels universally = resulted=20 in significantly smaller errors and internal inconsistencies, = independent of all=20 other factors. This was achieved at the expense of registration time, = which=20 increased by a factor roughly proportional to the increase in the number = of=20 voxels through interpolation. Because of this consistent and = significantly=20 poorer performance, registrations without interpolation to cubic voxels = were=20 excluded from the remainder of the analyses reported here. A convergence = threshold of 10^-4 for the RIU cost function = gave=20 significantly larger mean errors, maximum errors, and internal = inconsistencies=20 than more stringent thresholds of 10^-5 or = 10-6. Similarly, a convergence threshold of 10-1 for the SLS cost function resulted in = significantly=20 higher internal inconsistencies than more stringent thresholds of 10^-2 or 10^-6. These = least stringent=20 thresholds only modestly reduced registration times. The most stringent=20 thresholds resulted in longer registration times with no significant = improvement=20 in accuracy compared with intermediate thresholds.

The LS cost function produced = significantly larger=20 mean errors than the RIU cost function. No other significant differences = between=20 the RIU, SLS, and LS cost functions were identified when using stringent = or=20 moderate convergence thresholds. The SLS cost function tended to be = faster than=20 the RIU cost function. The original AIR algorithm (1)=20 produced a global mean error of 0.323 mm and a global maximum error of = 1.568=20 with a mean registration time of 180 s. The errors are similar to those = produced=20 by the new algorithm using the RIU cost function.

Increasing the final sampling interval = from every=20 voxel to every 27th voxel for the RIU and SLS cost function did not=20 significantly increase errors or internal inconsistencies and reduced=20 computation times by a factor of 3-5. Speed improvements without = significant=20 loss of accuracy were also achieved by using Marquardt-like rather than = full=20 Newton-type minimization. The reconciled mean transformations invariably = resulted in significantly smaller errors than the original = transformations.=20 Internal inconsistencies, reflecting deviations from these reconciled = mean=20 transformations, were always smaller than the errors based on the gold=20 standards. Larger errors were correlated with larger internal = inconsistencies=20 for any individual cost function, but the internal inconsistency = measures were=20 insensitive to the differences between cost functions.

The sites of focal differences in activity = resulting=20 from inflation or deflation of the balloons had no significant effect on = accuracy. The distributions of errors was invariably statistically=20 indistinguishable for registrations in which the balloons were in the = same state=20 as for registrations in which the balloons were in different = states.

Intrasubject MRI TOP=20

Visual inspection of the resampled MR = images from the=20 single subject confirmed that the registrations were qualitatively = correct.=20 Images of a transverse brain section through the anterior commissure = before and=20 after registration are shown in Fig.=20 2. The estimated maximum initial misalignment in the brain between = pairs of=20 images ranged from 0.26 to 3.80 mm. The estimated average initial = misalignment=20 across all brain regions and image pairs was 0.76 mm. Based on results = from the=20 SLS cost function, which were confirmed by direct inspection, the global = intensity scaling varied to a small extent from one image to the next = with the=20 largest pairwise differences being =E2=89=886%.

FIG. = 2. A=20 single transverse section from four of the eight image sets used = for MRI=20 validation. The first row shows the section before registration. = Subtle=20 differences due to misregistration can be seen and are highlighted = by=20 black arrows. The second row shows the same section after = registration of=20 the corresponding volumes using the LS cost function with sparse = sampling=20 at every 81st voxel, requiring <70 s per registration. Data = resampling=20 was performed using windowed sinc interpolation. All registration=20 strategies investigated resulted in images virtually = indistinguishable=20 from those on the second row.

Table=20 2 shows the average registration time and internal inconsistencies = of the=20 three different cost functions with and without smoothing and editing of = the=20 data.Tables=20 3-5 show discrepancies between results obtained when cost function,=20 smoothing, or editing was varied while keeping all other factors = constant.

TABLE = 2. Effects of cost function, smoothing filter, and = data=20 editing on MRI registration time and global internal = inconsistency=20

TABLE = 3. Global mean and maximum discrepancies between = MRI=20 registration results obtained with the three different cost = functions as a=20 function of smoothing filter and editing=20

TABLE = 4. Global mean and maximum discrepancies between = MRI=20 registration results obtained with different smoothing filters as = a=20 function of registration method (cost function and editing)=20

TABLE = 5. Global mean and maximum discrepancies between = results=20 obtained from edited versus unedited MRI data as a function of=20 registration method and Gaussian smoothing filter=20

None of the three cost functions produced=20 significantly better internal inconsistencies than either of the others. = The SLS=20 cost function was fastest and produced results within 32 =CE=BCm of the = LS cost=20 function. The RIU cost function was unexpectedly faster than the LS cost = function for unsmoothed images and gave results that always agreed with = the=20 other two cost functions to within 123 =CE=BCm. The mean discrepancy = between cost=20 functions was always=E2=89=A425 =CE=BCm for a given smoothing and = editing strategy.

Data smoothed with a 2 mm Gaussian filter = invariably=20 produced significantly smaller mean and maximum internal inconsistencies = than=20 unsmoothed data. This was true even for the RIU cost function, which = required=20 initialization with the optimal registration parameters obtained from = the=20 corresponding unsmoothed RIU registration when using smoothed data. = Without such=20 initialization, the RIU (and likewise the RIU-LM) cost function made = excessively=20 large steps at the first iteration and failed to find any way to improve = the=20 default initialization parameters. For unedited data (but not for edited = data),=20 even heavier smoothing with a 4 mm Gaussian filter produced additional=20 significant improvements in internal inconsistencies. However, when = smoothing=20 with an 8 mm Gaussian filter was investigated for the LS and SLS cost = functions,=20 internal inconsistencies were always worse than with 4 mm smoothing but = still=20 superior to those with unsmoothed data. Smoothing increased registration = time.Table=20 4 shows that smoothing altered registration results substantially. = Results=20 from unsmoothed images differed from those obtained with smoothed images = and=20 otherwise identical strategies by distances as large as 400 =CE=BCm, and = even global=20 mean discrepancies were at least 98 =CE=BCm.

To verify that these differences between = the smoothed=20 and unsmoothed results were not simply due to local minima, the = unsmoothed,=20 unedited registrations were repeated, using the corresponding optimal=20 registration parameters obtained with 4 mm Gaussian smoothing as the = starting=20 parameters for minimization. Likewise, the unedited registrations using = 4 mm=20 Gaussian smoothing were repeated, using the optimal parameters obtained = with=20 unsmoothed data as starting parameters. To eliminate sparse sampling as = a=20 potential confound, both sets of repeat registrations utilized full = sampling at=20 every voxel from the first iteration. Despite being initialized with a = different=20 set of parameters, the final results were unchanged. For all three cost=20 functions, global mean discrepancies between the original and repeat=20 registrations never exceeded 10 =CE=BCm.

The effects of editing the data to exclude = scalp,=20 skull, and dura varied depending on the smoothing strategy. In the = absence of=20 smoothing, editing led to significantly worse internal inconsistencies = for the=20 SLS and LS cost functions. In the presence of 2 mm Gaussian smoothing, = editing=20 led to significantly better internal inconsistencies for the SLS and RIU = cost=20 functions, and in the presence of 4 mm Gaussian smoothing, editing had = no=20 significant effect. Editing always reduced registration times, but by a = factor=20 of <2. Editing also substantially changed the registration results = with=20 global maximum discrepancies between edited and unedited results inTable=20 5 as large as 406 =CE=BCm and global mean discrepancies always = larger than 74=20 =CE=BCm.

Table=20 6 shows the effects of sampling intervals, convergence thresholds, = and=20 minimization procedures for each cost function using unsmoothed, = unedited data.=20 Baseline measurements for each cost function were performed using full=20 Newton-type minimization with a primary convergence threshold of zero, = forcing=20 as many iterations as the secondary termination criteria allow. = Different=20 sampling intervals and minimization strategies were then compared = against the=20 corresponding baseline values to compute the discrepancies shown in the = table.=20 The RIU convergence threshold of 10^-5 is = perhaps not=20 optimally stringent for MRI data since the results with this threshold = differed=20 by as much as 83 =CE=BCm from those obtained with a convergence = threshold of zero.=20 The SLS and LS thresholds of 1.0 gave results within 10 =CE=BCm of those = obtained=20 with a threshold of zero. For all three cost functions, sparser sampling = was=20 associated with a modest increase in global discrepancies from the = baseline, an=20 increase in global internal inconsistencies, and a substantial reduction = in=20 registration time. Marquardt-like minimization was also associated with = shorter=20 registration times, giving results almost identical to full Newton-type=20 minimization for the LS and SLS cost functions, but more divergent = results for=20 the RIU cost functions. Both minimization strategies proved several = times faster=20 than the strategy used in the original version of AIR (data not = shown).

TABLE = 6. The effects of minimization method, convergence = threshold,=20 and final sampling interval on MRI registration=20

With use of sinc interpolation to compute = the LS cost=20 function with a very stringent convergence threshold of 10-5 and a final sampling interval of every voxel = with=20 unedited, unsmoothed data resulted in a global mean internal = inconsistency of 20=20 =CE=BCm and a global maximum internal inconsistency of 63 =CE=BCm. The = distributions of=20 internal inconsistencies were not significantly better than with = comparable use=20 of trilinear interpolation and a less stringent convergence = threshold.

Within the ranges investigated, the = sampling interval=20 and minimization method made the least difference, the primary = convergence=20 threshold, the cost function, and the interpolation model had an = intermediate=20 effect, and smoothing and editing had the greatest impact on the final = results.=20 All methods were in general agreement with global maximum discrepancies=20 always<500 =CE=BCm. The results are compatible with, but not = conclusive proof of,=20 typical registration accuracies on the order of 75-150 =CE=BCm. It is = possible that=20 one of the strategies included here may be more accurate than the = others,=20 possibly with typical accuracies <10 =CE=BCm, but this clearly cannot = be the case=20 for all of the strategies simultaneously.

DISCUSSION TOP=20

The AIR 3.0 registration algorithm is = robust, fast,=20 accurate, and applicable to a diverse range of registration problems. = Except for=20 the difficulties with registration of smoothed MRI data using the RIU = cost=20 function, no registration failures were identified, and the algorithm=20 consistently achieved subvoxel accuracy. Reconciling internal = inconsistencies=20 among all possible pairwise registrations further improved registration=20 accuracy, and quantification of these internal inconsistencies provided = a basis=20 for evaluating registration accuracy in the absence of gold = standards.

Estimating Registration Accuracy TOP=20

When sufficiently accurate gold standards = are=20 available, quantification of registration accuracy is straight-forward. = For real=20 data from living humans, gold standards sufficient to quantify subvoxel = accuracy=20 are difficult to achieve, especially when even movements between the = brain and=20 the skull cannot be disregarded as negligible (11).=20 Accurate gold standards can be produced for phantoms or for simulated = data sets,=20 but the resulting assessments may be unrealistically optimistic since = these=20 methods generally fail to model all of the factors that contribute to=20 scan-to-scan variability in real human data. These factors include = noise,=20 non-rigid body movements due to respiratory or cardiac cycles (18),=20 artifacts due to rigid body movement during scan acquisition, = distortions due to=20 field inhomogeneities (2,19),=20 and inaccuracies in distance calibration (20).=20 For PET data, noise and spatial resolution are probably the main factors = limiting registration accuracy. The PET phantom models these particular = factors=20 reasonably well, and the accuracies on the order of 2 mm can reasonably = be=20 extrapolated to humans.

For MRI data, phantom studies and = simulations suggest=20 that registration accuracies of a fraction of a millimeter are possible, = but=20 this may have little bearing on real human data. Hajnal et al.(8,11)=20 have suggested that a statistical method(16)=20 applicable to least-squares minimization can be used to assess = registration=20 accuracy for real human data. However, the errors estimated by this = method are=20 only the errors in identifying the true minimum of the cost function. = Errors in=20 the underlying assumption that the images are correctly registered when = the cost=20 function is minimized are not evaluated. When applied to the MRI = least-squares=20 minimizations reported here, this method confirms that the true minimum = of the=20 cost function should be identified within 10 =CE=BCm, but the associated = internal=20 inconsistencies and discrepancies between strategies indicate that the = true=20 accuracy cannot be uniformly this good.

Even in the absence of gold standards, = internal=20 inconsistencies place an irrefutable limit on the degree of registration = accuracy that can be achieved since the true inaccuracy as measured = against any=20 gold standard is guaranteed to be at least as large as the internal=20 inconsistency. However, true accuracy can always be worse than the = internal=20 inconsistencies imply, so a strategy that performs worse in terms of = internal=20 consistency might nonetheless turn out to be the more accurate strategy = as=20 measured against independent gold standards. All other things being = equal, it=20 would seem prudent to prefer a method that produces significantly = smaller=20 internal inconsistencies over some less consistent method; but so long = as the=20 internal inconsistencies do not indicate intolerably poor accuracy, any = valid=20 advantage could provide sufficient justification for use of a less = consistent=20 method.

Choice of Minimization Procedure TOP=20

The performance of the multivariate = calculus-based=20 minimization procedures in AIR 3.0 matched and sometimes substantially = exceeded=20 the univariate minimization procedure used by the original AIR = algorithm. The=20 full Newton-type minimization procedure is a direct implementation of = the=20 theoretical foundation of all calculus- based minimization procedures = for image=20 registration (16).=20 Other minimization procedures, such as those used by Hajnal et al. (8),=20 Alpert et al.(13),=20 and Friston et al. (12)=20 are based on approximations to full Newton-type minimization, as is the=20 Marquardt-like approximation described here. All these approaches are = valid from=20 a theoretical standpoint (16),=20 and in many contexts the differences in speed and accuracy between them = are=20 probably unimportant. Use of the predicted change of the cost function = as the=20 primary termination criterion for minimization avoids the registration = failures=20 that can be seen in the presence of large movements when using a small, = fixed=20 number of iterations (12).=20 For the RIU cost function, a primary convergence criterion of 10-5 is generally sufficient for intrasubject = registration of=20 PET data, but a more stringent criterion might be better for MRI data. = For the=20 least-squares cost functions with 16 bit data (12 significant bits), a = primary=20 termination criterion around 1.0 is effective. For 8 bit data, a lower = value=20 around 0.1 may be appropriate. The ideal termination criteria may vary = as a=20 function of the nature of the data, so these values should be viewed as = only=20 approximate guidelines.

Choice of Cost Function and = Interpolation Model=20 TOP=20

The LS cost function performed = significantly worse=20 than the RIU cost function for PET data. This probably relates to the = fact that=20 the LS cost function was unable to take differences in global intensity = scaling=20 of the images into account. The SLS cost function, which includes an = explicit=20 intensity scaling factor, did not perform significantly worse than the = RIU cost=20 function for PET data. For MRI data, no significant differences among = cost=20 functions were found, except that the RIU cost function failed with = smoothed=20 data unless initialized with parameters already very close to the = correct=20 answer. The SLS cost function is probably the best overall choice for = MRI data=20 since it is insensitive to global intensity differences by design and = was robust=20 and generally the fasted in all contexts.

Snyder (14)=20 has suggested that the SLS cost function should be superior to the RIU = cost=20 function on theoretical grounds and has presented phantom PET data using = the LS=20 cost function in support of this suggestion. The discrepancy between = Snyder's=20 PET results and those presented here may be due to differences in = methodology=20 since Snyder kept the phantom's position stationary in all scans and=20 quantitatively scaled the image intensities to a consistent value. This = approach=20 does not simulate interpolation errors or positional distortions related = to=20 image acquisition, nor does it evaluate performance when image intensity = is=20 variable. Eberl et al. (21)=20 have also performed PET phantom studies and concluded that a cost = function=20 similar to the LS cost function is more accurate than the RIU cost = function.=20 However, interpretation of their results is clouded by the fact that = their=20 independent implementation of the RIU cost function produced errors far = larger=20 than have been reported previously(1)=20 or than were identified here.

Although sinc interpolation is = theoretically the ideal=20 interpolation function for resampling hand-limited MRI data when a fully = 3D=20 acquisition technique is used, the results presented here suggest that = trilinear=20 interpolation is adequate for registration. The results thus obtained = can then=20 be used to resample the data using sinc or chirp-z=20 interpolation to generate high quality final images.

Image Smoothing and Editing TOP=20

For PET registration, smoothing has been = shown=20 previously to improve registration accuracy of noisy images (1).=20 For the MRI data presented here, modest smoothing with a 2 to 4 mm = isotropic=20 Gaussian filter improved the internal consistency of registration = results.=20 Smoothing systematically changed the registration results more than = expected=20 based on the associated internal consistencies, showing that optimal=20 minimization of the cost function does not necessarily imply optimal=20 registration of the images. In the absence of gold standard MRI data = that would=20 allow definitive determination of which resolution gives the most = accurate=20 results, smoothing of MRI images with a 2 to 4 mm Gaussian filter is a=20 reasonable approach.

Hajnal et al. (11)=20 recommend routine editing of nonbrain structures for intrasubject = registration=20 since the brain can move with respect to the skull between image = acquisitions.=20 This is presumably less likely to be problematic during a single imaging = session=20 than with images acquired on separate days. If no relative movement has = actually=20 occurred, such editing effectively throws away sharply defined = boundaries that=20 could have served to improve registration accuracy. No evidence was = found to=20 suggest that movement of the brain relative to the skull occurred in MRI = data=20 sets used here, and editing sometimes significantly worsened internal=20 inconsistencies. Nonetheless, editing of the images is certainly = appropriate if=20 movements of the brain relative to the scalp and skull are known or = likely to be=20 present. Editing might also be helpful in situations where the rigid = body model=20 is significantly violated by imaging artifacts.

Spatial Transformation Model TOP=20

For intrasubject registration, a rigid = body model is=20 generally assumed to be the spatial transformation model of choice. = However, if=20 distances are poorly calibrated, a more general linear model with 7-11=20 parameters, depending on which distances are known to be inaccurate, may = give=20 better results. In data not presented here, use of a 12 parameter affine = model=20 to register the MRI data sets resulted in internal inconsistencies that = were=20 worse than those obtained with the rigid body model. Alternatively, = calibration=20 errors estimated from anatomic landmarks within the images can be used = to adjust=20 voxel dimensions before registration with a rigid body model(20).=20 It is advisable to correct significant spatial or intensity distorations = associated with magnetic susceptibility artifacts(19)=20 or with the use of local gradient coils(22)=20 before registration since these distorations are not readily modeled = with linear=20 scaling of spatial distortions or intensity. For extremely fast MR = techniques=20 such as echo planar imaging, it may be appropriate to consider the use = of=20 nonlinear spatial transformation models(6)=20 to compensate for distorations associated with respiratory and cardiac = cycles.=20 For slower imaging techniques, it should be kept in mind that, at best, = the=20 acquired images represent a temporal average of an object that is = subject to=20 constant small scale physiologic fluctuations. To the extent that this = temporal=20 averaging varies from one image to the next, the very notion of an = absolutely=20 correct spatial transformation model may be unjustified.

Registration Speed TOP=20

If reduced accuracy can be tolerated, = sparse sampling=20 of the data can be used to substantially increase registration speed. = For PET=20 data with anisotropic voxels, sparse sampling is a more accurate way to = increase=20 speed than omitting the step of interpolating the reference volume to = cubic=20 voxels. If initial misregistration is large, the secondary termination = criteria=20 may need to be adjusted when using sparse sampling to allow a sufficient = total=20 number of iterations.

Generality of Results TOP=20

Because of its generality, AIR is not = restricted to=20 any specific modality, species, or anatomic structure. The systematic=20 investigation of various factors that might alter speed or accuracy = reported=20 here provides insights that are applicable not only to AIR, but also to = other=20 similar registration algorithms. To the extent that the results may = depend on=20 the specific data sets investigated here, the general approaches = described for=20 validation are nonetheless broadly applicable. This is particularly true = of the=20 internal consistency measures, which require only that three or more = images be=20 available for pairwise registration. We hope that the methods and = results=20 described here will help users of a variety of registration methods = ensure the=20 accuracy of their results while maintaining practical levels of = performance.

Acknowledgment: This=20 work was supported by National Institute of Neurological Disorders and = Stroke=20 grant 1 K08 NS-01646, Department of Energy contract DE-FCO3-87ER60615, = NIH grant=20 5P01MH52176 to the International Consortium for Brain Mapping, the = Canadian=20 Medical Research Council SP-30, generous gifts from the Pierson-Lovelace = Foundation, the Ahmanson Foundation, the Tamkin Foundation, North Star = Fund, and=20 the Brain Mapping Medical Research Organization.

REFERENCES TOP=20

1. Woods RP, Cherry SR, Mazziotta JC. = Rapid=20 automated algorithm for aligning and reslicing PET images. J=20 Comput Assist Tomogr 1992;16:620-33.=20

[Context=20 Link]

2. Jiang AP, Kennedy DN, Baker JR, et = al. Motion=20 detection and correction in functional MR imaging. Hum Brain=20 Map 1995;3:224-35.=20

[Context=20 Link]

3. Woods RP, Mazziotta JC, Cherry SR. = MRI-PET=20 registration with automated algorithm. J Comput = Assist=20 Tomogr 1993;17:536-46.=20

[Context=20 Link]

4. Strother SC, Anderson JR, Xu XL, Liow = JS, Bonar=20 DC, Rottenberg DA. Quantitative comparisons of image registration = techniques=20 based on high-resolution MRI of the brain. J = Comput Assist=20 Tomogr 1994;18:954-62.=20

[Context=20 Link]

5. West J, Fitzpatrick JM, Wang MY, et = al.=20 Comparison and evaluation of retrospective intermodality image = registration=20 techniques.J Comput Assist Tomogr = 1997;21:554-66.=20

[Fulltext=20 Link] [CrossRef]=20 [Context=20 Link]

6. Woods RP, Grafton ST, Watson JDG, = Sicotte NL,=20 Mazziotta JC. Automated image registration: II. Intersubject validation = of=20 linear and nonlinear models. J Comput Assist = Tomogr=20 1998;22:155-67.=20

[Context=20 Link]

7. Rabiner LR, Schafer RW, Rader CM. The = chirp=20 z-transform algorithm and its application. Bell = System Tech=20 J 1969;48:1249-92.=20

[Context=20 Link]

8. Hajnal JV, Saeed N, Soar EJ, Oatridge = A, Young=20 IR, Bydder GM. A registration and interpolation procedure for subvoxel = matching=20 of serially acquired MR images. J Comput Assist = Tomogr=20 1995;19:289-96.=20

[Context=20 Link]

9. Unser M, Thevenaz P, Yaroslavsky L.=20 Convolution-based interpolation for fast, high-quality rotation of = images. IEEE Trans Image Proc 1995;4:1371-81.=20

[Context=20 Link]

10. Fraser D, Schowengerdt RA. Avoidance = of=20 additional aliasing in multipass image rotations. IEEE Trans=20 Image Proc 1994;3:721-35.=20

[Context=20 Link]

11. Hajnal JV, Saeed N, Oatridge A, = Williams EJ,=20 Young IR, Bydder GM. Detection of subtle brain changes using subvoxel=20 registration and subtraction of serial MR images. J Comput=20 Assist Tomogr 1995;19:677-91.=20

[Context=20 Link]

12. Friston KJ, Ashburner J, Frith CD, = Poline JB,=20 Heather JD, Frackowiak RSJ. Spatial registration and normalization of = images.Hum Brain Map 1995;3:165-89.=20

[Context=20 Link]

13. Alpert NM, Berdichevsky D, Levin Z, = Morris ED,=20 Fischman AJ. Improved methods for image registration. Neuroimage 1996;3:10-8.=20

[CrossRef]=20 [Context=20 Link]

14. Snyder AZ. Difference image versus = ratio image=20 error function forms in PET-PET realignment. In: Quantification of brain function using PET. San = Diego:=20 Academic Press, 1996;131-7.=20

[Context=20 Link]

15. Woods RP, Mazziotta JC, Cherry SR. = Automated=20 image registration. In: Uemura K, Lassen NA, Jones T, Kanno I, eds.Quantification of brain function: tracer kinetics and = image=20 analysis in brain PET. Amsterdam: Excerpta Medica, 1993;391-8.=20

[Context=20 Link]

16. Press WH, Teukolsky SA, Vetterling = WT, Flannery=20 BP.Numerical recipes in C: the art of scientific=20 computing. 2nd ed. Cambridge: Cambridge University Press, 1992.=20

[Context=20 Link]

17. Hoffman EJ, Cutler PD, Guerrero TM, = Digby WM,=20 Mazziotta JC. Assessment of accuracy of PET utilizing a 3D phantom to = simulate=20 the activity distribution of [¹⁸F]fluorodeoxyglucose=20 uptake in the human brain. J Cereb Blood Flow = Metab=20 1991;11:A17-25.=20

[Context=20 Link]

18. Poncelet BP, Wedeen VJ, Weisskoff = RM, Cohen MS.=20 Brain parenchyma motion: measurement with cine echo-planar MR = imaging.Radiology 1992;185:645-51.=20

[Context=20 Link]

19. Chang H, Fitzpatrick JM. A technique = for=20 accurate magnetic resonance imaging in the presence of field = inhomogeneities.IEEE Trans Med Imag 1992;11:319-29.=20

[Context=20 Link]

20. Freeborough PA, Woods RP, Fox NC. = Accurate=20 registration of serial 3D MR brain images and its application to = visualizing=20 change in neurodegenerative disorders. J Comput = Assist=20 Tomogr 1996;20:1012-22.=20

[Fullt= ext=20 Link] [CrossRef]=20 [Context=20 Link]

21. Eberl S, Kanno I. Fulton RR, Ryan A, = Hutton BF,=20 Fulham MJ. Automated interstudy image registration technique for SPECT = and=20 PET.J Nucl Med 1996;37:137-45.=20

[Context=20 Link]

22. Zeffiro TA, VanMeter JW, Eden GF, = Woods RP.=20 Signal gradient normalization improves head motion correction in = surface-coil=20 functional magnetic resonance imaging. Neuroimage=20 1996;3:S110.=20

[Context=20 Link]

Index Terms: Image registration; Magnetic = resonance=20 imaging; Emission computed tomography; Brain mapping

=C2=A9 Lippincott-Raven = Publishers