The present invention relates to aligning scanned images. In the exemplary embodiment provided in the present invention, two different intensity values of two images I1 and I2,T are modeled by random variables X1 and X2. For the purposes of the exemplary embodiment of the method, I1 corresponds to a reference (known) image and T is a spatial transformation applied to a study image I2 to produce a transformed image, I2,T.
For any given pixel location x, I2,T(x)=I2(T(x)). The value T belongs to an arbitrary class K of mappings ranging from rigid transformations to high-dimensional nonrigid deformations. Images are acquired on a finite sampling grid and values of I2,T are obtained from the study image I2 by interpolation. The joint density of the pair (X1, X2,T) is denoted, px
As opposed to the KL, MI and NMI evaluation techniques, the Earth Mover's Distance is used to allow for better evaluation of histogram similarity independently of quantization effects and deformations of the feature space. Earth Mover's Distance is referred to as a “cross-bin” distance as opposed to “bin to bin” distances such as used in KL. When applied to normalized histograms, Earth Mover's Distance is a discrete version of the Wasserstein metric in probability theory. This metric is the solution of an optimal mass transportation problem. The Earth Mover's Distance corresponds to the cost of optimally transporting one distribution into the other given a ground distance. The ground distance is typically derived from an Lp norm, for instance dp(x,y)=∥x−y∥p where (p=2 yields the Euclidean distance). More formally, given two normalized histograms h1 and h2, whose ith bins are centered at locations xi and yi, respectively, Earth Mover's Distance is defined by
Subject to the following constraints:
These values fi,j can be seen as elementary flows transporting elements of h1 and h2 (or h2 to h1). This approach generalizes to multidimensional histograms such as joint intensity distributions. When the ground distance is arbitrary, the optimal solution to the linear programming problem is given by a transportation simplex algorithm.
With this definition, two Earth Mover's distance based similarity measures are provided. The first measure is:
εa(I1,I2,T)=EMDL
This value is maximized to achieve registration. The second is used when a learned joint intensity distribution p0 is available.
εb(I1,I2,T)=EMDL
εb is minimized in order to align a pair of images. The implementation of these new measures raises several practical challenges. Computing Earth Mover's Distance for a pair of joint distribution and histograms involves solving a linear programming problem. The registration process requires multiple evaluations of Earth Mover's Distance on two dimensional histograms. To accomplish this, an algorithm providing the following 6 steps of pseudo-code allows for a compact solution.
Using this approach, the linear programming problem can be drastically simplified in order to get a much more efficient simplex algorithm. Computation times make the use of EMDL possible for registration purposes.
Referring to
εb(I1,I2,T)=EMDL
minimized to align the images as completely as possible.
Referring to
Experiments with two pairs of Magnetic Resonance brain images were performed. The first pair was acquired from a volunteer with a Siemens Magnetom Avanto 1.5T machine. One of the datasets was obtained with a T2-weighted HASTE sequence (matrix size 512×384, 22 slices, voxel size of 0.45×0.45×5 mm3) and the other with a T1 weighted Spin-Echo sequence (matrix size 192×192, 36 slices, voxel size of 1.2×1.2×3 mm3) A volunteer was asked to move his head after the first acquisition to simulate a misalignment of large amplitude. A pair of simulated T1 and T2 weighted MR images from a database are used to generate ground truth data and measure registration errors. The database had a volume size of 181×217×181 and isotropic voxels (1 mm3). Their intensity non-uniformity was set to 20%. The noise level of the T1 and T2 weighted images was set to 3% and 9% respectively.
In all registration experiments, the number of quantization levels was set to 16. This value was chosen empirically and provides a good compromise between registration accuracy and computational efficiency. Joint histograms were computed using partial volume interpolation. In the rigid registration experiments we applied a multi-resolution hill-climbing optimization strategy. With a 3 GHz microprocessor, the evaluation EMDL takes on average 1 ms for histograms of size 162. This computation time increases dramatically to 15 ms for a size of 322 and up to 313 ms for 642.
In the first experiment, the value of εa and NMI evolves when artificial rotations and translations are applied to one of the volumes. T1 and T2 weighted datasets which were perfectly aligned were used for analysis. Profiles of function εa with respect to rotation (
Like NMI, εa reaches a peak value at 0 and has a large capture range.
A first qualitative evaluation εa is conducted by rigidly aligning the real magnetic resonance data (T2 and T1 weighted scans with a large displacement). The configuration of the images before and after registration are shown in
In the next experiment, the average registration errors of εa− and NMI based rigid registration algorithms for 400 points evenly placed inside the T1− and T2− weighted datasets are computed. The results of the evaluation are presented in
In the final experiment, a prior distribution from aligned brain images is estimated and compared to the profiles of εb and KL with respect to artificial transformations when the intensity values of one of the volumes have been distorted. A gamma compression/expansion operator and an intensity offset are applied to after the intensity profile of the images (with y equal to 0.3 and 1.7, and an offset value equal to 400). These operations are performed before quantization of the data.
In
The present invention provides for an advantage compared to conventional technologies that the method efficiently aligns two image data sets. The alignment of the image data sets allows physicians the capability to not only evaluate a procedures progress with respect to two dimensions, but also in a third dimension, if needed as the alignment is highly accurate. This additional capability, allows for safer and less intrusive medical procedures for patients.
The method and apparatus of the present invention may be used with multiple types of data acquisition systems. An example of this are MRI systems. Other types of scanning may be used and are applicable to use with the method of the present invention. These types of scanning include, for example, x-ray, ultrasound or computed tomography (CT) scanning.
In the foregoing specification, the invention has been described with reference to specific exemplary embodiments thereof. It will, however, be evident that various modifications and changes may be made thereunto without departing from the broader spirit and scope of the invention as set forth in the appended claims. The specification and drawings are accordingly to be regarded in an illustrative rather than in a restrictive sense.
This is a United States non-provisional application of U.S. provisional patent application Ser. No. 60/836,596 filed Aug. 9, 2006 by Christophe Chefd'hotel and Guillaume Bousquet, the entirety of which application is incorporated by reference herein.
Number | Date | Country | |
---|---|---|---|
60836596 | Aug 2006 | US |