CROSS-REFERENCE TO RELATED APPLICATION

The present application claims priority to U.S. Provisional Patent Application No. 62/138,006, filed Mar. 25, 2015. The content of that patent application is hereby incorporated in its entirety.

GOVERNMENT RIGHTS

The subject matter disclosed herein was made with government support under award/contract/grant number HL 105212 awarded by the Department of Health and Human Services. The Government has certain rights in the herein disclosed subject matter.

TECHNICAL FIELD

The invention relates to a method for interactive non-uniformity correction and intensity standardization of magnetic resonance (MR) images using sample tissue regions specified for several different types of tissues.

BACKGROUND

Image non-uniformity and intensity non-standardness are two major hurdles encountered in human and computer interpretation and analysis of magnetic resonance (MR) images. The former phenomenon is caused by imperfections in the imaging device, specifically the magnetic field, and has been studied extensively, where many solutions have been proposed in the literature. See, for example, Sled JG, et al. â€œA nonparametric method for automatic correction of intensity non-uniformity in MRI data,â€ IEEE Trans. Medical Imaging, vol. 17, 87-97 (1998); Tustison, N. J., Avants, B. B., Cook, P. A., Zheng, Y., Egan, A., Yushkevich, P. A., Gee, J. C., â€œN4ITK: improved N3 bias correction,â€ IEEE Transaction on Medical Imaging, 29 (6), 1310-1320 (2010); and Vovk U, Pernus F, Likar B., â€œA review of methods for correction of intensity inhomogeneity in MRI,â€ IEEE Trans Med Imaging, 26(3), 405-421 (2007). Generally, the intensity non-uniformity component is assumed to be multiplicative or additive, i.e., the component is multiplicative or additive to an ideal image. Most frequently, the multiplicative model has been used as it is consistent with the inhomogeneous sensitivity of the receiver coil of the magnetic resonance imaging (MRI) scanner. The N3 and N4 algorithms described by Sled et al. and Tustison et al. are two popular correction approaches which have adopted the multiplicative model.

Although to a far lesser extent, solutions have also been presented for the second problem of intensity non-standardness. See, for example, NyÃºl LG, Udupa J. K., â€œOn standardizing the MR image intensity scale,â€ Magn. Reson. Med., 42 (6), 1072-1081 (1999); and NyÃºl LG, Udupa J. K., â€œNew variants of a method of MRI scale standardization,â€ IEEE, Trans. Medical Imaging, vol. 19, 143-150 (2000). This phenomenon causes the lack of a tissue specific numeric intensity meaning, even within the same MRI protocol, for the same body region, for images obtained on the same scanner, and for the same patient. It is also known that all non-uniformity correction methods introduce their own non-standardness as part of the correction process itself (see Madabhushi A, Udupa J. K., â€œInterplay between intensity standardization and inhomogeneity correction in MR image processing,â€ IEEE Trans. Medical Imaging, vol. 24, 561-576 (2005)). The interplay between the two processes, proper ways of handling the phenomena in combination, and the manner in which image segmentation performance improves substantially when these problems are addressed properly have also been demonstrated (see, for example, Bai W., Shi W., Ledig C., Rueckert D., â€œmulti-atlas segmentation with augmented features for cardiac MR,â€ Medical Image Analysis, 19(1), 98-109 (2015); Wu. S., Udupa J. K., Marinaki A., Weinstein S. P., Kontos D., â€œIntensity standardization in breast MR images improves tissue quantification,â€ Medical Imaging: Proc. SPIE, vol. 8668, 866822-1-866822-6 (2013); Ge Y, Udupa J. K., et al. â€œNumerical tissue characterization in MS via standardization of the MR image intensity scale,â€ Journal of Magn. Reson. Imag, vol. 12, 711-721 (2000); Zhuge Y, Udupa J. K., â€œMembership-based multiprotocol MR brain image segmentation,â€ Medical Imaging: Proc. SPIE, vol. 5032, 1572-1579 (2003); and Zhuge, Y., Udupa J. K., et al., â€œImage background inhomogeneity correction in MRI via intensity standardization,â€ Comp. Med. Imag. and Graph., 33(1),7-16 (2009)).

However, in many practical situations, automated methods fail to yield acceptable solutions for both problems. This is mainly because solutions for both problems require identifying regions representing the same tissue type for several different tissues, and the automatic strategies, irrespective of the approach, may fail in this task. It remains desirable to provide interactive strategies to overcome this problem in both phenomena wherein the required high level knowledge is provided by an operator to improve image quality substantially in such situations.

SUMMARY

Interactive non-uniformity correction and intensity standardization approaches are described that address the above and other needs in the art by providing improvements in the image based on MR images of the abdomen in an exemplary embodiment. These approaches will be referred to herein as iNC and iIS, respectively, for interactive non-uniformity correction and interactive intensity standardization.

The methods presented herein provide interactive strategies to overcome the problems noted above by requiring sample tissue regions to be specified for several different types of tissues. Interactive NC estimates the degree of non-uniformity at each voxel in a given image, builds a global function for non-uniformity correction, and then corrects the image to improve quality. Interactive IS includes two steps: a calibration step and a transformation step. In the calibration step, tissue intensity signatures of each tissue from a few subjects are utilized to set up key landmarks in a standardized intensity space. In the transformation step, a piecewise linear intensity mapping function is built between the same tissue signatures derived from the given image and those in the standardized intensity space to transform the intensity of the given image into standardized intensity. Preliminary results on abdominal T1-weighted and T2-weighted MR images of 20 subjects show that interactive NC and IS are feasible and can significantly improve image quality over automatic methods. Interactive IS for MR images combined with interactive NC can substantially improve numeric characterization of tissues.

In exemplary embodiments, the invention includes a system and method for providing interactive non-uniformity correction (NC) and interactive intensity standardization (IS). In an exemplary embodiment, the method includes the steps of:

specifying sample tissue regions for several different types of tissues;

estimating a degree of non-uniformity at each voxel in a sample image of the sample tissue regions;

building a global function for non-uniformity correction of the sample image;

correcting the sample image to improve image quality;

using tissue intensity signatures of each sample tissue from a few subjects to set up key landmarks in a standardized intensity space of the corrected sample image; and

building a piecewise linear intensity mapping function between the same tissue signatures derived from the corrected sample image and those in the standardized intensity space to transform the intensity of the corrected sample image into a standardized intensity.

In exemplary embodiments, the specifying step comprises selecting sample tissue that is normal or near normal, selecting a sample tissue region that is spatially well-distributed throughout the image volume, and selecting a sample tissue region that can be reliably sampled by an interactive operation. The specifying step further comprising the step of creating a binary mask file for each tissue type in the sample image by drawing tissue regions having only one tissue type.

In the exemplary embodiments, the estimating step comprises the steps of, for each image in a given set of images and a tissue region in each the image, estimating the degree of non-uniformity at each voxel in the tissue region by comparing a value of the voxel to a mean intensity of the tissue region to estimate either an additive of a multiplicative factor, and fitting a second degree polynomial to the estimated non-uniformity factors in the tissue region. The step of building a global function comprises combining functions fitted to the second degree polynomial from all tissue regions in each image.

In the exemplary embodiments, determining the tissue intensity signatures of each sample tissue includes finding mean intensities corresponding to each tissue region and performing a piecewise linear transformation or curve fitting to the mean intensities.

The scope of the invention also includes a system including a processor that executes stored instructions for executing the steps of the method. The above and other characteristic features of the invention will be apparent from the following detailed description of the invention.

BRIEF DESCRIPTION OF THE DRAWINGS

The present application is further understood when read in conjunction with the appended drawings. In the drawings:

FIG. 1 illustrates the iIS operation where the abscissa represents the mean intensities of the reference tissues in a test image and the ordinate denotes the standardized reference intensities.

FIG. 2 illustrates a display of a sample slice from the images of three subjects. Top row: D. Middle row: D_{aNC, aIS}. Bottom row: D_{iNC, iIS}.

DETAILED DESCRIPTION OF ILLUSTRATIVE EMBODIMENTS

Certain specific details are set forth in the following description with respect to FIGS. 1-2 to provide a thorough understanding of various embodiments of the invention. Certain well-known details are not set forth in the following disclosure, however, to avoid unnecessarily obscuring the various embodiments of the invention. Those of ordinary skill in the relevant art will understand that they can practice other embodiments of the invention without one or more of the details described below. Also, while various methods are described with reference to steps and sequences in the following disclosure, the description is intended to provide a clear implementation of embodiments of the invention, and the steps and sequences of steps should not be taken as required to practice the invention.

Materials and Methods

Image Data: Abdominal MR images of acceptable quality from 20 subjects (between 26-85 years) with normal appearing livers as per a board-certified radiologist's (DAT) review were utilized for illustrating the new iNC and iIS approaches. All MR images had been previously acquired on 1.5T MRI scanners. In-phase T1-weighted gradient recalled echo (TR=150 msec, TE=4 msec) and non-fat suppressed heavily T2-weighted fast spin echo (TRËœ490 msec, TE=90 msec) images were utilized for this study. Scene size is 512Ã—432Ã—20-30 and voxel size is 0.78Ã—0.78Ã—8.00 mm³.

Notations: NCâ€”non-uniformity correction; ISâ€”intensity standardization; aNCâ€”automated NC by the method described by Zhuge et al.; N3â€”NC method of Sled et al.; aISâ€”automated IS by the method described by Nyul et al.; iNC, iISâ€”interactive NC and interactive IS, respectively, by the new approach presented herein; Dâ€”a given set of images; D_câ€”the set of images resulting from applying a NC method c to the images in D; and D_c,sâ€”the set of images resulting from applying an IS method s to the images in D_c. If an IS method s is applied directly on the original images in D, the resulting set will be denoted by D. In the inventor's study, câˆˆ{aNC, N3, iNC} and sâˆˆ{aIS, iIS}. t1, . . . ,tnâˆ’n different types of tissues. D^t1, . . . , D^tnâ€”the sets of binary images corresponding to the images in D where the binary images represent masks indicating sample tissue regions for t1, . . . ,tn, respectively.

The basic idea behind both iNC and iIS approaches is to first decide on the different tissue types that can be used as reference and to subsequently delineate sample regions for these tissues in the images in D. The selection of tissue type should be guided by three criteria: (i) the tissue is normal or near normal, (ii) the tissue region is spatially well-distributed throughout the image volume, and (iii) the tissue region can be reliably sampled by an efficient interactive operation such as painting with broad brush strokes. The sampled tissue regions are subsequently utilized for performing both iNC and iIS operations.

For each image in the given set D of images, a binary mask file for each of n tissue types t1, . . . , tn is created by drawing tissue regions guided by the above criteria. The mask does not need to be the exact delineation of the whole boundary of the tissue region, but the marked regions should be pure in the sense of containing only one tissue type. Any combination of quick interactive segmentation tools can be used for this purpose such as painting by broad brush strokes. The inventors perform this and all other operations described herein by using the open source software system CAVASS (see Grevera G., Udupa J. K., et al., â€œCAVASS: A computer-assisted visualization and analysis software system,â€ Journal of Digital Imaging, vol. 20, 101-118 (2007)). The inventors employ the same set of mask data sets D^t1, . . . , D^tnfor both iNC and iIS processes as described below. In other words, user interaction is needed only once and not separately for each process.

Interactive Non-Uniformity Correctionâ€”iNC

For iNC, the inventors are given the original image set D and the associated binary image sets D^t1, . . . , D^tn. For each image I in D and the tissue region R_icorresponding to ti in I, the inventors estimate the degree of non-uniformity at each voxel in R_ias traditionally done by comparing the voxel value to the mean intensity in R_ifor estimating either an additive or a multiplicative factor. The inventors fit a second degree polynomial to the estimated non-uniformity factors in R_i. In this process, the inventors treat each tissue region as an independent entity. Subsequently, the inventors combine the fitted functions from all tissue regions into a global function. The image I is then corrected using the combined global function. For the input set D, the output image set produced by the iNC operation will be D_iNC.

Interactive Intensity Standardizationâ€”iIS

The iIS operation consists of two stepsâ€”a one-time calibration step and a transformation step that is performed for each image to be standardized. To reap the full benefits of standardization, as recommended by Madabhushi, the inventors perform iIS after iNC. Otherwise, if iIS is performed first followed by iNC, as demonstrated by Madabhushi, residual non-standardness introduced by the (any) NC method will be left over. That is, the inventors perform iIS on the set D_iNCresulting from applying iNC on D. Note that iIS can also be applied directly on D, or D_cfor any câˆˆ{aNC, N3, iNC}, or even D_aIS. Below D_cwill be used to denote any input image set for the iIS operation.

The purpose of the calibration step is to estimate parameters of the standard intensity scale, as shown in FIG. 1. Some sample images from the input set D_care used for this purpose. The parameters of the standard scale are the â€œsignaturesâ€ Î¼_t1, . . . , Î¼_mof the different reference tissues. In this application, the inventors take Î¼_t1, . . . , Î¼_mto be the mean of the intensities in the n tissue regions R_iwithin the sample images of D_c. Instead, other tissue-specific signatures may also be employed. Additionally two other parameters are usedâ€”p₁which denotes the lowest meaningful intensity, and p₂that indicates the largest meaningful intensity. As in Nyul et al., the inventors take p_ito be the smallest intensity over all sample images and p₂to be the intensity at the 99.8 percentile within the sample images. By observing thousands of MR images under different protocols and body regions, it has been found by Nyul et al. that there is usually a distinct cut off at 99.8 percentile beyond which the intensities represent mostly outlier values.

In the transformation step, given any image IâˆˆD_c, it is standardized by first finding its tissue signatures, in this case the mean intensities m_t1, . . . , m_tncorresponding to its n tissue regions, and then performing a piecewise linear transformation, as illustrated in FIG. 1. Instead of the piecewise linear mapping, any curve fitting method may also be used. The set of images that result will be, as per the above notation, D_c,iIS. Note that the transformation extends beyond (p₂, Î¼₂) along the positive axes directions and beyond (p₁, Î¼₁) along the negative axes directions following the corresponding linear segments. In this manner, no input intensities are lost and a 1:1 mapping can be guaranteed.

To sum up, user interaction is needed by the proposed approach only once for each image for performing both iNC and iIS operations. As pointed out earlier, many combinations of the NC and IS operations and their iterative executions are possible, although they are not explored herein.

Results

The inventors use the data set D described in the previous section for testing the iNC and iIS approaches. For the calibration step of iIS, the inventors used 5 images from the set D_iNC. For comparison, the inventors also applied the automatic methods aNC and aIS in that order to D to yield the non-uniformity corrected set D_aNCand corrected and standardized image set D_{aNC, aIS}. The inventors used two tissue types (n=2), fat and muscle, as reference for the correction and standardization operations.

Note that many different combinations of operations can be carried out and compared for the quality of the final images. Following the above notation, the combinations which are compared in this initial study are as follows. For NC: D, D_N3, D_aNC, and D_iNC. For IS: D, D_aNC,aIS, and D_iNC,iIS.

To assess the performance of an NC operation câˆˆ{aNC, N3, iNC}, the inventors estimate the percent coefficient of variation (% cv) of intensity within each specific tissue region in each image I in the resulting image set D_cand then describe the mean and standard deviation of this % cv value over all images in this set. Percent cv values within a tissue region in an image I in D_cgive an indication of the amount of residual non-uniformity left over in I as attested by the tissue region. For this evaluation, in addition to the fat and muscle tissue regions used as reference, the inventors use the liver object region as an independent additional test tissue region. To facilitate this analysis, the inventors segmented the liver 3D region in all images by using the user-driven live wire method. The image sets have been selected carefully to make sure that the liver region is normal in every image.

For assessing the performance of an IS operation s in conjunction with an NC operation c in the resulting image set D_c,s, the mean intensity is estimated in each tissue region in each image I in D_c,s, and the mean and standard deviation of this mean intensity are described over all images in the set. This is a sensible approach since non-standardness is a trait of not a given image but of how the same tissue property changes over a subject population. The inventors use the same three tissue regions for IS evaluation as for NC.

Table 1 lists the mean and standard deviation of the % cv values for NC operations on T1-weighted in-phase images for the sets D, D_N3, D_aNC, and D_iNC. The p-values resulting from separate t-tests comparing D_N3, D_aNC, and D_iNCwith D are also listed in the table. Similarly, Table 2 lists % cv and p-values for the T2-weighted images. Tables 3 and 4 analogously display the mean and standard deviation of mean tissue intensities for T1- and T2-weighted images illustrating the performance of the IS operations.

TABLE 1Results for non-uniformity correction. Mean and standard deviationof % cv values for T1-weighted in-phase MR images.p-p-p-LivervalueFatvalueMusclevalueD17.66 (3.05)â€”16.82 (4.48)â€”â€‚17.27 (15.66)â€”D_N318.45 (2.86)0.1511.53 (5.59)0.0â€‚12.64 (4.89)0.08D_aNC17.50 (7.31)0.4619.85 (9.70)0.0615.16 (6.26)0.25D_iNC17.35 (5.26)0.3614.01 (5.87)0.0111.76 (3.40)0.03

TABLE 2Results for non-uniformity correction. Mean and standarddeviation of % cv values for T2-weighted MR images.p-p-p-LivervalueFatvalueMusclevalueD33.05 (9.47)â€”18.40 (4.97)â€”36.25 (7.58)â€”D_N333.44 (8.22)0.35â€‚9.65 (3.32)0.0â€‚36.79 (7.33)0.2â€‚D_aNC33.64 (9.42)0.4617.17 (6.99)0.06â€‚45.48 (18.94)0.25D_iNC30.83 (9.04)0.3613.38 (3.24)0.0134.89 (6.95)0.04

TABLE 3Intensity standardization. Mean and standard deviation oftissue mean intensities for T1-weighted in-phase MR images.DD_{aNC, aIS}D_{iNC, iIS}LiverFatMuscleLiverFatMuscleLiverFatMuscleMean325.21575.42273.551460.82534.641248.321713.332567.201278.11intensitySD191.24282.72133.92161.47264.65294.91183.1187.0224.38cv %58.8149.1348.9611.0610.4423.6210.693.391.91

TABLE 4Intensity standardization. Mean and standard deviationof tissue mean intensities for T2-weighted MR images.DD_{aNC, aIS}D_{iNC, iIS}LiverFatMuscleLiverFatMuscleLiverFatMuscleMean128.40678.08117.51531.332475.21373.10647.862483.35497.19intensitySD56.32290.8969.79161.50193.04169.78166.6138.7720.96cv %43.8742.9059.3930.407.8045.5125.721.564.22

FIG. 2 displays one slice of the T1-weighted images from each of three subjects from the sets D, D_{aNC, aIS}, and D_{iNC, iIS}. The images from D shown in the first row are difficult to display at a fixed window gray level setting because of the intensity non-standardness present among subjects. The second row displays images from the set D_{aNC, aIS}where automatic NC and IS have been applied. The images clearly demonstrate an improvement when displayed at a fixed gray level window setting. The last row shows corresponding images for the same three subjects drawn from D_{iNC, iIS}resulting from interactive NC and IS. The fixed gray level window display seems to show the best portrayal of each type of tissue in the three subjects at about the same brightness level, clearly indicating the improved numeric characterization achieved for the tissues across subjects.

From Tables 1 and 2, the inventors observe that iNC generally performs better than the other two automatic NC methods. More importantly, the quality of the images obtained after iNC and iIS operations is substantially better than the quality of the images after the automated operations aNC and aIS for all tissue regions.

Conclusions

The inventors have presented an interactive strategy for correcting image intensity non-uniformities and non-standardness. The interaction required for both operations is only once in the form of regions of two or more tissue regions painted in each given image. Non-standardness correction has an initial one-time calibration step. The inventors presented an evaluation of the two methods comparing them to some automated methods from the literature on T1- and T2-weighted abdominal MR images. The overall quality improvement achieved by the two operations applied in tandem seems to be substantial in terms of the tissue specificity achieved by the numeric pixel values.

More extensive evaluations of other combinations of the NC and IS operations on other types of images from different body regions are underway. One disadvantage of the methods is the human interaction required on a per image basis. The inventors are investigating ways of reducing this interaction to a minimum where one can still achieve substantial quality advantage for the resulting images. The inventors are also studying severe cases of non-uniformity that can be handled satisfactorily by the proposed methods where other automatic methods may fail.

It will be appreciated that the methods described herein may be implemented in software that operates on a processor that executes instructions stored in a memory component. The processor may include a standardized processor, a specialized processor, a microprocessor, or the like. The processor may execute instructions including, for example, instructions for implementing the method as described herein. On the other hand, the memory component stores the instructions that may be executed by the processor. The memory component may include a tangible computer readable storage medium in the form of volatile and/or nonvolatile memory such as random access memory (RAM), read only memory (ROM), cache, flash memory, a hard disk, or any other suitable storage component. In one embodiment, the memory component may be a separate component in communication with a processor, while in another embodiment, the memory component may be integrated into the processor. Such non-transitory memory components may be used as a computer readable storage device to store the instructions for implementing the methods and software features described herein.

Those skilled in the art also will readily appreciate that many additional modifications are possible in the exemplary embodiment without materially departing from the novel teachings and advantages of the invention. Accordingly, any such modifications are intended to be included within the scope of this invention as defined by the following exemplary claims.