Applications of Automatic Anatomy Recognition In Medical Tomographic Imagery Based On Fuzzy Anatomy Models

Abstract

A computerized method of providing automatic anatomy recognition (AAR) includes gathering image data from patient image sets, formulating precise definitions of each body region and organ and delineating them following the definitions, building hierarchical fuzzy anatomy models of organs for each body region, recognizing and locating organs in given images by employing the hierarchical models, and delineating the organs following the hierarchy. The method may be applied, for example, to body regions including the thorax, abdomen and neck regions to identify organs.

Claims

1. A computerized method, comprising the steps of: building a fuzzy anatomy model of a body region of interest from an existing set of patient images for the body region; obtaining an image of a particular patient body region of interest, wherein obtaining the image comprises mapping slice locations associated with the particular patient body region of interest from a scanner coordinate system to slice locations in a standardized anatomic space; identifying where a volume-to-area correlation becomes maximal for one or more body regions of the patient; using automatic anatomy recognition (AAR) to recognize objects in the particular patient body region of interest and to output spatial information associated with the recognized objects; and delineating, using automatic anatomy recognition (AAR) and based at least on the spatial information associated with the recognized objects, the recognized objects to determine one or more contours or boundaries of the recognized objects.

2. A computerized method as in claim 1, wherein building a fuzzy anatomy model of the body region of interest comprises gathering image data from patient image sets, formulating precise definitions of each body region and organ and delineating them following said definitions, building hierarchical fuzzy anatomy models of organs for each body region, and selecting from the hierarchical fuzzy anatomy models the fuzzy anatomy model of the body region of interest.

3. A computerized method as in claim 2, wherein using AAR to recognize and delineate objects in the body region of interest comprises recognizing and locating organs in the image by employing the hierarchical models, and delineating objects in the body region of interest following a hierarchy of the hierarchical fuzzy anatomy models.

4. A computerized method as in claim 2, wherein object size and positional relationships of objects are specifically encoded into a hierarchy of the hierarchical fuzzy anatomy models and subsequently exploited to recognize and delineate the objects in the body region of interest.

5. A computerized method as in claim 2, further comprising automatically determining an optimal hierarchy for the body region that will yield a best recognition and delineation results.

6. A computerized method as in claim 2, wherein using AAR to recognize objects starts from large, well-defined objects and proceeds down a hierarchy of the hierarchical fuzzy anatomy models in a global to local manner.

7. A computerized method as in claim 1, wherein using AAR to recognize and delineate objects in the body region of interest comprises creating a fuzzy model-based version of an Iterative Relative Fuzzy Connectedness (IRFC) delineation algorithm including an affinity function and a seed specification by integrating fuzzy model constraints into a delineation algorithm used to delineate the objects in the body region of interest.

8. A computerized method as in claim 1, wherein using AAR to recognize and delineate objects in the body region of interest comprises applying a one-shot method to roughly segment an organ and to exclude extraneous information in the image.

9. A computerized method as in claim 1, wherein using AAR to recognize and delineate objects in the body region of interest comprises applying a thresholded optimal search to further refine recognition in the image before recognizing and delineating an object in the body region of interest.

10. A computerized method as in claim 1, wherein the body region of interest and the one or more body regions includes one or more of a pelvis, a thorax, an abdomen, or a neck region.

11. A computerized method as in claim 1, further comprising estimating body fat of the patient by measuring an area of fat from a single anatomic slice at a site of the maximum volume-to-area correlation and using the single anatomic slice as a marker.

12. A computerized method as in claim 1, further comprising providing slice localization in a patient's image by linear mapping whereby slice locations are linearly mapped from all patients so that a superior-most and inferior-most anatomic slice locations match in a longitudinal direction for all patients.

13. A computerized method as in claim 1, further comprising providing slice localization in a subject's image by nonlinear mapping whereby in addition to a superior-most and an inferior-most location, a plurality of key landmark locations chosen in a longitudinal direction also match for all subjects.

14. A computerized method as in claim 1, wherein the image is an image of a patient's thorax, further comprising automatically localizing lymph node stations in the body region of interest, comprising the steps of: gathering image data from patient image sets; manually delineating on the patient image sets various lymph node stations using a standard definition of the lymph node stations; building at least one hierarchical fuzzy anatomy model of the lymph node stations from the delineations on the patient image sets; and automatically locating the lymph node stations on the image of the patient's thorax using the at least one hierarchical fuzzy anatomy model.

15. A computerized method as in claim 14, wherein the standard definition of the lymph node stations is based on IASLC standard definitions.

16. A computerized method as in claim 14, wherein automatically locating the lymph node stations further comprises using object recognition algorithms to locate the lymph node stations.

17. A computerized method as in claim 16, wherein the object recognition algorithms comprise a thresholded optimal search algorithm and automatically locating the lymph node stations further comprises using the thresholded optimal search algorithm to refine object pose by an optimal search based on thresholding a test image.

18. The computerized method as in claim 1, wherein the standardized anatomic space comprises locations of imaging slices associated with anatomic landmarks in a direction perpendicular to an imaging plane of the imaging slices.

19. A system, comprising: a database that stores image data from patient image sets; an imaging device that obtains an image of a particular patient body region of interest, wherein slice locations associated with the particular patient body region of interest are mapped from a scanner coordinate system to slice locations in a standardized anatomic space; a memory storing computer instructions; and a processor that processes said computer instructions to build a fuzzy anatomy model of the body region of interest from existing patient image sets, uses automatic anatomy recognition (AAR) to recognize and delineate objects in the particular patient body region of interest, and provides contours of delineated objects, wherein the processor further processes said computer instructions to identify where a volume-to-area correlation becomes maximal for one or more body regions of the patient.

20. A device comprising: one or more processors; and memory storing instructions that, when executed by the one or more processors, cause the device to: build a fuzzy anatomy model of a body region of interest from an existing set of patient images for the body region; obtain an image of a particular patient body region of interest, wherein the image is obtained based on mapping slice locations associated with the particular patient body region of interest from a scanner coordinate system to slice locations in a standardized anatomic space; identify where a volume-to-area correlation becomes maximal for one or more body regions of the patient; use automatic anatomy recognition (AAR) to recognize objects in the particular patient body region of interest and to output spatial information associated with the recognized objects; and delineate, using automatic anatomy recognition (AAR) and based at least on the spatial information associated with the recognized objects, the recognized objects to determine one or more contours or boundaries of the recognized objects.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0032] The present application is further understood when read in conjunction with the appended drawings. For the purpose of illustrating the subject matter, there are shown in the drawings exemplary embodiments of the subject matter; however, the presently disclosed subject matter is not limited to the specific methods, devices, and systems disclosed. In addition, the drawings are not necessarily drawn to scale. In the drawings:

[0033] FIG. 1 illustrates a schematic representation of the AAR schema showing the three main steps of model building, object recognition, and object delineation.

[0034] FIG. 2 illustrates the hierarchies devised for the three body regions studied, where H is a hierarchy, represented as a tree, of the objects in a body region

[0035] FIG. 3 illustrates organs from one training set for each body region displayed via surface rendering. For each row, objects in one picture are listed as { . . . }. Top row: Thorax. 3rd picture: {RPS, TB, E, AS, VS, PC}. Middle row: Abdomen. 3rd picture: {Ask, Lvr, LKd, IVC, AIA, Spl, SAT, Ms1}. Bottom row: Neck. 5th picture: {Mnd, Tng, NP, OP, Ad, FP, Tns1}.

[0036] FIG. 4 illustrates volume renditions of fuzzy models of objects in different combinations for the three body regions. For each row, objects in one picture are listed as { . . . }. Top row: Thorax. 5th picture: {LPS, AS, TB}. Middle row: Abdomen. 3rd picture: {ASk, Lvr, LKd, RKd, AIA, IVC, Spl}. Bottom row: Neck: 5th picture: {Mnd, Tng, NP, OP, Ad, FP}.

[0037] FIG. 5 illustrates volume renditions of fuzzy models created without (Rows 1 and 3) and with (Rows 2 and 4) orientation alignment for several non-sparse (Rows 1 and 2) and sparse (Rows 3 and 4) objects. Row 1: PC, RPS, LKd, Lvr. Row 3: AS, E, AIA, IVC, TB.

[0038] FIG. 6 illustrates sample recognition results for the Thorax for the alignment strategy shown in Equation (10). Cross sections of the model are shown overlaid on test image slices. Left to right: TSkn, TSk, LPS, TB, RPS, E, PC, AS, VS.

[0039] FIG. 7 illustrates sample recognition results for the Abdomen for the alignment strategy shown in Equation (10). Cross sections of the model are shown overlaid on test image slices. Left to right: ASkn, ASk, SAT, Lvr, RKd, LKd, Spl, Msl, AIA, IVC.

[0040] FIG. 8 illustrates sample recognition results for the Neck for the alignment strategy shown in Equation (10). Cross sections of the model are shown overlaid on test image slices. Left to right: NSkn, FP, Mnd, NP (note that NP is a combination of nasal cavity and nasopharynx), Ad, OP, RT, LT, Tng, SP.

[0041] FIG. 9 illustrates the recognition results for the Mill data set DS4 (right) with the model cross-section overlaid on test image slices for ASkn and SAT.

[0042] FIG. 10 illustrates sample delineation results for the Thorax. Left to Right: TSkn, TMS, LPS, AS, RPS, PC, TB, E.

[0043] FIG. 11 illustrates sample delineation results for the Abdomen. Left to Right: ASkn, SAT, Lvr, SAT, RKd, LKd, Spl, Msl, AIA.

[0044] FIG. 12 illustrates sample delineation results for the Neck. Left to Right: NSkn, FP, NP, OP, RT, LT, Tng, SP, Ad.

[0045] FIG. 13 illustrates sample delineation results for DS4 for ASkn (left) and SAT (right).

[0046] FIG. 14 illustrates a schematic representation of the approach of standardized anatomic space in an exemplary embodiment for quantifying abdominal fat.

[0047] FIG. 15 illustrates calibration to create a standard scale. M standard locations selected on three patients are shown on the right. They are mapped linearly on to the standard scale shown (thick) on the left.

[0048] FIG. 16 illustrates non-linear mapping from patient space to standard anatomic space.

[0049] FIG. 17 illustrates anatomic locations of slices in B=the Abdominal Region for 50 subjects. The abscissa shows subject numbers, and the ordinate indicates the extent of B in different subjects in the cranio-caudal direction in terms of the vertebral bodies.

[0050] FIG. 18 illustrates surface renditions of the skeletal components in B of some of the subjects who showed wide variation in vertebral positions in FIG. 17 including subjects 4, 7, 11, 12, 14, 20, 22, 23, 33, 34, 37, and 48.

[0051] FIG. 19 illustrates correlation values from linear mapping (top row) and non-linear mapping (bottom row) for SAT (left) and VAT (right). The vertical axis shows the correlation value. The horizontal axis shows the location of image slices (Slice 1 is at the inferior most position). Some key landmark positions are indicated along the horizontal axis in the bottom row.

[0052] FIG. 20 illustrates anatomic locations (marked with ‘*’) of maximum correlation between (single) slice area and volume (SAT on left, VAT on right, linear method in top row, non-linear method in bottom row). The horizontal axis shows subject numbers, and the vertical axis shows anatomic location from T7 to L5.

[0053] FIG. 21 illustrates anatomic locations where maximum correlation is achieved for the two mapping methods. The top two rows are for SAT (left) and VAT (right) by using linear mapping; and the bottom two rows represent SAT and VAT for non-linear mapping. The 1st and 3rd rows are from the same subject, and the 2nd and 4th rows are from another subject. The spine is used as a reference to show the slice locations (as a white line) of maximum correlation.

[0054] FIG. 22 illustrates correlation curves from non-linear mapping derived from a calibration data set different from that used in FIG. 19. Left: SAT; Right: VAT. The vertical axis shows the correlation value. The horizontal axis shows the location of image slices (Slice 1 is at the inferior most position). Some key landmark positions are indicated along the horizontal axis.

[0055] FIG. 23 illustrates correlation curves from linear and non-linear mapping. Rows 1 and 2: SAT (left) and VAT (right) results for linear mapping using 2 and 3 slices. Rows 3 and 4 are similarly for non-linear mapping. The vertical axis shows the correlation values and the horizontal axis shows the combination number for the different combinations among all possible choices of 2 and 3 slices. The oscillations seen are due to the systematic pattern of multiple slice number combinations.

[0056] FIG. 24 illustrates examples of the mediastinal lymph node stations delineated on axial CT slices from one study (Top row: Stations 1, 2R, 2L, 3p, 3a; Bottom row: Stations 4R, 4L, 5, 6, 7).

[0057] FIG. 25 illustrates hierarchy HA for automatic recognition, where HA has the thoracic skin (tskin) object as root and trachea and bronchi (tb) as offspring of tskin and parent of all stations.

[0058] FIG. 26 illustrates volume rendered displays of the generated fuzzy models showing stations 1, 2 (2R, 2L), 3p, 4 (4R, 4L), and 7 and fuzzy models of the right pleural space (rps) and trachea & bronchi (tb) for anatomic context (Top row (L to R): [rps, stations 1, 2R, 2L, 3p, 4R, 4L, and 7], [stations 1, 2R, 2L, 3p, 4R, 4L, and 7], [stations 2R, 2L, 3p, 4R, 4L, and 7]; Bottom row (L to R): [rps, stations 3p and 7], [station 1 and tb], [rps, stations 2R, 2L, 3p, 4R, 4L, and 7]).

[0059] FIG. 27 illustrates recognition results for three stations by overlaying the cross section of the fuzzy model at recognition on a slice of the test image. The true delineation of the stations is also displayed for reference on the corresponding slice.

[0060] FIG. 28 illustrates a generalized anatomic hierarchy for the thoracic region.

[0061] FIG. 29 illustrates the use of a superset instead of a precise boundary when modeling sparse objects such as the trachea and bronchi or esophagus.

[0062] FIG. 30 illustrates a system in accordance with the invention for providing radiation therapy planning in patients receiving radiation treatment.

DETAILED DESCRIPTION OF ILLUSTRATIVE EMBODIMENTS

[0063] Certain specific details are set forth in the following description with respect to FIGS. 1-30 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.

[0064] In the detailed description to follow, the AAR methodology for automatic anatomy recognition for different body regions will be described. Then, examples will be provided for implementing the AAR methodology to quantify abdominal fat, provide automatic localization of IASLC-defined mediastinal lymph node stations, and radiation therapy planning in exemplary embodiments.

Building Fuzzy Model of Body Region

Notations and Overall Approach

[0065] As described herein, the AAR approach has five unique characteristics: (1) direct hierarchical codification of the prior object geographic and geometric relationship information; (2) a “what-you-see-is-what-you-get” entirely digital fuzzy modeling strategy; (3) hierarchical object recognition strategies that go from a broader gestalt to narrower specifics in locating objects; (4) demonstrated generality of applicability of the same approach to different organ systems, body regions, and modalities; and (5) adaptability of the system to different applications. The AAR approach of the invention is graphically summarized in FIG. 1. As illustrated, the body is divided into body regions B.sub.1, . . . , B.sub.K. Models are built for each specific body region ∈{B.sub.1, . . . , B.sub.K} and each population group G (whatever way G is defined). Throughout this description, and G are treated as variables, and each body region is considered separately and independent of other body regions. The three main blocks in FIG. 1 correspond to model building, object recognition, and object delineation. A fuzzy model FM() is built separately for each of the L objects in , and these models are integrated into a hierarchy chosen for . The output of the first step is a fuzzy anatomic model FAM(, G) of the body region for group G. This model is utilized in recognizing objects in a given patient image I of belonging to G in the second step. The hierarchical order is followed in this process. The output of this step is the set of transformed fuzzy models FM.sup.T() corresponding to the state when the objects are recognized in I. These modified models and the image I form the input to the third step of object delineation which also follows the hierarchical order. The final output is in the form of delineated objects O.sub.1.sup.D, . . . , O.sub.L.sup.D, where each is a binary image.

[0066] The following notation is used herein. G: the population group under consideration. : the body region of focus. O.sub.1, . . . , O.sub.L: L objects or organs of (such as esophagus, pericardium, etc. for =Thorax). ={I.sub.1, I.sub.N}: the set of images of for G from N subjects which are used for model building and for training the parameters of the AAR algorithms. : the binary image representing the true delineation of object in the image I.sub.n∈.sup.b={: 1≤n≤N & 1≤,≤L} is the set of all binary images used for model building. FM(): Fuzzy model of object derived from the set of all binary images ={: 1≤n≤N} of . FAM(,G): Fuzzy anatomy model of the whole object assembly in with its hierarchyFM.sup.T(): Transformed (adjusted) FM() corresponding to the state when is recognized in a given patient image I. : Delineation in I represented as a binary image. Any image I will be represented by a pair I=(C,f), where C denotes a 3D rectangular array of voxels, and f is a mapping f: C.fwdarw. where is a set of integers (except when dealing with fuzzy sets, which are also expressed as images for computational purposes, in which case is a set of real numbers) denoting the image intensities. For any binary image J=(C,f.sub.b), PAS(J) will be used to denote the principal axes system derived from the set X of voxels of J with value 1. PAS(J) is described by the geometric center of Xand the eigenvectors derived from X via principal component analysis.

[0067] The following description will follow the schematic of FIG. 1. Table 1 below lists brief anatomic definitions of all objects from all three body regions considered herein.

TABLE-US-00001 TABLE 1 Anatomic definitions of organs Thoracic objects Acronym Definition of object. Thoracic TSkn The outer boundary of the thoracic skin skin (arms excluded). The interior region constitutes the entire thoracic body region. The inferior boundary is defined to be 5 mm below the base of the lungs and the superior boundary is defined to be 15 mm above the lung apices. Thoracic TSk All skeletal strictures contained in the skeleton thoracic body region, including the spine, ribs, sternum, and the portions of the scapulae and clavicles that are inside the body region. Respiratory RS Grouping of RPS, LPS, and TB. system Right lung RPS The outer boundary of the right lung along the right pleura. Left lung LPS The outer boundary of the left lung along the left pleura. Trachea TB The outer boundary of the trachea and and bronchi bronchi from the superior thoracic trachea to the distal main stem bronchi. Internal IMS Grouping of PC, E, AS, and VS. media- stinum Pericardial PC Region within the boundary of pericardial region sac. The superior aspect is defined by the branching of the main pulmonary artery. Esophagus E The outer boundary of the esophagus from the superior aspect of thorax to the level of gastric cardia. Arterial AS The outer boundary of the ascending aorta, system aortic arch, descending thoracic aorta, pulmonary arteries, innominate artery, proximal left common carotid artery, and proximal left subclavian artery. The superior aspect is defined by the branching of the innominate artery. Venous VS The outer boundary of the superior system vena cava, right and left brachiocephalic veins, and azygos vein. Abdominal Objects Acronym Definition of object Abdominal ASkn The outer boundary of the abdominal skin. The skin interior region constitutes the entire abdominal body region. The superior boundary is defined by the superior aspect of the liver. The inferior boundary is defined by the bifurcation of the abdominal aorta into the common iliac arteries. Abdominal Ask All skeletal structures contained in the skeleton abdominal body region, including lumbar spine and portion of the inferior ribs within the body region. Soft tissue ASTs Grouping of Kd, Spl, Msl, AIA, IVC. Kidneys Kd Grouping of RKd and LKd. Right RKd The outer boundary of the right kidney. kidney All external blood vessels are excluded. Left kidney LKd The outer boundary of the left kidney. All external blood vessels are excluded. Spleen Spl The outer boundary of the spleen. All external blood vessels are excluded. Muscle Msl The outer boundaries of the abdominal musculature, including the rectus abdominis, abdominal oblique, psoas, and paraspinal muscles. Abdominal AIA The outer boundary of the abdominal aorta. aorta The superior and inferior slices of AIA are the same as those of the abdominal region. Inferior IVC The outer boundary of the inferior vena cava. vena cava The superior and inferior slices of IVC are the same as those of the abdominal region. Liver Lvr The outer boundary of the liver. The intrahepatic portal veins and hepatic arteries are included in this region. Fat Fat Grouping of SAT and VAT Sub- SAT Adipose tissue in the subcutaneous cutaneous region in the abdomen. adipose tissue Visceral VAT Adipose tissue internal to the adipose abdominal musculature. tissue Neck objects Acronym Definition of object Head and NSkn The outer boundary of the head and Neck skin neck skin, where the interior region constitutes the entire head and neck body region. The superior boundary is defined by a level 6.6 mm above the superior aspect of the globes. The inferior boundary is defined by a level 6.6 mm inferior to the inferior aspect of the mandible. Air and A&B Grouping of Mnd and Phrs. Bone Mandible Mnd The outer boundary of the mandible. Pharynx Phrx Grouping of NP and OP. Naso- NP The outer contour of the nasal and pharyngeal nasopharyngeal air cavity, extending airway to the inferior aspect of the soft palate. Oro- OP The outer contour of the oropharyngeal pharyngeal air cavities, extending from the inferior airway aspect of the soft palate to the superior aspect of the epiglottis. Fat pad FP The outer boundary of the parapharyngeal fat pad. Neck soft NSTs Grouping of Tnsl, Tng, SP, Ad. tissues Palatine Tnsl Grouping of RT and LT. tonsils Right RT The outer boundary of the right palatine tonsil. palatine tonsil Left LT The outer boundary of the left palatine tonsil. palatine tonsil Tongue Tng The outer boundary of the tongue. Soft palate SP The outer boundary of the soft palate. Adenoid Ad The outer boundary of the adenoid tissue. tissue

Gathering Image Database for and G

[0068] The basic premise of the AAR approach is that the fuzzy anatomic model of for G should reflect near normal anatomy. Consequently, the cleanest way of gathering image data for model building will be to prospectively acquire image data in a well-defined manner from subjects in group G who are certified to be near normal. Such an approach would be expensive and may involve radiation exposure (in case of CT imaging). For developing the concepts and testing the feasibility of AAR, therefore, the inventors have taken a vastly less expensive and simpler approach of utilizing existing human subject image data sets. For the thoracic and abdominal body regions, a board certified radiologist (co-author DAT) selected all image data sets (CT) from the health system patient image database in such a manner that the images appeared radiologically normal for the body region considered, with exception of minimal incidental focal abnormalities such as cysts, small pulmonary nodules, etc. Images with severe motion/streak artifacts or other limitations were excluded from consideration. For these two body regions, the population groups considered have an age range of approximately 50-60 years. This age range was selected to maximize the chances of finding sufficient number of near normal images. For the neck body region, the inventors utilized image data (MR1) previously acquired from normal subjects for the study of pediatric upper airway disorders. Gin this instance is female subjects in the age range of 7-18. The modeling schema is such that the population variables can be defined at any desired “resolution” in the future and the model can then be updated when more data are added.

[0069] Some organs in are better defined in a slice plane different from the slice plane used for imaging others. For example, for =neck, the best plane for slice imaging is sagittal for tongue and soft palate, while for the upper airways and other surrounding organs, axial slices are preferred. The AAR methodology automatically handles organs defined in images with different orientations of digitization by representing image and object data in a fixed and common scanner coordinate system of reference.

Delineating Objects of in the Images in the Database

[0070] There are two aspects to this task—forming an operational definition of both and the organs in in terms of their precise anatomic extent, and then delineating the objects following the definition. These considerations are important for building consistent and reliable models, and, in the future, if similar efforts and results for body-wide models are to be combined, exchanged, and standardized.

[0071] Definition of body regions and objects: Each body region is defined consistently in terms of a starting and ending anatomic location. For axial slice data, these locations are determined in terms of transverse slice positions. For example, for =Thorax, the body region is considered to extend axially from 5 mm below the base of the lungs to 15 mm above the apex of the lungs. Arms are not included in this study. For other orientations of slice planes in slice imaging, the same definitions are applied but translated into other planes. Similarly, each object included in is defined precisely irrespective of whether it is open-ended—because it straddles body regions (for example, esophagus)—or closed and contained within but is contiguous with other objects (for example, liver with hepatic portal vein, common hepatic artery, and bile duct). For each body region, the inventors have created a document that delineates its precise definition and the specification of the components and boundaries of its objects. This document is used as a reference by all involved in generating data sets for model building. These definitions are summarized in Table 1 above.

[0072] Each body region is carved out manually, following its definition, from the data sets gathered for it. In the notation herein, denotes the resulting set of such standard images that precisely cover as per definition. The inventors assume the scanner coordinate system, SCS, as a common reference system with respect to which all coordinates will be expressed.

[0073] Delineation of objects: The objects of are delineated in the images of adhering to their definition, by a combination of methods including live wire, iterative live wire (Souza, A., and Udupa, J. K. , 2006. Iterative live wire and live snake: New user-steered 3D image segmentation paradigms, SPIE Medical Imaging. SPIE, pp. 1159-1165), thresholding, and manual painting, tracing and correction. To minimize human labor and to maximize precision and accuracy, algorithms in terms of a proper combination of these methods and the order in which objects are delineated are devised first, all of which operate under human supervision and interaction. For illustration, in the abdomen, to delineate subcutaneous adipose tissues (SAT) as an object, the skin outer boundary ASkn (as an object) is first segmented by using the iterative live wire method. Iterative live wire is a version of live wire in which once the object is segmented in one slice, the user commands next slice, the live wire then operates automatically in the next slice, and the process is continued until automatic tracing fails when the user resorts to interactive live wire again, and so on. Subsequently, the interface between the subcutaneous and visceral adipose compartments is delineated by using also the iterative live wire method. Once these two object boundaries are delineated, the subcutaneous and visceral components are delineated automatically by using thresholding and morphological operations. On MR images, the same approach works if background non-uniformity correction and intensity standardization (Nyul, L.G., Udupa, J. K., 1999. On standardizing the MR image intensity scale. Magnetic Resonance in Medicine 42, 1072-1081) are applied first to the images in . If direct delineation by manual tracing or even by using live wire is employed, the process would become complicated (because of the complex shape of the adipose and visceral compartments) and much more labor intensive.

[0074] Because of the enormity of this task, a number of trainees, some with medical and biomedical but some with engineering background, were involved in accomplishing this task. All tracings were examined for accuracy by several checks—3D surface renditions of objects from each subject in various object combinations as well as a slice-by-slice verification of the delineations overlaid on the gray images for all images. The set of binary images generated in this step for all objects is denoted by .sup.b={: 1≤n≤N & 1≤<L}. The set of binary images generated just for object is denoted by ={: 1≤n≤N}.

Constructing Fuzzy Object Models

[0075] The Fuzzy Anatomy Model FAM(, G) of any body region for group G is defined to be a quintuple:

FAM(,G)=(H,,ρ,λ,η).   (1)

Briefly, the meaning of the five elements of FAM(, G) is as follows. H is a hierarchy, represented as a tree, of the objects in □ as illustrated in FIG. 2. is a collection of fuzzy models, one model per object in . ρ describes the parent-to-offspring relationship in H over G. λ is a set of scale factor ranges indicating the size variation of each object over G. η represents a set of measurements pertaining to the objects in . A detailed description of these elements and the manner in which FAM(, G) is derived from and are presented below.

[0076] Hierarchy H: This element describes the way the objects of are considered ordered anatomically as a tree structure. This order currently specifies the inclusion of an offspring object O.sub.k anatomically in the parent object O.sub.l though other arrangements are possible for H. While each has its own hierarchy, itself forms the offspring of a root denoting the whole body, WB, as shown in FIG. 2. FIG. 2(a) illustrates the hierarchy for whole body WB, while FIG. 2(b) illustrates the hierarchy for the Thorax, where TSkn: Outer boundary of thoracic skin as an object; RS: Respiratory System; TSk: Thoracic Skeleton; IMS: Internal Mediastinum; RPS, LPS: Right & Left Pleural Spaces; TB: Trachea & Bronchi; E: Esophagus; PC: Pericardium; AS, VS: Arterial & Venous Systems. FIG. 2(c) illustrates the hierarchy for the Abdomen, where ASkn: Outer boundary of abdominal skin; ASk: Abdominal Skeleton; Lvr: Liver; ASTs: Abdominal Soft Tissues; SAT & VAT: Subcutaneous and Visceral Adipose Tissues; Kd: Kidneys; Spl: Spleen; Msl: Muscle; AIA: Aorta and Iliac arteries; IVC: Inferior Vena Cava; RKd & LKd: Right and Left Kidneys. FIG. 2(d) illustrates the hierarchy for the Neck, where NSkn: Outer boundary of skin in neck; A&B: Air & Bone; FP: Fat Pad; NSTs: Soft Tissues in neck; Mnd: Mandible; Phrx: Pharynx; Tnsl: Tonsils; Tng: Tongue; SP: Soft Palate; Ad: Adenoid; NP & OP: Nasopharynx and Oropharynx; RT & LT: Right and Left Tonsils.

[0077] An object that is exactly a union of its offspring will be referred to herein as a composite object. Examples: RS, Fat, Kd, etc. Note that none of the skin objects is a composite object since the full body region inside the skin is not fully accounted for by the union of the offspring objects. The notion of composite objects is useful in combining objects of similar characteristics at a higher level of the hierarchy, which may make object recognition (and delineation) more effective. Thin tubular objects will be called sparse objects: TB, E, AS, VS, AIA, IVC, Phrx, NP, and OP. Compact, blob-like objects will be referred to as non-sparse: TSkn, RS, IMS, LPS, RPS, PC, ASkn, Fat, SAT, VAT, Lvr, Spl, Kd, RKd, LKd, NSkn, FP, NSTs, Tnsl, Tng, SP, Ad, RT, and LT. Some objects are a hybrid between these two types, consisting of both features. Examples: TSk, Ask, ASTs, A&B, and Mnd.

[0078] Fuzzy model set : The second element (in the description of FAM(, G) represents a set of fuzzy models, ={FM(): 1≤<L}, where FM() is expressed as a fuzzy subset of a reference set ⊂Z.sup.3 defined in the SCS; that is, FM()=(). The membership function (v) defines the degree of membership of voxel v∈ in the model of object . Ideally, for any , 1≤<L, the inventors would like the different samples of in different subjects to differ by a transformation involving translation, rotation, and isotropic scaling. The idea behind the concept of the fuzzy model of an object is to codify the spatial variations in form from this ideal that may exist among the N samples of the object as a spatial fuzzy set, while also retaining the spatial relationship among objects in the hierarchical order.

[0079] Given the training set of binary images of object , the inventors determine and FM() for as follows. The inventors permit only such alignment operations, mimicking among the members of that are executed precisely without involving search and that avoid the uncertainties of local optima associated with optimization-based full-fledged registration schemas. In this spirit, the inventors handle the translation, rotation, and scaling components of in the following manner.

[0080] For translation and rotation, for each manifestation of in , the inventors determine, within SCS, the principal axes system PAS() of . Subsequently, all samples are aligned to the mean center and principal axes (though orientation alignment is often not necessary). The scale factor estimation is based on a linear size estimate (in mms) of each sample of and resizing all samples to the mean size. The size of in is determined from √{square root over (e.sub.1+e.sub.2+e.sub.3)}, where e.sub.1, e.sub.2, and e.sub.3 are the eigenvalues corresponding to the principal components of in

[0081] After aligning the members of via a distance transform is applied to each transformed member for performing shape-based interpolation (Raya, S. P., Udupa, J. K., 1990. Shape-Based Interpolation of Multidimensional Objects. IEEE Transactions on Medical Imaging 9, 32-42; and Maurer, C. R., Qi, R. S., Raghavan, V., 2003. A linear time algorithm for computing exact Euclidean distance transforms of binary images in arbitrary dimensions. IEEE Transactions on Pattern Analysis and Machine Intelligence 25, 265-270). The distances are averaged over all members and converted through a sigmoid function to obtain the membership values and subsequently FM().

[0082] Parent-to-offspring relationship ρ: This element describes the parent-to-offspring spatial relationship in H for all objects in . Since each object O.sub.k has a unique parent, this relationship is represented by ρ={ρ.sub.k: 1≤k≤L} (which also encodes WB to body region relationships, although this is not taken into account in the current implementation). For each O.sub.k, ρ.sub.k codifies the mean position as well as the orientation relationship between O.sub.k and its parent over N samples. The inventors adopt the convention that ρ.sub.i denotes the relationship of the root object of relative to SCS. Let G be the geometric center of in . Then, the mean positional relationship between and O.sub.k is considered to be the mean of the vectors in the set {GC.sub.n,k-G1≤n≤N}. To find the mean orientation .sub.k, the inventors make use of the eigenvectors , , and of the shape of in estimated over all N samples. The inventors take an average of each , over N samples for i=1, 2, 3. However, for some n and i, may be more than 90 degrees from the average, in which case is replaced by − while simultaneously replacing by − for some j different from i so as to keep the system right-handed. The inventors then recalculate the average, and repeat until the eigenvector is within 90 degrees of the average. Then, starting from either the first or the third eigenvector, whichever has the eigenvalue farther from the second, the inventors normalize and make the others orthogonal to it. is then taken to be the transformation that aligns the eigenvector system of the parent with that mean orientation. This method guarantees a robust orientation estimate despite the 180-degrees switching property of eigenvectors.

[0083] In order not to corrupt ρ.sub.k by the differences in size among subjects, before estimating ρ.sub.k, the parent and all offspring objects O.sub.k of are scaled with respect to the center G of as per a common scale factor, estimated for via the method described above. The reasoning behind this scaling strategy is that an object and its entire offspring should be scaled similarly to retain their positional relationship information correctly.

[0084] Scale range λ: The fourth element λof FAM(, G) is a set of scale factor ranges, λ={=[]: 1≤≤L}, indicating the size variation of each object over its family . This information is used in recognizing in a given image to limit the search space for its pose, as explained in the description of recognizing objects below.

[0085] Measurements η: This element represents a set of measurements pertaining to the object assembly in . Its purpose is to provide a database of normative measurements for future use. This element also serves to improve knowledge about object relationships (in form, geographical layout, etc. in ) and thence in constructing better hierarchies for improving AAR, as explained briefly below.

[0086] There are several parameters related to object recognition (next section) and delineation (see below), some of which are image modality specific (They are identified by .sup.m and in the next section and σ.sub.Ψo, m.sub.Ψo, m.sub.ΨB, m.sub.ΨO and σ.sub.ΨB in the following section.) The values of these parameters are also considered part of the description of η. The definition of these parameters and the process of their estimation are described at relevant places in the following sections for ease of reading, although their actual estimation is done at the model building stage.

[0087] The fuzzy anatomy model FAM(, G) output by the model building process is used in performing AAR on any image I of for group G as described in the following sections.

Recognizing Objects

[0088] The process of what is usually referred to as “segmenting an object in an image” may be thought of as consisting of two related phenomena—object recognition (or localization) and object delineation. Recognition is a high-level process of determining the whereabouts of the object in the image. Given this information for the object, its delineation is the meticulous low-level act of precisely indicating the space occupied by the object in the image. The design of the entire AAR methodology is influenced by this conceptual division. The inventors believe that without achieving acceptably accurate recognition it is impossible to obtain good delineation accuracy. The hierarchical concept of organizing the objects for AAR evolved from an understanding of the difficulty involved in automatic object recognition. Once good recognition accuracy is achieved, several avenues for locally confined accurate delineation become available, as discussed in the next section. The goal of recognition in AAR is to output the pose (translation, rotation, and scaling) of FM(), or equivalently the pose-adjusted fuzzy model FM.sup.T(), for each in a given test image I of such that FM.sup.T() matches the information about present in I optimally.

[0089] The recognition process proceeds hierarchically as outlined in the procedure AAR-R presented below. In Step R.sub.1, the root object is recognized first by calling algorithm R-ROOT. The inventors assume that the field of view in I fully encloses the root object. For the hierarchies shown in FIG. 2, the root object is the skin outer boundary which is typically more-or-less, although not perfectly, fully included within the imaging field of view. Proceeding down the tree represented by H in the breadth-first order, other objects are recognized by calling algorithm R-OBJECT. The latter makes essential use of the parent fuzzy model and the parent-to-offspring relationship ρ encoded in FAM(, ).

TABLE-US-00002 Procedure AAR-R Input:  An image I of    , FAM(   , G). Output:  FM.sup.T (   ),    = 1, . . . , L. Begin R1. Call R-ROOT to recognize the root object in H; R2. Repeat R3.  Find the next offspring O.sub.k to recognize in H (see text), R4.  Knowing FM.sup.T (   ), ρ.sub.k, and λ.sub.k, call R-OBJECT to recognize O.sub.k; R5. Until all objects are covered in H; R6. Output FM.sup.T (   ),    = 1, . . . , L; End

[0090] Two strategies are described here for each of algorithms R-ROOT and R-OBJECT. The first, a global approach, does not involve searching for the best pose. The inventors call this the One-Shot Method since the model pose is determined directly by combining the prior information stored in FAM(, ) and information quickly gathered from the given image I. The one-shot method is used as initialization for a more refined second method called Thresholded Optimal Search.

One-Shot Method

[0091] A threshold interval Th.sub.1 corresponding to the root object O.sub.1 is applied to I followed by a morphological opening operation to roughly segment O.sub.1 to produce a binary image J. The purpose of the morphological operation is to exclude as much as possible any significant extraneous material, such as the scanner table and patient clothing, from J. Then the transformed model FM.sup.T(O.sub.1) is found by applying a transformation .sup.m to FM(O.sub.1). .sup.m is devised to express the mean relationship between the roughly segmented O.sub.1 and the true segmentation of O.sub.1 represented in the binary images I.sub.n,1∈.sup.b. The estimation of .sup.m is done at the model building stage of AAR as mentioned above. To determine .sup.m, similar thresholding and morphological operations are performed on each gray image I.sub.n, in the training set to obtain a rough segmentation of O.sub.1, denoted J.sub.n,1, in I.sub.n. The relationship between this rough segmentation J.sub.n,1 and the true segmentation I.sub.n,1 of O.sub.1 in .sup.b is found as a transformation that maps PAS(J.sub.n,1) to PAS(I.sub.n,1). The mean, denoted .sup.m, of such transformations over all training images is then found.

[0092] Once the root object O.sub.1 is recognized, the poses for other objects in I in the hierarchy H are determined by combining (in the sense of composition) .sup.m with the parent to offspring relationship information stored in pk for each parent-offspring pair. The transformed models FM.sup.T() are then found from this information.

Thresholded Optimal Search

[0093] [This is a strategy to refine the results obtained from the one-shot method. Its premise is that the overall image intensity of the objects in can be characterized by threshold intervals such that, at the model's pose corresponding to the best match of the model with an underlying object in the given test image I, the mismatch between the thresholded result and the model is minimal. For MR images for this approach to make sense, it is essential to correct for background intensity non-uniformities first followed by intensity standardization (Nyul and Udupa 1999).

[0094] Suppose that at the model building stage, the optimal threshold interval for each object has already been determined automatically from the training image set. It will be explained below as to how this is accomplished. Then, at the recognition stage, the threshold for is fixed at this learned value . Starting from the initial pose found by the one-shot method, a search is made within the pose space for an optimal pose p* of the fuzzy model over I that yields the smallest sum of the volume of false positive and false negative regions, where the model itself is taken as the reference for defining false positive and negative regions. Specifically, let FM.sup.P() denote the fuzzy model of at any pose p, expressed as an image, and let J denote the binary image resulting from thresholding I at . Then:

[00001]$\begin{array}{cc}{p}^{*}\begin{array}{c}\in \arg \min \\ p\end{array}\left(.Math.{\mathrm{FM}}^P\left(O_{\ell }\right)-J.Math.+.Math.J-F{M}^P\left(O_{\ell }\right).Math.\right).& \left(2\right)\end{array}$

Since arg min is a set, “∈” means one of the values chosen from the set is assigned to p*. Image subtraction here is done in the sense of fuzzy logic, and |x| denotes the fuzzy cardinality of x, meaning that it represents the sum total of the membership values in x. The search space to find p* is limited to a region around the initial pose. This region is determined from knowledge of ρ.sub.k and its variation and the scale factor range λ.sub.k. For the positional vector, the inventors search in an ellipsoid with its axes in the coordinate axis directions and with length four times the standard deviation of the corresponding coordinate. When searching in orientation space, the inventors search in an ellipsoid with its axes in the direction of the eigenvectors of the rotation vector distribution (covariance matrix) and with length four times the square root of the corresponding eigenvalue. (A rotation vector has magnitude equal to the angle of rotation and direction along the axis of right-handed rotation. The rotation referred to is the rotation of required to bring it into coincidence with .) For the scale factor, the inventors search in an interval of size four times the standard deviation of the scale factor.

[0095] Determining The at the model building stage: To estimate , the inventors run a rehearsal of the recognition method described above as follows, essentially for attempting to learn the recognition process. Imagine the inventors already built and estimated ρ and λ. Suppose that the recognition process is now run on the training images. Since the optimal threshold is not known but the true segmentations are known, the idea behind this learning of the recognition process is to test recognition efficacy for each of a number of threshold intervals t and then select the interval that yields the best match of the model with the known true segmentations for each . That is, if J.sub.n(t) is the binary image resulting from thresholding the training image I.sub.n at t, then:

[00002]$\begin{array}{cc}T{h}_{\ell }\in \begin{array}{c}\arg \min \\ p,t\end{array}\underset{n}{.Math.}.Math.\left(J_n\left(t\right)\times F{M}^p\left(O_{\ell }\right)\right)-I_{n,\ell }.Math.+.Math.I_{n,\ell }-\left(J_n\left(t\right)\times F{M}^p\left(O_{\ell }\right)\right).Math..& \left(3\right)\end{array}$

Here, × denotes fuzzy intersection. In words, the optimal threshold is found by searching over the pose space over all training data sets and all thresholds the best match between the true segmentation of with the result of thresholding I.sub.n restricted to the model. In the present implementation, 81 different values of the intervals are searched (9 for each end of the interval). The 9 positions for the lower end are the 5.sup.th, 10.sup.th, . . . , 45.sup.th percentile values of the cumulative object intensity histogram determined from the training image set. Similarly, for the upper end, the positions are 55.sup.th to 95.sup.th percentile values.

[0096] To summarize, the thresholded optimal search method starts the search process from the initial pose found by the one-shot method. It uses the optimal threshold values determined at the training stage for each object and finds the best pose for the fuzzy model of in the given image I by optimally matching the model with the thresholded version of I. The only parameters involved in the entire recognition process are the thresholds one threshold interval per object, and .sup.m. Their values are automatically determined in the model building stage from image and binary image sets and and they become part of the model FAM(, G) itself.

Delineating Objects

[0097] Once the recognition process is completed and the adjusted models FM.sup.T() are output for a given image I of , delineation of objects is performed on I in the hierarchical order as outlined in the procedure AAR-D. As in recognition, in Step D1, the root object is first delineated by calling D-ROOT AAR-D then proceeds in the breadth-first order to delineate other objects by calling D-OBJECT.

TABLE-US-00003 Procedure AAR-D Input:  An image I of    , FAM(   , G) FM.sup.T(   ),    = 1, . . . , L; Output: .sup.   .sup. D,    = 1, . . . , L. Begin D1 Call D-ROOT to recognize the root object in H; D2. Repeat D3.  Traverse H and find the next offspring O.sub.k to deliniate in H; D4.  Knowing delineation of    , call D-OBJECT to delineate O.sub.k in  I; D5. Until all objects are covered in H; D6. Output   .sup. D,    = 1, . . . , L; End

[0098] For D-ROOT and D-OBJECT, the inventors have chosen an algorithm from the fuzzy connectedness (FC) family in view of the natural and intimate adaptability of the FC methods to prior information coming in the form of fuzzy sets. In particular, since the inventors focus on the problem of delineating one object at a time, for both Steps D1 and D4, the inventors have selected the linear-time Iterative Relative FC (IRFC) algorithm of Ciesielski, K. C., Udupa, J. K., Falcao, A. X., Miranda, P. A. V., 2012. Fuzzy Connectedness Image Segmentation in Graph Cut Formulation: A Linear-Time Algorithm and a Comparative Analysis. Journal of Mathematical Imaging and Vision 44, 375-398, for separating each object from its background. The novel adaptations are in incorporating fuzzy model information into the IRFC formulation and in making the latter fully automatic. These modifications are described below.

Fuzzy Model-Based IRFC (FMIRFC)

[0099] There are two aspects that need to be addressed to fully describe the FMIRFC algorithm: affinity function and seed specification. Affinity is a local concept indicating the degree of connectedness of voxels locally in terms of their spatial and intensity nearness. In the FC family, this local property is grown into a global phenomenon of object connectedness through the notion of path strengths.

[0100] Affinity function: The FC framework (Udupa and Samara sekera 1996, Ciesielski et al. 2012) is graph-based. An ordered graph (C, α) is associated with the given image I=(C,f) where α is an adjacency relation on C such as 6-, 18-, or 26-adjacency. Each ordered pair (c, d) of adjacent voxels in α is assigned an affinity value k(c, d) which constitutes the weight assigned to are (c, d) in the graph. To each path π in the graph (or equivalently in I) in the set of all possible paths π.sub.a,b between two voxels a and b of C, a strength of connectedness K(π) is determined, which is the minimum of the affinities along the path. The connectivity measure K*(a, b) between a and b is then defined to be K*(a,b)=max{K(π): π∈πa,b}. The notion of connectivity measure can be generalized to the case of “between a set A and a voxel b” by a slight modification: K*(A,b)=max{K(π): π∈π.sub.a,b & a∈A}. By using a fast algorithm to compute K*(A, b), the machinery of FC allows a variety of approaches to define and compute “objects” in images by specifying appropriate affinity functions and seed sets. In particular, in IRFC, two seed sets A.sub.O and A.sub.B are indicated for an object O and its background B, respectively. Then the object indicated by A.sub.O is separated optimally from the background indicated by A.sub.B by an iterative competition in connectivity measure between A.sub.O and every voxel c∈C and between A.sub.B and c. In published IRFC methods, A.sub.O and A.sub.B are specified usually with human interaction.

[0101] In FMIRFC, affinities κ.sub.o(c, d) and κ.sub.B(c, d) for O and B are designed separately. Subsequently, they are combined into a single affinity κ by taking a fuzzy union of κ.sub.o and κ.sub.B. Each of κ.sub.o and η.sub.B has three components. The description below is for κ.sub.o. The same applies to κ.sub.B.

κ.sub.o(c,d)=w.sub.1Ψ.sub.0(c,d)+w.sub.2φ.sub.0(c,d)+w.sub.2γ.sub.0(c,d).   (4)

Here, Ψ.sub.0(c, d) represents a homogeneity component of affinity, meaning, the more similar image intensities f(c) and f(d) are at c and d, the greater is this component of affinity between c and d. As commonly done in the FC literature, the inventors set

[00003]$\begin{array}{cc}{\psi }_O\left(c,d\right)=\exp \left[-{\left(f\left(c\right)-f\left(d\right)\right)}^2/2{\sigma }_{{\psi }_O}^2\right],& \left(5\right)\end{array}$

where σ.sub.Ψo is a homogeneity parameter that indicates the standard deviation of intensities within object O. φ.sub.o(c, d), the object feature component, on the other hand, describes the “degree of nearness” of the intensities at c and d to the intensity m.sub.φo expected for the object O under consideration. Denoting the standard deviation of object intensity by σ.sub.Ψ0, this nearness is expressed by:

[00004]$\begin{array}{cc}{\phi }_o\left(c,d\right)=\exp \left[-(\max \left\{{\left(f\left(c\right)-m_{{\phi }_o}\right)}^2,{\left(f\left(d\right)-m_{{\phi }_o}\right)}^2\right\}/2{\sigma }_{{\phi }_o}^2\right].& \left(6\right)\end{array}$

The third component γ.sub.o incorporates fuzzy model information into affinity by directly taking the larger of the two fuzzy model membership values μ.sub.o(c) and μ.sub.o(d) at c and d for the object:

γ.sub.o(c,d)=max{μ.sub.o(c),μ.sub.o(d)}.   (7)

Finally, a combined single affinity κ on I is constructed by

κ(c,d)=max{κ.sub.o(c,d),κ.sub.B(c,d)}.   (8)

[0102] The weights in Equation (4) are chosen equal and such that they add up to 1. The homogeneity parameter is set equal for object and background (σ.sub.Ψo=σ.sub.ΨB) and estimated from uniform regions in the training images (after leaving out high gradient regions), as commonly done in the FC literature (Saha and Udupa 2001). The remaining parameters (σ.sub.Ψo, σ.sub.ΨB, m.sub.Ψo, m.sub.ΨB) are estimated automatically from the training data sets from the knowledge of O and B regions for each object.

[0103] Seed specification: Seed sets A.sub.O and A.sub.B are found by a joint criterion of a threshold for image intensity and for model membership for each of O and B. The threshold interval Th.sub.o for O is the same as the one used for recognition, namely The threshold interval Th.sub.B for background is a union of similar threshold intervals for the background objects. (In principle, all objects other than O can be considered to be background objects of O; however, in practice, only the anatomic neighbors of O matter.) The only new parameters are Th.sub.O.sup.M and Th.sub.B.sup.M used as model thresholds for indicating A.sub.O and A.sub.B , respectively. These parameters are used as follows:

TABLE-US-00004 Algorithm FMIRFC Input'.  Image I of    , FAM(   , G), FM.sup.T(   ) at recognition. Below, assume O =    . Output     . Begin FC1. Determine background     of O; FC2. Retrieve affinities κ.sub.O and κ.sub.B from FAM(   , G); FC3. Compute combined affinity κ; FC4. Retrieve thresholds Th.sub.O, Th.sub.B, Th.sub.O.sup.M, and Th.sub.B.sup.M from FAM(   , G), and determine seed sets A.sub.O and A.sub.B in I via (9); FC5. Call the IRFC delineation algorithm with κ, A.sub.O, A.sub.B, and I as arguments; FC6. Output image    .sup.D returned by the IRFC algorithm; End

A.sub.O={v∈C:f(v)−Th.sub.O & μ.sub.O(v)∈Tho.sup.M},

A.sub.B={v∈C:f(v)∈Th.sub.B & μ.sub.B(v)∈Th.sub.B.sup.M}.   (.sup.9)

In an exemplary implementation, Th.sub.O.sup.M is fixed at [0, 0.9] and [0, 0.5] for non-sparse and sparse objects, respectively, and ThB.sub.M is set to [0, 0]. The FMIRFC algorithm is summarized above.

Illustrations, Experimental Results and Discussion

Image Data

[0104] The data sets used for the three body regions are summarized in Table 2. Data sets DS1 and DS2 are from CT and are selected from a hospital patient image database, and were verified to be of acceptable quality and radiologically normal, with exception of minimal incidental focal abnormalities, in the body regions for which they are chosen. Note the typical clinical resolution for pixel size (˜1 mm) and slice spacing (5 mm) in these data sets and hence the challenge for object recognition and delineation. The goal in focusing on these data was to challenge the AAR system to perform on typical clinical data sets. DS3 is from an on-going research project investigating the association of Polycystic Ovary Syndrome with Obstructive Sleep Apnea in obese pediatric female subjects (Arens, R., Sin, S., Nandalike, K., Rieder, J., Khan, U. I., Freeman, K., Wylie-Rosett, J., Lipton, M. L., Wootton, D. M., McDonough, J. M., Shifteh, K., 2011. Upper Airway Structure and Body Fat Composition in Obese Children with Obstructive Sleep Apnea Syndrome. American Journal of Respiratory and Critical Care Medicine 183, 782-787). It consists of both axial and sagittal acquisitions and a mix of T1- and T2-weighted images. DS1-DS3 represent the three body regions for which the hierarchy of organs was depicted in FIG. 2. DS4 (Wagshul, M. E., Sin, S., Lipton, M. L., Shifteh, K., Arens, R., 2013. Novel retrospective, respiratory-gating method enables 3D, high resolution, dynamic imaging of the upper airway during tidal breathing. Magn Reson Med. Doi:10.1002/mrm.24608), however, is used for testing the ability of the AAR method to rapidly prototype an application by using existing models for the same body region. In this case, models built from DS2 from CT are deployed on DS4 from MRI.

TABLE-US-00005 TABLE 2 Summary of data sets used in the experiments. Data Body Group Number of Image Identifier Region     G (age) Subjects N Modality Imaging Protocol Details Image Information DS1 Thorax 50-60 50 CT Contrast-enhanced, axial, 512 × 512 × 51-69, male normal breath-hold 0.9 × 0.9 × 5 mm.sup.3 DS2 Abdomen 50-60 50 CT Contrast-enhanced, axial, 512 × 512 × 51-69, male normal breath-hold 0.9 × 0.9 × 5 mm.sup.3 DS3 Neck 8-17 15 MRI T-2 weighted, axial & T1- & 400 × 400 × 35-50, male & normal T2-weighted sagittal. T2: 0.5 × 0.5 × 3.3 mm.sup.3 female TR/TE = 8274.3/82.6 msec, T1: TR/TE = 517.7/7.6 msec DS4 Abdomen 8-17 14 MRI T2-weighted, axial. 400 × 400 × 45-5, male & 6 normal, TR/TE = 1556.9/84 msec 0.7 × 0.7 × 6 mm.sup.3 female 8 obese patients

[0105] In all data sets, any extra slices falling outside the body region as per definition are removed manually first. Note the variation in the size of the body region in Table 2 (expressed roughly as slice spacing×number of slices). In the case of MRI, the resulting images are processed, first to suppress background non-uniformities and subsequently to standardize the image intensities (Nyul and Udupa 1999). Standardization is a post-acquisition image processing technique which significantly minimizes the inter-subject and intra- and inter-scanner image intensity variations for the same tissue and achieves tissue-specific numeric meaning for MR images. It has been shown to significantly improve the accuracy of delineation algorithms (Zhuge, Y., Udupa, J. K., 2009. Intensity standardization simplifies brain MR image segmentation, Computer, Vision and Image Understanding 113, 1095-1103). It is done separately for each MRI protocol and body region. For DS1 and DS2, one half of the image data sets were used for model building, which included the estimation of the parameters of the recognition and delineation algorithms (.sup.m, , σ.sub.ΨO, m.sub.φO, m.sub.φB, σ.sub.φO, and σ.sub.φB), and the remaining data sets were used for testing the methods. For DS3, the train-test sets were set up as 11 and 4, and this was repeated 30 times for different choices of 11 and 4 data sets. For DS4, all data sets were used for testing, and model building was based on one half of the data sets in DS2. This provided an interesting scenario for the challenge for the AAR method, in that, models built from normal CT data sets for one patient group were used for performing AAR on MRI data sets from normal subjects and patients from another group.

Model Building

[0106] In FIG. 3, the organs defined in the image of one of the subjects employed in model building are displayed for each body region in different combinations of the organs. The inventors have examined all data sets under DS1-DS3 in this manner which has helped in a proper understanding of the organ relationships. This is important for devising effective hierarchies, recognition strategies, and delineation algorithms.

[0107] FIG. 4 displays fuzzy models FM() of objects in various combinations for the three body regions. Since the volumes are fuzzy, they are volume rendered by using an appropriate opacity function. Note that although the models appear blurred, they portray the overall shape of the objects they represent and the object relationships. From consideration of the difficulties in model building, recognition, and delineation, the inventors divided objects in the body into sparse, non-sparse, and hybrid groups. Sparse objects pose special challenges for recognition and delineation, stemming mostly from difficulties in model building. Variations in the form, shape, and orientation of sparse objects cause them to overlap far less, or often not at all, compared to non-sparse objects, when forming the model by gathering fuzzy overlap information. In other words, the models tend to diffuse or become too fuzzy. For example, in AS (thorax), the descending aortic portion extends from superior to inferior. However, this part is often either bent from the vertical or is crooked, and the pattern of the brachiocephalic and subclavian arteries arising from the aortic arch is different. If the variation is just in orientation only, then aligning by orientation may produce sharper models. But the issue is not one of producing less fuzzy models but of building models that have the right/correct amount of fuzziness so that the recognition process will be least misguided by the model. This dilemma of the disconnection between model building and recognition is common to all model/atlas-based methods and is the real challenge in automatic recognition of sparse and hybrid objects.

[0108] To study the effect of orientation alignment, FIG. 5 displays models created without and with orientation adjustment, for several sparse as well as non-sparse objects from all three body regions. The volume renditions were created with exactly the same settings for each object for its two versions of models. Orientation adjustment does not produce any dramatic difference in the models created, although close scrutiny reveals that the model definition improves slightly; examine especially LPS, AIA, AS, and Lvr.

[0109] Relating to the fifth element η of FAM(, G), Tables 3-5 show correlations among objects in their size for the three body regions. Object size is determined as explained above. As may be expected, bilateral organs, such as LPS and RPS, LKd and RKd, and LT and RT, are strongly correlated in size. That is, their sizes go together, whatever way they may be related to the subject's body size. There are also other interesting strong, poor (or no), and even weak negative, correlations, as highlighted in the tables; for example, TSk with RS and RPS; VS with TB, PC, and E; ASkn with ASTs, SAT and Msl; ASTs with SAT and Msl; Msl with SAT; NSkn with A&B; Ad with NSkn, FP, NP, and SP. Although the inventors have not explored the utility of such information herein, those skilled in the art will appreciate that this and other information will be useful in devising hierarchies more intelligently than guided by just anatomy, and hence in building better FAM(, G).

TABLE-US-00006 TABLE 3 Size correlation among objects of the Thorax. TSkn RS TSk IMS RPS TB LPS PC E AS VS TSkn 1 RS 0.76 1 TSk 0.76 0.93 1 IMS 0.48 0.76 0.71 1 RPS 0.6 0.92 0.88 0.75 1 TB 0.06 0.41 0.5 0.56 0.59 1 LPS 0.64 0.93 0.87 0.74 0.96 0.57 1 PC 0.47 0.51 0.45 0.65 0.28 0.11 0.3 1 E 0.42 0.65 0.56 0.58 0.72 0.58 0.78 0.18 1 AS 0.44 0.53 0.49 0.71 0.54 0.24 0.51 0.35 0.35 1 VS 0.3 0.31 0.35 0.34 0.34 0.09 0.34 −.01 0.05 0.42 1

TABLE-US-00007 TABLE 4 Size correlation among objects of the Abdomen. ASkn ASk ASTs Lvr SAT Msl Spl RKd LKd AIA IVC ASkn 1 Ask 0.68 1 ASTs 0.9 0.8 1 Lvr 0.61 0.48 0.58 1 SAT 1 0.69 0.92 0.61 1 Msl 0.91 0.79 0.99 0.63 0.94 1 Spl 0.62 0.43 0.61 0.51 0.65 0.62 1 RKd 0.53 0.64 0.57 0.61 0.51 0.6 0.34 1 LKd 0.53 0.56 0.52 0.51 0.49 0.54 0.34 0.87 1 AIA 0.6 0.85 0.7 0.27 0.58 0.68 0.49 0.51 0.5 1 IVC 0.32 0.58 0.47 0.29 0.32 0.46 0.3 0.38 0.36 0.67 1

Object Recognition

[0110] Results for recognition are summarized in FIGS. 6-8 and Tables 6-9 for the different body regions. FIGS. 6-8 and Tables 6-8 illustrate recognition results for the three body regions for the best set up involving orientation adjustment selectively for different objects. The alignment strategy was as follows for the different objects in these results. [0111] Non-sparse & hybrid objects: RS, LPS, RPS, TMS, TSk, Ask, Kd, Spl, Msl, LKd, RKd, A&B, FP, NSTs, Mnd, Tnsl, Tng, SP, Ad, RI, LT—no orientation alignment. [0112] Sparse objects: TB, E, AS, VS, AIWA, IVC, Phrx, NP, OP — orientation alignment by all axes.

Equation (10)

[0113] The recognition accuracy is expressed in terms of position and size. The position error is defined as the distance between the geometric centers of the known true objects

TABLE-US-00008 TABLE 5 Size correlation among objects of the Neck. NSkn A&B FP Mnd NP OP Tng SP Ad LT RT NSkn 1 A&B 0.89 1 FP 0.76 0.81 1 Mnd 0.75 0.96 0.83 1 NP 0.39 0.12 −.06 −.12 1 OP 0.63 0.59 0.44 0.54 0.14 1 Tng 0.83 0.75 0.76 0.66 0.19 0.65 1 SP 0.5 0.27 0.23 0.14 0.46 0.26 0.37 1 Ad −.2 0.61 −.19 0.1 −.29 −.06 −.07 −.19 1 LT 0.61 0.56 0.58 0.48 0.28 0.5 0.64 0.25 −.1 1 RT 0.61 0.56 0.58 0.48 0.28 0.5 0.64 0.25 −.1 1 1
in and the center of the adjusted fuzzy model FM.sup.T(). The size error is expressed as a ratio of the estimated size of the object at recognition and true size. Values 0 and 1 for the two measures, respectively, indicate perfect recognition. Note in FIGS. 6-8 that the model bleeds into adjacent tissue regions with some membership value since it is fuzzy. This should not be construed as wrong segmentation. The main issue is if the model placement via recognition is accurate enough to yield good delineation. Similarly and due to the slice visualization mode, sparse object components may appear to be missed or to present with low membership values.

[0114] Although the inventors have not conducted extensive experiments to test all possible arrangements for orientation alignment for non-sparse and sparse objects, the inventors generally found that orientation adjustment for non-sparse objects does not improve recognition results. In some cases like PC, it may actually lead to deterioration of results. In the inventors' experience, the set up in Equation (10) turned out to be an excellent compromise from the viewpoint of accuracy of results and efficiency. For comparison, recognition results for the thorax are demonstrated in Table 9 with no orientation adjustment for any object in both model building and recognition.

TABLE-US-00009 TABLE 6 Recognition results (mean, standard deviation) for Thorax for the strategy in (10). (“Mean” excludes VS.) TSkn RS TSk IMS LPS TB RPS E PC AS VS Mean Location 3.9 5.5 9.0 5.6 6.3 11.6 10.4 9.8 8.6 10.7 31.8 8.1 Error (mm) 1.5 2.3 5.0 3.5 3.1 5.0 4.7 4.8 5.0 5.4 12.0 4.0 Size error 1.0 0.99 0.96 0.95 0.97 0.91 0.98 0.9 0.95 1.01 0.77 0.96 0.01 0.02 0.05 0.05 0.03 0.06 0.04 0.14 0.05 0.08 0.06 0.06

TABLE-US-00010 TABLE 7 Recognition results (mean, standard deviation) for Abdomen for the strategy in (10). ASkn SAT ASk Lvr ASTs Kd Spl Msl AIA IVC RKd LKd Mean Location 5.9 20.2 11.7 7.2 7.2 10.6 11.6 7.7 8.2 8.7 11.3 7.3 9.8 Error (mm) 3.4 8.5 7.9 5.4 3.0 9.8 13.9 3.6 2.8 7.2 11.6 7.4 7 Size error 1.0 0.97 0.96 0.93 1.0 0.94 1.2 1.01 1.1 1.15 0.97 0.93 1.01 0.02 0.03 0.06 0.07 0.02 0.09 0.19 0.03 0.13 0.1 0.1 0.08 0.07

TABLE-US-00011 TABLE 8 Recognition results (mean, standard deviation) for Neck for the strategy in (10). NSkn A&B FP NSTs Mnd Phrx Tnsl Tng SP Ad NP OP RT LT Mean Location 3 7.8 4.2 4.8 12.5 10.4 2.8 4.9 5.1 1.8 11.1 10 2.9 2.3 5.96 Error (mm) 1.2 3.8 2.1 2.1  3.7  4.5 1.8 2.8 1.8 0.8  6.8  8.7 2.2 2.1 1.96 Size error 1 0.9 1 0.92  0.74  0.8 1 1.02 0.93 0.9  0.65  0.74 0.92 0.9 0.93 0.01 0.03 0.03 0.06  0.05  0.04 0.1 0.06 0.24 0.12  0.07  0.2 0.11 0.12 0.04

TABLE-US-00012 TABLE 9 Recognition results for Thorax with no orientation alignment. (“Mean” excludes VS.) TSkn RS TSk IMS LPS TB RPS E PC AS VS Mean Location 3.9 5.5 9 5.6 6.3 8 10.4 14.2 8.6 8.1 33.6 8.0 Error (mm) 1.5 2.3 5 3.5 3.1 6.5 4.7 10.5 5 7.5 15.1 4.9 Size error 1.01 0.99 0.96 0.95 0.97 0.83 0.98 0.85 0.95 0.99 0.77 0.95 0.01 0.02 0.05 0.05 0.03 0.08 0.04 0.12 0.05 0.08 0.06 0.05

[0115] Size error is always close to 1 for all body regions and objects. Generally, recognition results for non-sparse objects are excellent with a positional error of mostly 1-2 voxels. Note that for DS1 and DS2, voxels are quite large. In particular, since recognition results do not improve much with finer discretization of the model but only increase computation for recognition, the inventors constructed models with isotropic voxels of side equal to one half of the largest dimension of the voxels in the original data. Thus, for DS1 and DS2, the model voxels are of size 2.5×2.5×2.5 mm.sup.3. The inventors observed that the positional accuracy within the slice plane is better than across slices. In other words, errors listed in the tables are mostly in the third dimension in which voxel size is large. Orientation adjustment improves recognition somewhat for some sparse objects, but has negligible effect for non-sparse objects, at least in the thorax.

[0116] The recognition results for the MRI data set DS4 are demonstrated in FIG. 9 and Table 10. Again, since the model is fuzzy, it will encroach into adjacent tissue regions with some membership value. Since the goal here was just to measure subcutaneous adiposity, the hierarchy was simplified as shown in FIG. 9. Again the position error is 1-2 voxels. These results are particularly noteworthy since they are generated by using the models built from image data sets acquired from a different modality, namely CT, and for a different group with an age difference of about 40 years and with a different gender. This underscores the importance of understanding the dichotomy between recognition and delineation. Recognition is a high-level and rough process which gives anatomic context. The models do not have to be detailed attempting to capture fine details. Obtaining the anatomic context is a necessary step for achieving accurate delineation. It is important to note here that for the cross modality operation to work in this manner, the MR image intensities should be standardized (Nyul and Udupa 1999).

TABLE-US-00013 TABLE 10 Recognition accuracy for the objects shown in FIG. 9 ASkn SAT Position 4.6 12.97 Error 2.5 5.3 (mm) Size Error 1.01 1 0.05 0.03

Object Delineation

[0117] Sample delineation results are displayed in FIGS. 10-13 for DS1-DS4. Delineation accuracy statistics for these data sets, expressed as false positive and false negative volume fractions (FPVF, FNVF) as well as mean Hausdorff distance (HD) between the true and delineated boundary surfaces, are listed in Tables 11-14. The HD measure is defined as the mean over all test subjects of the median of the distances of the points on the delineated object boundary surface from the true object boundary surface.

[0118] Delineation results for VS (Thorax) are not presented since the recognition accuracy for VS is not adequate for reliable delineation. It is noted that the delineation of 21 non-sparse objects achieves a mean FPVF and FNVF of 0.02 and 0.08, respectively, and a mean HD of 0.9 voxels, which are generally considered to be excellent. Six sparse objects also achieve good delineation outcome, with the above mean measures reading 0.05, 0.15, and 1.5, respectively. However, sparse objects VS, E, IVC, Mnd, and NP pose challenges for effective delineation. Often, even when their recognition is effective, it is difficult to guarantee placement of seed sets A.sub.O and A.sub.B appropriately within and outside these objects because of their sparse nature. In DS3 (MR images of neck), it is very difficult to properly delineate Mnd, NP, and OP because of their poor definition in the image. To test the effectiveness of the models created from these data (DS3) in segmenting the same objects on CT data of a group of three different pediatric subjects, the inventors devised a simple hierarchy with NSkn as the root and with Mnd, NP, and OP as its offspring objects. The delineation results obtained for these four objects were excellent, with a mean FPVF of 0, 0.01, 0, and 0.02, and mean FNVF of 0.01, 0.01, 0.02 and 0.1, respectively.

TABLE-US-00014 TABLE 11 Delineation results for Thorax (mean and standard deviation). TSkn RS TSk IMS LPS RPS E PC TB AS FPVF 0.02 0.0 0.19 0.03 0.01 0.01 0.0 0.01 0.01 0.01 0.02 0.0 0.05 0.01 0.03 0.02 0.0 0.00 0.00 0.00 FNVF 0.05 0.06 0.13 0.07 0.04 0.04 0.49 0.09 0.16 0.17 0.06 0.04 0.07 0.07 0.02 0.02 0.19 0.06 0.14 0.17 HD 3.6 1.24 10.6 6.2 2.9 2.1 3.1 3.5 5.2 5.3 (mm) 4.5 0.42 2.4 1.8 8.8 4.7 0.87 1.3 1.8 2.5

TABLE-US-00015 TABLE 12 Delineation results for Abdomen (mean & standard deviation). ASkn ASk Lvr ASTs SAT RKd LKd Spl Msl AIA FPVF 0.01 0.06 0.04 0.12 0.05 0.00 0.01 0.0 0.13 0.01 0.00 0.01 0.02 0.05 0.03 0.00 0.01 0.0 0.03 0.0 FNFV 0.05 0.14 0.1 0.15 0.12 0.13 0.1 0.13 0.09 0.13 0.08 0.09 0.05 0.09 0.02 0.04 0.02 0.03 0.08 0.03 HD 1.7 6.9 5.3 1.74 1.6 2.4 5.4 6.8 2.5 5.6 (mm) 2.7 1.5 1.6 1.0 0.8 1.1 4.8 6.0 1.1 1.8

TABLE-US-00016 TABLE 13 Delineation results for Neck (mean & standard deviation). NSkn FP Mnd NP OP RT LT Tng SP Ad FPVF 0.0 0.0 0.01 0.01 0.0 0.01 0.01 0.02 0.01 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.01 0.01 0.0 FNFV 0.0 0.1 0.49 0.32 0.2 0.06 0.06 0.02 0.08 0.07 0.01 0.05 0.08 0.2 0.02 0.02 0.01 0.01 0.01 0.04 HD 2.8 0.83 3.3 3.8 7.6 3.3 3.2 8.4 8.03 2.2 (mm) 0.06 0.53 0.56 1.01 2.4 0.62 1.4 1.92 4.0 0.3

TABLE-US-00017 TABLE 14 Delineation results for DS4. ASkn SAT FPVF 0.0 0.06 FNVF 0.03 0.01 HD (nun) 1.7 3.9

Comparison With a Non-Hierarchical Approach

[0119] To study the effect of the hierarchy and the knowledge encoded in it on recognition, Table 15 lists the recognition performance of a non-hierarchical approach. The results are shown for Thorax wherein each object is recognized on its own by using the same fuzzy models FM() as used in the hierarchical AAR system. The initial pose for search is taken to be the center of the image and search range covers roughly the whole body region with the scale factor range the same as that for the hierarchical approach. In comparison to the hierarchical approach (Tables 6 and 9), it is clear that non-hierarchical recognition performance is much worse.

TABLE-US-00018 TABLE 15 Recognition results for Thorax: non-hierarchical approach (mean & standard deviation). TSkn RS TSk IMS LPS TB RPS E PC AS VS Mean Location 10.5 12.9 21.1 27.7 91.4 53.3 72.3 42.4 45.5 23.1 82.2 43.8 Error (mm) 9.5 13.1 21.8 9.8 10.8 20.9 12.9 34.5 12.5 15.2 33.8 17.7 Size error 1.0 1.01 0.96 0.92 0.8 0.82 0.8 0.86 0.9 0.97 0.81 0.9 0.02 0.09 0.08 0.07 0.09 0.06 0.07 0.14 0.06 0.11 0.08 0.08

Computational Considerations

[0120] Program execution times are estimated on a Dell computer with the following specifications: 4-core Intel Xeon 3.6 GHz CPU with 8 GB RAM and running the Linux-jb18 3.7.10-1.16 operating system. Mean computational times for the AAR steps are listed in Table 16. Model building includes the construction of fuzzy models and the estimation of ρ, λ, and all parameters related to recognition and delineation, including the optimal threshold parameters . This latter step takes about 12 seconds per object. As seen from Table 16, each of the three main operations takes under 1 minute per object. Among these operations, only the time for model building depends on the number of training data sets, while recognition and delineation are independent of this factor. On average, model building times per object per training data set for Thorax, Abdomen, and Neck are, respectively, 1.4 sec, 1.7 sec, and 1 sec. In statistical atlas based methods, the computational time for image registration becomes the bottleneck. The calculation taking Elastix as a representative registration tool kit (Klein, S., Staring, M., Murphy, K., Viergever, M. A., Pluim, J. P. W., 2010. Elastix: a toolbox for intensity based medical image registration, IEEE Transactions on Medical Imaging, 29, 196-205) indicates that the creation of a single atlas for each of the 11 objects of the Thorax at a reduced image resolution of 2.5×2.5×2.5 mm.sup.3 for the 25 training data sets of DS1 would take about 23.5 hours compared to 6.4 min for the AAR system. The time per object for recognition and delineation can also take several minutes for these methods. Even with 100 data sets for training and 15 objects in a body region, the total time needed for the AAR model building step would be about 40 minutes, whereas atlas building may take days to complete especially when multi-atlas strategies are used.

TABLE-US-00019 TABLE 16 Mean computational time in seconds per object for different operations and body regions. Operation Thorax Abdomen Neck Model building 35 42 24 Object recognition 30 46 6 Object delineation 47 56 24

Comparison With Other Methods

[0121] The publications reporting works that are directly related to the invention in spirit are (Baiker et al. 2010, Chu et al. 2013, Criminisi et al. 2013, Lu et al. 2012, Linguraru et al. 2012, Okada et al. 2008, Zhou et al. 2012). In Table 17, a comparison to the AAR system of the invention is provided based on the results reported in these works. It is noted that a quantitative grading/understanding of the methods is impossible since the data sets used, acquisition protocols and resolutions, considered objects, training and test data set subdivisions, cross validation strategies, and computing platforms are all different in these methods. Interestingly, a commonality among them is that they all focused on CT image data sets.

TABLE-US-00020 TABLE 17 A comparison with the current methods from the literature that are related to the inventive methods. Unknown and irrelevant entries are indicated by “~”. Region Training- overlap to-test (Dice, Voxel data Location Jackard size pro- error Index Method Objects (mm.sup.3) portion (mm) (JI), etc.) Lu et al. Prostate, ~ × ~ × 141 to 2.4 to 4.2 ~ 2012 bladder, 0.8 to 5 47, 4-fold rectum Linguraru Liver, (0.5 to 27 to 1, 0.8 to 1.2 90.9% to et al. spleen, 0.9).sup.2 × 1 28-fold 94.8% 2012 kidneys to 5 Okada et Liver, 0.7 × 20 to 8 (for liver) 88% al. 2008 vena cava, 0.7 × 2.5 1.5 to 2.8 Gallbladder Chu et al. Liver, (0.55 to 90 to 10, ~ 56% 2013 spleen, .82).sup.2 × 10-fold (pan- pancreas, 0.7 to 1 creas- kidneys (estimated) JI) to 95.2% (liver- Dice) Criminisi 26 anatomic (0.5 to 308 to 82 9.7 to 19.1 ~ et al. structures 1).sup.2 × 1 (mean for 2013 in the torso to 5 each structure) Zhou et 12 organ (0.6 to 300 to 6 to 14 ~ al. 2012 regions in 0.7).sup.3 1000 for mode thorax, locations abdomen, pelvis Baiker et Brain, (0.332).sup.3 MOBY ~ 47% to al. 2010 heart, kidneys, atlas, 73% lungs, liver, 26 skeleton datasets

[0122] Among these methods, (Chu et al. 2013, Linguraru et al. 2012, Lu et al. 2012, Okada et al. 2008) comprise one group wherein the body region of focus was the pelvis or abdomen, with 3-5 objects considered for segmentation. They all employ an object localization step, which is achieved either through an atlas (Chu et al. 2013, Linguraru et al. 2012, Okada et al. 2008), statistical shape models (Okada et al. 2008), or machine learning techniques (Lu et al. 2012), and subsequently a delineation step that uses graph cuts (Chu et al. 2013, Linguraru et al. 2012), information theory (Lu et al. 2012), and MAP or ML estimation (Chu et al. 2013, Okada et al. 2008). In the second group (Criminisi et al. 2013, Zhou et al. 2012), the aim is only to locate the objects via machine learning techniques. The third group is constituted by (Baiker et al. 2010), the only work that considered body-wide organs, but in mice, using a kinematic model of the skeletal joints to localize objects relative to different skeletal components.

[0123] It is observed that, for the same objects (liver, kidneys, and spleen), the results expressed herein are comparable to, often better than, the current results from literature, especially considering the 5 mm slice spacing and the equal training-to-test data set proportion for the evaluation. Those skilled in the art will appreciate that the development of a general AAR system that can be readily applied and adapted to different body regions, multitudes of organs, and modalities has not yet been demonstrated in the literature.

Observations

[0124] The methods presented herein provide a general body of techniques for automatic anatomy recognition and delineation whose principles are not tied to any specific body region, organ system, or imaging modality. The inventors took a fuzzy approach for building the models and attempted to harness as much specific anatomic information as possible to be embedded into the fuzzy anatomic model. The inventors demonstrated the generality of the approach by examining the performance of the same AAR system on three different body regions using CT and MR image data sets. The inventors also illustrated the potential of the system for rapid prototyping by demonstrating its adaptability to a new application on a different modality (DS4). The system is set up to operate fully automatically. All image modality-specific parameters needed—threshold intervals for objects in for recognition and affinity parameters for delineation—are estimated automatically from the training data sets. When a new application is sought at a modality different from those considered in the anatomy model FAM(, G), a few sample segmentations of the objects of interest and the matching images are needed for relearning these image intensity-related parameter values (specifically, and the affinity parameters). All other modality-independent aspects of the model do not need retraining. In the case of MRI, images from each separate MRI protocol have to be standardized for image intensity so that setting up these parametric values becomes sensible. Separation of modality-independent from dependent aspects, organization of objects in a hierarchy, encoding object relationship information into the hierarchy, optimal threshold-based recognition learning, and fuzzy model-based IRFC are novel and powerful concepts with consequences in recognition and delineation, as demonstrated herein.

[0125] While the above strengths of the AAR system of the invention are quite unique, the system has some limitations at present. First, the inventors have not studied the performance of the system on patient images that contain significant pathology. However, it is noted that DS4 indeed includes image data sets of patients who are obese. Note also that these image data sets are from a very different age and gender group and on a different imaging modality from those used to build FAM(, G). The inventors believe that it is important to make the system operate satisfactorily on normal or near-normal images before testing it on images with diverse pathologies. Second, the accuracy is inadequate for some sparse objects for recognition (VS, IVC) and delineation (E, Mnd, NP). Also, the inventors have not considered herein other important and challenging sparse objects such as the adrenal glands, pancreas, and the spinal cord. If recognition is inadequate, delineation will become unacceptable because it becomes impossible to appropriately initialize the delineation process and to exploit the model for making up for missing boundary information in the image in delineation. When these cases were closely examined, it became clear that there are fundamental challenges in the model building stage itself for sparse objects. Generally, the inventors found that sparse objects have much greater variation than their non-sparse counterparts in form, topology, and geographic layout, compared to their size. As an example, consider AS and VS (Thorax). The descending aortic portion of AS is often straight and directed vertically downward while in some subjects it may be inclined, curved, or even tortuous, with other portions, especially the aortic arch, not varying much. The branching pattern of the left and right brachiocephalic veins and the course of the azygos vein in VS also vary considerably. In view of such difficulties, the inventors have come to the realization that sparse objects should not be modeled directly from their precise shape information in the binary image set .sup.b, instead only their rough super form (such as a minimal super set that subsumes such variations) should be utilized in model building. The inventors are exploring the use of rough sets (Maji, P., Pal, S. K., 2012. Rough-Fuzzy Pattern Recognition: Applications in Bioinformatics and Medical Imaging, John Wiley & Sons, Inc. New York) for this purpose.

[0126] The AAR methodology seems to have definite computational advantages over atlas-based approaches. Further, in atlas-based methods, it is perhaps much more challenging to incorporate the extensive object-level knowledge that the AAR approach exploits at various stages for recognition and delineation. These incorporations constitute highly non-linear and discontinuous phenomena which are effected in intensity, geometric, and topological spaces. The kinematic model employed by Baiker et al. 2010 is a good analogy of how one may encode object relationships via a model that are difficult to emulate through continuous and smooth image/atlas deformations.

[0127] The problem of automatically determining the body region following the definition of within the given data set was not explored herein. As demonstrated in (Chen et al. 2012), it is possible to determine the slices delimiting a body region B automatically based on slice profiles. Furthermore, the information about the relationship between and WB can also be encoded into the hierarchy as illustrated in FIG. 2(a) for each .

[0128] The use of composite objects often leads to better recognition accuracy. This is because the multiple objects contained in a composite object offer tighter constraints in recognition search. The aspect of how objects can be grouped to achieve optimum recognition results needs further investigation. A related topic is how to device optimal hierarchies for a given body region. The hierarchies the inventors have considered so far are anatomically motivated. Perhaps there are “optimal” hierarchies from the view point of achieving the best recognition (and hence, delineation) results. In such an investigation, matters of how objects should be grouped as well as ordered in the hierarchy can both be addressed simultaneously using graph optimization techniques.

[0129] The inventors have set up the AAR-R and AAR-D procedures in a general way. Recognition and delineation algorithms other than those described herein can be used independently for R-ROOT and R-OBJECT and for D-ROOT and D-OBJECT within the same hierarchical set up. Similar to composite object recognition, delineation done simultaneously for multiple objects, unlike the one-object-at-a-time approach of AAR-D, may improve overall accuracy.

[0130] Computationally, there are three expensive operations in the AAR system—image interpolation, distance transform, and the delineation algorithm (FMIRFC). To make recognition and delineation operate in practical time in a clinical setting, implementations of these operations will desirably be sped up. Toward this goal, these operations may be implemented in a GPU. GPU implementations of some fuzzy connectedness algorithms have already been published (Zhuge, Y., Cao, Y., Udupa, J. K., Miller, R. W., 2011. Parallel fuzzy connected image segmentation on GPU. Medical Physics 38 (7), 4365-4371; and Zhuge, Y., Ciesielski, K. C., Udupa, J. K., Miller, R. W., 2013. GPU-based relative fuzzy connectedness image segmentation. Medical Physics 40 (1), 011903-011903-10).

[0131] Finally, along the lines of the study underlying DS4, those skilled in the art will appreciate that AAR system may be adapted to several clinical applications. Some of the avenues the inventors are currently exploring for the proposed AAR approach are delineated in the examples provided in more detail below. For example, the AAR techniques described above may be modified for quantifying abdominal fat through the use of standardized anatomic space, providing automatic localization of IASLC-defined mediastinal lymph node stations, and providing radiation therapy planning in exemplary embodiments as described in detail below.

Example 1: Optimization of Abdominal Fat Quantification on CT Imaging through use of Standardized Anatomic Space

Notations and Overall Approach

[0132] As noted above, it is desired to answer three questions related to fat quantification that have not been addressed in the literature. How does one ensure that the slices used for correlation calculation from different subjects are at the same anatomic location? At what single slice anatomic location do the areas of SAT and VAT estimated from a single slice correlate maximally with the corresponding volume measures? Are there combinations of multiple slices (not necessarily contiguous) whose area sum correlates better with volume than does single slice area with volume? The techniques of the invention adapted the AAR methodology described above to address these questions.

[0133] Let V(B, Q, G) denote the set of all possible 3D images of a precisely-defined body region , taken as per a specified image acquisition protocol Q, from a well-defined group of subjects G. For example, B may be the abdominal region, which is defined by its superior bounding plane located at the superior most aspect of the liver and its inferior bounding plane located at the junction where the abdominal aorta bifurcates into common iliac arteries. Variable Q may be CT imaging with a specified set of acquisition parameters, and G may denote normal male subjects in the age range of 50 to 60 years. The reason for relating all analysis to a specified set V(, Q, G) is that, it may not be possible to generalize the conclusions drawn about fat distribution when one changes some of the variables associated with V, especially patient group G. The inventors denote by V the set of images available for the study, which is assumed to be a representative subset of V(, Q, G). Let P be an image in V of some subjects of his body region . The inventors view I.sup.s as a set of ns axial slices

I.sup.s={S.sub.1.sup.s, . . . , S.sub.n.sub.s.sup.s}.

Since I.sup.s is an image of , S.sub.1.sup.s and S.sub.n.sub.s.sup.s represent anatomic planes bounding . It is assumed that they correspond to the superior and inferior bounding planes, P.sub.H.sup.s, and P.sub.D.sup.s of of subject s, respectively. All locations and coordinates are assumed to be specified with respect to a fixed Scanner Coordinate System (SCS) for all subjects. If the acquired images have extra slices, it is assumed that they have been removed to satisfy this condition. Note that if I.sup.s and I.sup.t are images in V(, Q, G) of two subjects s and t, then the number of slices n.sub.s and n.sub.t representing in the two subjects may not be equal. Similarly, the same numbered slices in I.sup.s and I.sup.t may not correspond to the same anatomic axial location in subjects s and t. Suppose one discretizes the anatomic axial positions in from P.sub.H.sup.s and P.sub.D.sup.s into L anatomic locations l.sub.1.sup.s, . . . , l.sub.L.sup.s such that l.sub.1.sup.s and l.sub.L.sup.s always correspond to P.sub.H.sup.s and P.sub.D.sup.s, respectively, for all subjects s. For example, position l.sub.i.sup.s may correspond to the location of an axial plane passing through the middle of the body of the L1 lumbar vertebra of subjects, in this case, l.sub.i.sup.t represents an axial plane at the same anatomic location for subject t. Locations l.sub.1.sup.s, . . . , l.sub.L.sup.s may be also thought of as representing anatomic landmarks labeled l.sub.1, . . . , l.sub.L. In the above example, l.sub.i, is the name of the landmark associated with location l.sub.i.sup.s. The inventors denote these anatomic landmarks by the ordered set AL={l.sub.1, . . . , l.sub.L}. In order to perform volume to area correlation analysis correctly, the inventors need to first assign a correct label from the set AL to every slice in every image I.sup.s in V. Since it is customary to use the vertebral column as reference for specifying homologous anatomic locations, the inventors will follow this same approach herein. It is noted however that the methods are general and can use any other reference system for locations. The inventors think of anatomic landmarks to be defined in a Standard Anatomic Space (SAS), and the process of assigning labels from AL to slices in any given image I.sup.s as a mapping from SCS to SAS.

[0134] The overall approach to seek answers to the three questions posed above is depicted schematically in FIG. 14. The four steps involved are described below in detail.

Segmenting SAT and VAT Regions in Images in V

[0135] In this first step (FIG. 14, step 10), SAT and VAT regions in the images in V are automatically segmented by modifying the automatic anatomy recognition (AAR) system. As described above, the AAR system operates by creating a fuzzy anatomy model for a body region 13 and subsequently using this model to recognize and delineate organs in . The fuzzy anatomy model consists of a hierarchical arrangement of the organs of , a fuzzy model for each included organ, and organ relationships in the hierarchical order. The modification for purposes of this application includes considering just three objects—skin boundary, SAT, and VAT, with skin as the root object and SAT and VAT as its offspring objects, in place of all 10-15 major organs of that are otherwise included in the model. The rest of the processes remain the same as the recognition and delineation methods described above. If V is a set of MR images, then image background non-uniformity correction and intensity standardization will have to be performed before applying AAR-based segmentation.

Assigning Landmark Labels l.SUB.1., . . . , l.SUB.L .to Slices in each Image P

[0136] Ideally, once the set AL of anatomic landmarks is selected, one could identify manually the anatomic location, and hence the landmark label, to be assigned to each slice S.sub.i.sup.s of each image I.sup.s of V. Such a manual approach can be realized on CT imagery as follows. First segment the vertebral column in I.sup.s, create 3D surface renditions of the column, and interactively indicate in this display the axial locations. The inventors use a visualization software package (e.g., CAVASS) for selecting locations quantitatively precisely on shaded surface renditions. In MR images, however, this approach will be more difficult since segmentation of the vertebral column is challenging. Since this manual approach is labor intensive, the inventors explored two alternative approaches—linear and non-linear, and compared them to the manual approach. In all approaches, the input is the set V of images and the result is a mapping that indicates the anatomic location (label) associated with each slice of each image of V (FIG. 14, step 20).

Linear Approach

[0137] This approach assumes that, once it is guaranteed that the bounding planes P.sub.H.sup.s and P.sub.D.sup.s, and hence slices S.sub.1.sup.s and S.sub.n.sub.s.sup.s of image I.sup.s, correspond to locations l.sub.1.sup.s and l.sub.L.sup.s, respectively, then anatomic locations corresponding to slices S.sub.2.sup.s and S.sub.n.sub.s.sub.-1.sup.s can be found by linearly mapping the n.sub.s slices from S.sub.1.sup.s and S.sub.n.sub.s.sup.s to L slices U.sub.1.sup.s and PU.sub.L.sup.s via linear interpolation for any subjects. Note that L can be less than, or greater than, or equal to n.sub.s. The only requirement on L is that it should be at least 2. Of course, it is possible that in this approach the mapped location of a slice U.sub.i.sup.s may not match the true location of landmarks. To implement this approach, the inventors first identify the data set in V whose domain is the smallest in the longitudinal direction in terms of the number of slices, take this number to be L, and then linearly interpolate all other data sets to yield this number of slices. The inventors assume the slice locations of this data set to correspond approximately to landmarks l.sub.1, . . . , l.sub.L. Mapping of anatomic locations to all other data sets is now established by the correspondence of the same numbered slices. In this manner, for any given slice number for any subject, the corresponding (linearly mapped) slice numbers for all other subjects are identified. A drawback of the linear approach is that non-linearities in the relationships among anatomic locations used as landmarks in the longitudinal direction cannot be accounted for.

Non-Linear Approach

[0138] In the linear approach, two anatomic landmarks l.sub.1 and l.sub.L were employed to anchor the first and the last slice of B and to predict the anatomic location of all other slices. Generally such a linear mapping does not yield locations that are sufficiently close to actual anatomic locations (as demonstrated below) of landmarks. This deficiency can be overcome by non-linear mapping. In this approach, in addition to I.sub.1 and l.sub.L, other key anatomic landmarks are used to refine mapping. The method consists of two stages—calibration and transformation.

[0139] The purpose of the calibration stage is to learn any non-linearities that may exist in the relationships among anatomic locations. (Here “learning” does not have the same meaning as “training” widely used in machine learning.) Typically, the inventors select M<L key anatomic landmarks, denoted by m.sub.1, . . . , m.sub.M, from among l.sub.1, . . . , l.sub.L. In the present embodiment, the inventors selected the mid-points (in the vertical direction) of the vertebral bodies from T11 to L4 as key landmarks (so M=6). Next, these key landmarks are identified manually in a set T⊂V(, Q, G) of images. For any image I.sup.s in T, the locations of these key landmarks for subjects will be denoted by m.sub.1.sup.s, . . . , m.sub.M.sup.s. A standard anatomic scale is then determined to be of length which is the largest of the lengths from P.sub.H.sup.s and P.sub.D.sup.s over all data sets in T. Locations m.sub.1.sup.s, . . . , m.sub.M.sup.s for every data set in T are then mapped linearly on to the standard scale (see FIG. 15) and the mean positions μ.sub.1, . . . . , μ.sub.M of the key points on the standard scale over all mapped data sets of T are computed. The mapping from SCS to SAS is subsequently determined to be the piece-wise linear function that maps m.sub.1.sup.s, . . . , m.sub.M.sup.s to μ.sub.1, . . . , μ.sub.M as depicted in FIG. 16.

[0140] In the transformation stage (FIG. 16), given any image I.sup.s, first the locations of its anatomic landmarks m.sub.1.sup.s, . . . , m.sub.M.sup.s are identified. Then the mapping function from SCS to SAS determined in the calibration stage is used to determine the label to be assigned to each slice S.sub.i.sup.s of l.sup.s. The algorithm, called SAS, for mapping to the standardized anatomic space, summarized herein, is straightforward and requires no special data structures or optimization in implementation.

TABLE-US-00021 Algorithm SAS Input: Two disjoint sets of images T and V, T ⊂ (   , Q, G), V ⊂ (   , Q, G); AL; { m.sub.1, ..., m.sub.M}. Output A mapping form SCS to SAS; the set V of images with a label assigned to each slice of each image of V. Begin Calibration Stage C1. Determine standard scale and identify key landmarks m.sub.1, ..., m.sub.M in each image in T; C2. Map key landmarks linearly to standard scale; C3. Estimate mean locations μ.sub.1,... μ.sub.M of key landmarks to standard scale; Transformation Stage T1. For each image I.sup.s of V and for each of its slices, determine its key landmark locations m.sub.1.sup.s, ..., m.sub.M.sup.s; T2. Find the mapping of these locations as per SCS to SAS function; T3. Based on this mapped value assign label to each slice of l.sup.s;

Experimental Results and Discussion

Image Data

[0141] Variables G and Q defining V(, Q, G) for the experiments were as follows. Contrast-enhanced abdominal CT image data sets from fifty 50-60 year-old male subjects with an image voxel size of 0.9×0.9×5 mm.sup.3 were utilized in the study. The subjects were radiologically normal with exception of minimal incidental focal abnormalities. The abdominal body region was defined in the same way for the 50 subjects, with P.sub.H.sup.s located at the superior most aspect of the liver and P.sub.D.sup.s corresponding to the point of bifurcation of the abdominal aorta into common iliac arteries. Of the 50 data sets, 5 were used for calibration (constituting T) and the rest (constituting V) were used for testing.

[0142] To illustrate the anatomic variability that exists among subjects, FIG. 17 plots schematically the locations of the mid points of vertebral bodies in the cranio-caudal (vertical) direction for all subjects considered in the study. The top and the bottom of the vertical line drawn for each subject indicate the extent of in relation to the vertebral bodies. For example, in subject numbered 50 (the right-most location on the abscissa), the abdominal region starts from roughly the T11 vertebra and ends at the L5 vertebra. The locations of both the top-most and bottom-most slices have significant variability in terms of anatomic correspondence as seen in FIG. 17. To further illustrate qualitatively the variability in the layout of the vertebrae among subjects, FIG. 18 displays surface renditions of the skeletal components in for some of the 50 subjects who show wide variation in FIG. 17.

Correlation Analysis

[0143] To study the nature of the volume-to-area correlation, the inventors analyzed the relationship between 3D volume and area estimated from a single slice as well as summed up areas estimated from 2 and 3 slices where the slices were selected at all possible locations and not necessarily contiguously situated (FIG. 14, step 30).

Correlation With Single Slice

[0144] The inventors considered 34 subjects for correlation analysis by selecting those subjects whose body region covered vertebrae from T10 to L4 as a common/overlap region among the 50 subjects. The reason for this decision is to guarantee that the body region of the subjects for calculating correlation will be in the same anatomic range in SAS. Some subjects for whom slices start from T12 or even higher positions as shown in FIG. 18 are not selected. For all 34 subjects, six spinal landmarks were selected from T10 to L3 as the mid points of the respective vertebral bodies. Although the method is illustrated by using 6 landmarks here, this number can be set to any value greater than or equal to 2 and any other landmarks can also be used.

[0145] In order to study how correlation may vary for different anatomic slice locations, FIG. 19 (using data from 34 subjects) displays the correlation values as a curve for different slice locations for SAT and VAT by using both linear and non-linear mappings. Some key landmark positions are indicated along the horizontal axis in the bottom row of FIG. 19. The number of slices for linear and non-linear mapping is different and as such there is not much meaning in comparing the slices for the two methods by numbers. This is due to the fact that for linear mapping, the mapped slices are found by mapping the total number of slices to the same number of (smallest/largest) slices for every subject. For non-linear mapping on the other hand, because of the fact that the distance between successive landmarks is allowed to be different for different subjects, the total number of slices in SAS may not be the same as that of linear mapping. To examine how the location of maximum correlation may vary across subjects (FIG. 14, step 40), FIG. 20 displays the anatomic landmark locations at which maximum correlation occurred for SAT and VAT for the two methods for different subjects. FIG. 21 demonstrates the anatomic locations where maximum correlation is achieved by the two methods for two sample subjects, where slice images and their locations in the anatomic space for both SAT and VAT are shown with reference to the 3D rendered spine.

[0146] The following observations may be made from FIGS. 19-21. The maximum correlation for VAT derived from non-linear mapping achieved a 10% higher value than that of linear mapping. It is also clear that the maximum correlation occurs at different anatomic locations for SAT and VAT. With non-linear mapping, the maximum correlation occurs at T12 for SAT and L3 for VAT as shown in the bottom figures in FIG. 19. The maximum correlation also occurs at different anatomic locations for the different mapping methods. For SAT and VAT, the locations of maximum correlation derived from linear mapping vary considerably among subjects. The goal of correlation calculation is to find one or more anatomic locations which are optimal for estimating abdominal fat distribution. The slices from all subjects used for correlation calculation are expected to be at the same anatomic location. Correlation testing will have no meaning if every subject has the slice considered for correlation calculation at a different location. The site of maximum correlation derived from non-linear mapping has more precision than from simple linear mapping. In the bottom row of FIG. 20, there is a small variation of locations over all subjects which is less than 5.0 mm (same as slice spacing), implying that the slice localizations are anatomically very precise in the SAS (see Table 18). Table 19 shows the correlation derived from a single slice at different true anatomic locations where the maximum correlation for VAT occurs at L3-L4 and at T12-L1 for SAT. The results are similar to the correlation calculated from the non-linear method as shown in the bottom row of FIG. 20. Again, the anatomic location of maximum correlation for SAT is different from that of VAT.

TABLE-US-00022 TABLE 18 Correlation coefficients and slice location variation (in mm) for linear and non-linear mapping techniques. Correlations shown are maximum values. Single slice Linear Non-linear at L4-L5 mapping mapping SAT VAT SAT VAT SAT VAT Correlation 0.74 0.87 0.89 0.81 0.88 0.92 Location — — 17.80 15.70 4.38 2.63 variation (mm)

TABLE-US-00023 TABLE 19 Correlation with single slice at different true anatomic locations Anatomic Correlation with single slice slice location SAT VAT T10-T11 0.85 0.79 T11-T12 0.87 0.81 T12-L1 0.88 0.90 L1-L2 0.88 0.89 L2-L3 0.85 0.82 L3-L4 0.76 0.92 L4-L5 0.74 0.87

[0147] Examining the top two rows derived from linear mapping for SAT and VAT in FIG. 21, it is observed that the anatomic locations of maximum correlation are significantly different for different subjects. Yet, for the non-linear mapping, for both SAT and VAT, the anatomic locations of maximum correlation are much closer even though they come from different subjects. In particular, for the non-linear mapping, there appears to be relatively constant high correlation values for SAT in the lower thoracic/upper abdominal region and for VAT in the lower abdomen.

[0148] To test the sensitivity of the results to the choice of the calibration data set, FIG. 22 shows the single-slice correlation curves for SAT and VAT derived from a different set of randomly chosen 5 calibration data sets. The maximum correlation achieved was again 0.88 and 0.91 for SAT and VAT, respectively, which are the same as the results in FIG. 19, and the slice locations where fat volume maximally correlated with fat area are also the same. The curves are remarkably similar, as expected, except for some minor differences at the ends of the curves. Compared with linear mapping and earlier methods, one advantage of the proposed non-linear mapping approach is to guarantee that the slice where maximum correlation occurs is at the same anatomical location irrespective of patient-to-patient anatomical variability.

Correlation With Multiple Slices

[0149] To address the question as to whether single slice or multiple (contiguous or non-contiguous) slices yield better area-to-volume correlation, the inventors calculated the correlation by using multiple slices with both linear and non-linear mappings.

TABLE-US-00024 TABLE 20 Correlation by using one or more slices per subject, where correlations shown are maximum values for the two mapping techniques. Multiple slices 1 slice 2 slices 3 slices Linear mapping SAT 0.89 0.90 0.91 VAT 0.81 0.84 0.85 Non-linear mapping SAT 0.88 0.88 0.89 VAT 0.92 0.95 0.95 Slices at L4-L5 SAT 0.74 0.74 0.75 VAT 0.87 0.88 0.89

[0150] Table 20 lists the maximum correlation achieved by using one, two, and three slices per subject with linear and non-linear mapping. The correlation derived from the slice at the L4-L5 junction is also listed for comparison since this location is most commonly used. For this case, the choice of 2 and 3 slices is such that the slices are contiguous and they are as close to the L4-L5 junction as possible. Note that non-linear mapping with multiple slices achieved the highest correlation. Table 21 lists anatomic locations where maximum correlation is achieved for the two methods. Slice locations are shown after mapping (linear and non-linear) where the maximum correlation is achieved. The values listed in the table are the slice numbers in the volume files (where number 1 indicates the bottom slice of the abdominal region and larger numbers are located closer to the top of the abdominal region). The anatomic locations in SAS are also listed for the non-linear mapping. The locations of maximum correlation for SAT and VAT are again different, and the multiple slices achieving maximum correlation are not contiguous. For non-linear mapping, the sites in the standardized anatomic space where maximum correlation is achieved are also listed in Table 21. One possible explanation for the findings is that discontinuous slice location combination may allow for a more representative sampling of the average fat area per slice across the abdominal region (vs. the scenario where all the slices are from contiguous slices through the abdomen).

TABLE-US-00025 TABLE 21 1 Slice 2 Slices 3 Slices Linear SAT 36 27, 53 26, 52, 54 mapping VAT 25 4, 27 3, 25, 29 Non-linear SAT 33 22, 37 22, 33, 35 mapping (T12) (L1, L2, T11) (L1-L2, T12, T11-T12) VAT 8 2, 36 8, 16, 36 (L3-L4) (L3, L4, T11) (L3-L4, L1-L2, T11)

[0151] FIG. 23 shows the correlation curves when multiple slices are used for correlation calculation. Here, only the results using 2 and 3 slices are shown since single slice results have been shown above. From FIGS. 19 and 23 one skilled in the art will observe that higher maximum correlation can be achieved when more slices are used for correlation calculation. The maximum correlation for VAT derived from non-linear mapping is substantially (10%) greater than that from linear mapping as shown in Table 20. Especially note that the correlation curves for VAT go down considerably for certain combinations of 2 and 3 slices for linear mapping. This is due to the fact that the multiple slices are in fact from much different anatomic locations among subjects as in the single-slice case, while the non-linear method performs much better since slices are selected anatomically more accurately in the standardized anatomic space.

Observations

[0152] Correlation analysis to determine the optimal anatomic slice locations in the abdomen for estimating body fat has not previously been performed. The inventors have found that the optimal anatomic slice locations for single-slice SAT and VAT estimation are not the same, contrary to common assumption. This result is important since these fat components may have different effects upon the pathophysiology of different disease processes. Use of multiple slices can achieve higher correlation than use of a single slice. The optimal locations of slices in this latter case are not contiguous. Experimental results on 50 abdominal CT image data sets showed that the standardized anatomic space created through non-linear mapping of slice locations achieves better anatomic localization than linear mapping. The method of the invention can be extended with greater or fewer landmarks than those adopted herein. The method has been illustrated by using CT image data sets, though the inventors will continue to explore the applicability of this method on MR image data sets in the future.

[0153] Overall, one skilled in the art will appreciate the following from the above detailed description of the methods of abdominal fat quantification in accordance with the invention: [0154] 1. The maximum area-to-volume correlation achieved is quite high, suggesting that it may be reasonable to estimate body fat by measuring the area of fat from a single anatomic slice at the site of maximum correlation. However, the site of maximum correlation and the degree of correlation itself may both depend on the particular patient group or disease condition studied. This disclosure focused on (near) normal male subjects in the age group of 50-60 years. [0155] 2. The site of maximum correlation is not at L4-L5 as commonly assumed, but is more superiorly located at T12-L1 for SAT and at L3-L4 for VAT. Furthermore, the optimal anatomic locations for SAT and VAT estimation are not the same, contrary to common assumption. [0156] 3. It is important to make sure that the slices for different subjects are selected at the same anatomic locations for correlation analysis. These locations seem to vary non-linearly from subject to subject, at least for the population (G) and body region number () considered herein. The standardized space mapping achieves this consistency of anatomic localization by accurately managing non-linearities in the relationships among landmarks. The dependence of VAT on the precision of anatomic localization seems to be far greater than that of SAT, perhaps due to the complex shape of the distribution of VAT compared to SAT.

[0157] Multiple slices achieve greater improvement in correlation for VAT than for SAT. The optimal locations of slices are not contiguous.

[0158] The method of abdominal fat quantification in accordance with the methods of the invention thus help one skilled in the art to find optimal location(s) of slices for any given patient group and body region utilizing the data sets under any given image modality. Once the optimal locations are determined in the manner demonstrated herein, actual acquisition of images at precisely those locations in clinical practice can be implemented without much difficulty by making appropriate changes to the scan protocol, for example by marking off plane locations on scout views.

[0159] One drawback of the described methods is that it is difficult to implement on MR images since it is quite challenging to segment vertebral bodies in MR images. However, if certain features to tag anatomic locations reliably can be identified on slice images, then the method can be implemented in a straightforward manner.

Example 2: Automatic Localization of IASLC-Defined Mediastinal Lymph Node Stations on CT Images using fuzzy Models

Materials and Methods

[0160] The image data sets utilized for this study were 45 routine contrast-enhanced chest CT examinations collected from the patient image database of the University of Pennsylvania health system. Subjects were male patients with an average age of 54.7±3.9 years. The images included were considered to be radiologically near normal by a board-certified radiologist (DAT). The CT examinations had been performed on 16 or 64 multi-detector row CT scanners (Siemens Medical Solutions, Malvern, Pa.) during a full inspiratory breath-hold and during the venous phase of enhancement following intravenous contrast administration. Each examination consisted of an average of 60 axial slices covering the entire thorax, with a pixel size of 0.77 mm×0.77 mm and a slice spacing of 5 mm. For model building and training, image data sets from 23 subjects were used, and the remaining 22 image data sets were utilized for testing the recognition method.

[0161] The AAR lymph node (AAR-LN) system utilizing the methods of the invention is comprised of two parts: model building and automatic recognition of the stations.

Model Building

[0162] The building of fuzzy models consists of three main processes: (a) gathering the image database, (b) delineating each nodal station on each 3D image to indicate the 3D region occupied by the station, and (c) constructing fuzzy models. Part (a) has been described in the aforementioned articles of Udupa et al. Part (b) describing how the regions were defined and how the delineation was implemented is presented below. Part (c) utilizes the same algorithms as used previously by the inventors for fuzzy modeling of anatomic organs, which take as input the delineated binary images of all samples of all organs for the population and output a hierarchical fuzzy anatomy model of the body region as described above with respect to the AAR system. In the present embodiment, the organs are replaced by the different mediastinal lymph node stations.

Delineating Lymph Node Stations

[0163] Each lymph node station was defined consistently according to the thoracic anatomic landmarks, generally as a 3D hexahedral object. The location and orientation of the faces of the hexahedron were defined in the anterior (A), posterior (P), right (R), left (L), superior (S), and inferior (I) aspects. An illustration of the CT images and the delineated lymph node stations is shown in FIG. 24. The stations were manually outlined according to a modified version of the IASLC definitions by trained specialists under the supervision of a board-certified radiologist (DAT). The modifications of the IASLC were necessary to make the specifications fully defined in terms of image features so as to facilitate its implementation and for handling situations where boundaries were ill-defined by a 2009 proposal for a lymph node map for lung cancer (Ichimura et al., “Caudal border of level 2R in the new international lymph node map for lung cancer,” Journal Thoracic Oncology, Vol. 5 No. 4, April, 2010, pp. 579-80 (2010)). Station 1 was adapted to the lower cervical region in the study as the full cervical region was not completely imaged on the standard chest CT examinations. In Table 22 below, the inventors present the regions that define the nodal stations, and their respective anatomical landmarks as used in the inventors' implementation.

TABLE-US-00026 TABLE 22 Description of the mediastinal lymph node stations via boundary definitions of the hexahedron, with the planes that define the hexahedral faces in the anterior (A), posterior (P), right (R), left (L), superior (S), and inferior (I) aspects. Region Definition Region Definition Station 1 A: Plane Station 4R A: Anterior wall Partial of between Right of superior vena cava Low anterior Lower or origin of left Cervical/ clavicles (*) Para- common carotid artery Supra- P: Anterior tracheal (whichever is clavicular aspect spine (*) Nodes more anterior) (*) Nodes R: Lateral tip of P: Posterior wall right transverse of trachea (*) process of R: Right pleural sac (*) spine (*) L: Left wall of trachea L: Lateral tip of S: Axial level/ left transverse plane where the inferior process aspect of left of spine (*) brachiocephalic S: Lung (innominate) vein apices (*) crosses anterior to the I: Superior left side of the aspect of trachea (*) manubrium (*) I: Inferior aspect of horizontal portion of azygos vein Station 2R A: Anterior wall Station 4L A: Anterior wall of Right of superior vena Left superior vena cava Upper cava or origin of Lower or origin of left Para- left common Para- common carotid artery tracheal carotid artery tracheal (whichever is Nodes (whichever is Nodes more anterior) (*) more anterior) P: Posterior wall (*) of trachea (*) P: Posterior wall R: Left wall of trachea (*) of trachea R: Right L: Oblique plane pleural sac (*) along left wall of S: Superior aortic arch (*) aspect of S: Superior aspect manubrium (*) of aortic arch I: Superior I: Superior aspect aspect of of left main aortic arch pulmonary artery Station 2L A: Anterior wall Station 5 A: Anterior wall Left of superior vena Subaortic of superior vena cava Upper cava or origin of Nodes or origin of left Para- left common common carotid artery tracheal carotid artery (whichever is Nodes (whichever is more anterior) (*) more P: Posterior aspect anterior) (*) of aortic arch (*) P: Posterior wall R: Oblique plane of trachea (*) along left wall of R: Left wall aortic arch (*) of trachea L: Left L: Left pleural sac (*) pleural sac (*) S: Inferior S: Superior aspect of aortic arch aspect of I: Superior manubrium (*) aspect of left main I: Superior pulmonary artery aspect of aortic arch Station 3a A: Posterior Station 6 A: Anterior aspect Pre- wall of sternum Para- of aortic arch (*) vascular P: Anterior wall aortic P: Posterior aspect Nodes of superior vena Nodes of aortic arch (*) cava or origin of R: Oblique plane left common along left wall of carotid artery aortic arch (*) (whichever is L: Left more pleural sac (*) posterior) (*) S: Superior aspect R: Lateral tip of of aortic arch a thoracic spinal I: Inferior aspect right transverse of aortic arch process + 2 cm (*) L: Lateral tip of a thoracic spinal left transverse process + 2 cm (*) S: Superior aspect of manubrium (*) I: Carina Station A: Posterior Station 7 A: Plane along 3p wall of trachea Sub- anterior wall of Retro- P: Anterior carinal mainstream bronchi tracheal aspect of Nodes P: Plane along Nodes spine (*) anterior aspect of R: Right wall the thoracic spine of trachea (*) R: Left wall of L: Left wall right bronchi (*) of trachea (*) L: Right wall of S: Superior left bronchi (*) aspect of S: Carina manubrium (*) I: Inferior aspect I: Carina of bronchus intermedius (on right) and superior aspect of left lower lobe bronchus (on left) (whichever is the most inferior) (*) (*) Modifiedfom or added to the original IASLC definitions [1] for consistency and for enabling actual computer implementation.

Constructing Fuzzy Models

[0164] The fuzzy anatomic model, FAM(B), of the stations in the thoracic body region B, is defined to be a quintuple FAM(B)=(H, M, ρ, λ, η), where H is a hierarchical order considered in FAM(B) for the organs and stations in B; M is a set fuzzy models FM(), one model for each organ/station ; ρ describes the position and orientation relationships, and their variations, between each offspring object O.sub.k and its parent in the hierarchy; λ is a family {1≤≤L}, where each expresses the variation in scale factor (size) of organ/station over its population; and η represents the statistics of a set of measurements pertaining to the organ/station assembly in B. Details on model building are described by above in the description of the AAR system.

Automatic Recognition of Organs and Lymph Node Stations

[0165] The inventors have experimented with automatic recognition, with hierarchy H.sub.A. The hierarchy is displayed in FIG. 25; H.sub.A has the thoracic skin (tskin) object as root and trachea and bronchi (tb) as offspring of tskin and parent of all stations. For the recognition process, the inventors used two different methods: thresholded optimal search and the one-shot method as described above in the description of the AAR system. As described therein, the thresholded optimal search algorithm refines object pose by an optimal search based on thresholding the test image. The one-shot method uses knowledge embedded in FAM(B) to recognize the different objects. The inventors combined both methods in a mixed strategy, with part of recognition done by using thresholded optimal search (for the objects tskin, tb and rs) and others by using the one-shot strategy (stations).

Evaluation

[0166] To evaluate recognition performance, the inventors used two metrics: distance error and scale ratio. Distance error is the distance between the geometric centers of the known true object and the fuzzy model at the time of recognition. Scale ratio is the size of the estimated model divided by the size of the true object. Object size is measured by the square root of the sum of the cigenvalues where the eigenvalues result from principal component analysis of the object.

Results

[0167] FIG. 26 depicts volume renditions of the fuzzy models of several stations and organs constructed from 23 of the 45 image data sets. The remaining 22 data sets were used for testing the recognition process.

[0168] In FIG. 27, the inventors display recognition results for three stations by overlaying the cross section of the fuzzy model at recognition on a slice of the test image. The true delineation of the stations is also displayed for reference on the corresponding slice.

[0169] In Table 23, the inventors present the recognition results for the different hierarchies. The results shown are the mean position error (in mm) for each object and scale ratio, respectively. In Table 23, one can observe that the anatomic organs (tskin, tb and rs) have mean errors from 5.27 to 7.13 mm. These values are comparable to the 5 mm spacing between slices, and are excellent results for automatic recognition. For recognition of the nodal stations, the location error ranges from 6.91 mm to 34.24 mm. Noting that the voxel size is limited by slice spacing, these errors expressed in terms of voxels are in the range of 1 to 6 voxels. Some stations such as #4R and #4L have a mean error of 7.55 mm and 6.91 mm, which is around 1 voxel. For station localization, the inventors believe that these results are excellent. The estimated scale ratio in the automatic recognition is perfect for tskin (1.00) and slightly underestimated for tb (0.90). The scale ratio for the nodal stations follows the ratio of the parent object, where the stations have a size ratio of around 0.90.

TABLE-US-00027 TABLE 23 Mean position error and standard deviation (in mm) and mean scale ratio. Station Metric Tskin Tb 1 2R 2L 3a 3p 4R 4L 5 6 7 Distance 5.2 7.13 17.02 11.73 18.51 11.79 9.73 7.55 6.91 13.95 34.24 13.77 error (std) (2.51) (5.20) (7.06) (6.73) (9.18) (7.70) (6.10) (2.88) (3.47) (5.96) (7.69) (4.55) [mm] Scale ratio 1.00 0.90 0.81 0.90 0.85 0.87 0.88 0.88 0.85 0.85 0.87 0.95

Observations

[0170] The automatic recognition of lymph node stations is a challenging problem. The proposed approach involves two main steps of AAR: fuzzy model building and object recognition. Definition of the lymph node stations was performed consistently and refined for computational specificity in order to obtain reliable computational fuzzy models. The recognition method was based on the one-shot and thresholded optimal search algorithms described above in the description of the AAR system. The results presented here are preliminary but indicate the potential utility of the AAR approach, originally designed for recognizing and delineating anatomy in medical imagery, adapted for automated lymph node station definition and localization. The results indicate that localization of thoracic nodal stations within 1-3 voxels is feasible for a majority of the tested stations. Considering the ambiguity (and fuzziness) that exists in the expert perception of the nodal stations, these results seem to be excellent.

[0171] The hierarchy used with the fuzzy anatomy model can impact the results considerably. The inventors are investigating this observation, especially keeping in mind the challenges offered by nodal stations 2R, 2L, and 6. Understanding the relationship among the objects of the thoracic region and lymph node stations as to which relationships are less variable over the patient population may allow the inventors to devise optimal hierarchies with improved recognition accuracy in the future.

[0172] While the IASLC formulation is helpful to standardize interpretation and reporting of lymph node disease conditions, it leaves the radiologist with the arduous task of following the detailed specifications and finding the lymph node stations on images subjectively, particularly since not all boundaries (superior, inferior, anterior, posterior, right, and left) are precisely defined for each of the lymph node stations. At present, to the inventors' knowledge, no method exists to assist radiologists during image interpretation to automatically identify and indicate IASLC lymph node stations on cross-sectional imaging. The inventors believe that when the AAR-LN system is implemented on dedicated workstations, the acceptance of the IASLC standard and the consistency of its interpretation will be greatly facilitated, not just for CT imaging assessment, but potentially also for magnetic resonance imaging (MRI), positron emission tomography (PET) imaging, and hybrid modalities such as PET/CT and PET/MRI for imaging assessment of the thorax. The inventors note, however, that other standards besides the IASLC standard may be readily used by those skilled in the art without departing from the teachings of the invention.

Example 3: Use of AAR for Radiotherapy Planning

[0173] The AAR approach described above is a significant departure from state-of-the-art segmentation methodologies in the following key considerations: [0174] (1) Fuzzy modeling: All reported model-based approaches have a statistical framework, none taking a fuzzy approach. Fuzzy set concepts have been used extensively otherwise in image processing and visualization. The AAR approach allows bringing anatomic prior information in an all-digital form into graph-based object delineation and recognition algorithms (such as fuzzy connectedness and graph cuts) without having to make “continuous” assumptions on matters such as shapes and random variables and their nature, density distribution, and independence, etc. They also allow capturing information about uncertainties at the patient level (e.g., blur, partial volume effects) and population level, and codification of this information within the model. [0175] (2) Prior information encoding: The AAR methodology takes a novel approach to gathering and encoding in an optimal hierarchical manner very detailed prior information about individual objects and their relationship to other objects, more globally as well as in their neighborhood. This information has a direct bearing on object recognition/localization and delineation. This codification permits gradual refinement starting from conspicuous global information to more subtle local details. [0176] (3) Numerousness of objects, generality of AAR: The AAR approach can handle all major organs in a body region, and is applicable to different body regions, and even image modalities. This generality allows rapid prototyping for a new application for the same body region, involving the same or different image modality. The extensive object-level knowledge that the AAR approach exploits at various stages for recognition and delineation is challenging to incorporate in prior art atlas-based methods and methods that rely on smooth spatial deformations and transformations.

[0177] Using the AAR methodology described herein for radiotherapy planning will shift current clinical practice by facilitating advanced RT methods such as IMRT/PBRT. The clinical benefits of adaptive planning have been well established in improving tumor radiation dosing and reducing dose to normal adjacent structures for head and neck and lung cancers. Despite promising findings showing less toxicity and greater sparing of organs at risk with adaptive planning for use of PBRT for non-small cell lung carcinoma (NSCLC), there is limited momentum in the field to use adaptive planning since it remains highly labor-intensive and impractical to deliver. The AAR approach described herein can be used to dramatically reduce contouring time, making adaptive re-planning much more feasible, facilitating its wide-spread clinical use.

Step 1: Build Hierarchical Fuzzy Anatomy Models (FAMs) of the Region of Interest

[0178] As noted above, a Fuzzy Anatomy Model of a body region B for a group G, FAM(B, G), is defined as a quintuple: FAM(B, G)=(H, , ρ, λ, η). H here denotes a hierarchy, represented as a tree, of the objects in B. As noted above, AAR uses anatomic hierarchies, an example of which is shown in FIG. 28. Methods are used in accordance with the methods of the invention to automatically construct a hierarchy that is optimal for object recognition, as explained below. Variable ={FM(O.sub.k): k=1, . . . , L} is a set of fuzzy models, with one model FM(O.sub.k) per each object O.sub.k. For constructing FM(O.sub.k) the inventors first convert all existing contour data into 3D binary images. Then the inventors follow the above-described approach of codifying as fuzzy membership values the manner in which the different samples of O.sub.k vary spatially over group G from an ideal homothetic affine transformation, while also retaining the spatial relationship among objects in the hierarchical order. Variable ρ describes the parent-to-offspring relationship in H over G in both position and orientation. Variable λ is a set of scale factor ranges indicating the size variation of each object O.sub.k over G. Variable η represents a set of measurements describing the morphology, shape, and image intensity properties of the objects in B over G. FAM(B, G) is then enhanced from this basic form in two ways as described below.

[0179] Optimal hierarchies: Consider a directed graph whose nodes are the 16 thoracic objects of Table 24 and in which every pair of objects is connected by a directed arc. Suppose a cost is assigned to each directed arc (x, y) that indicates how undesirable it is to have object x as the parent of y in the hierarchy. An optimal spanning tree for the graph will then express the desire of arranging objects in an optimal hierarchy. One possible form for cost is as a function of the spatial and intensity nearness of x and y—very near on both counts should imply high cost. If x is much larger in size than y and if both are similar in image intensity, then x is far less influenced by y during y's localization than vice versa and hence asymmetry of the cost function.

TABLE-US-00028 TABLE 24 Acronyms used for thoracic objects; see FIG. 28 SB Skin RS Respiratory TS Thoracic VS Venous Boundary System = Skeleton System LPS Left LPS + H Heart (main Pleural RPS + TB E Esophagus trunks Space BP Brachial AS Arterial only) RPS Right Plexus System SC Spinal Pleural S Stomach (main Cord Space IMS Internal trunks LBP Left TB Trachea & Mediastinum only) Brachial Bronchi Plexus RBP Right Brachial Plexus

[0180] Sparse object modeling: Sparse objects (such as trachea & bronchi—TB, esophagus—E) are challenging to model because small subject-to-subject variations in their form (e.g., the descending aorta portion of the arterial system (AS) being wavy versus straight) can lead to ill-defined models. This is true for all model-based approaches. The inventors propose to employ rough sets to construct fuzzy models FM(O.sub.k) for all sparse objects O.sub.k. The basic idea is as follows. Instead of using the precise boundary of the object in each training sample, a superset containing the region of the object is used (see FIG. 29). The supersets are defined carefully so that they are as small as possible, their variability over G is as small as possible, and they do not interfere with other nearby objects. The rest of the fuzzy model building process then will proceed as above using the supersets instead of the real object contours.

Step 2: Algorithms for the Automatic Localization and Delineation of Thoracic Objects

[0181] The AAR object recognition/localization procedure uses the model built in Step 1 and follows the procedure for AAR-R set forth above. Procedure AAR-R proceeds hierarchically in H. It finds the root object first (by calling R-ROOT) and subsequently localizes the offspring objects (by calling R-OBJECT) in a breadth-first order by combining optimally the prior information in FAM(B,G) with the information in the given image I. Different algorithms may be chosen for R-ROOT and R-OBJECT. Several such algorithms may be used, including the one-shot method, ball-scale strategy, Fisher Linear Discriminant, and Thresholded Optimal Search method. Among these, the Thresholded Optimal Search method described above is the top performer.

[0182] Methods for delineating objects require the above recognition results as the starting point. The inventors have chosen to base delineation methods on the Iterative Relative Fuzzy Connectedness (IRFC) engine since IRFC is extremely fast, completing segmentation of very large objects (5123) within 5 seconds. It is provably robust to seed set location and size, requires very few seeds, and is naturally adaptable to prior information encoded in the form of fuzzy models. Lastly, even when an object has poorly defined boundaries, it produces effective segmentations when co-object components are properly identified through seeds. Four extensions to the basic algorithm may also be provided as follows: (i) For each object, its topological neighbors will be determined and this anatomic information, encoded in FAM(B,G), is used for automatically specifying seeds for object and co-object components from knowledge of the recognition results for these components. (ii) For each object and its neighbors, their component tissues are identified and their characteristics (mean and standard deviation), determined from the training data sets, are encoded into FAM(B,G). This information is used for automatically and optimally specifying the fuzzy affinity, which forms the core of the IRFC (and any FC) engine. (iii) A fuzzy model component of affinity is designed to bring prior information seamlessly into the IRFC framework. (iv) The spatial, intensity-homogeneity-based, tissue-specific, and the model components of affinity are integrated into TRFC for a globally optimal separation of object and co-object components.

[0183] For segmenting images corresponding to different treatment fractions, the effectiveness of propagating recognition results from the previous fraction to the next for initializing recognition and delineation processes is also considered since the changes from fraction to fraction are usually not as much as variations found in the same object, with or without abnormalities, over a patient population.

[0184] The AAR methodology has the following strengths and unique features which set it apart from state-of-the-art object localization and delineation methods: (1) It starts with a precise and consistent definition of each body region and all major objects in each body region, which is then implemented in the AAR system. (2) Its generality has been demonstrated on four body regions—neck, thorax, abdomen, and brain—involving about 40 objects and on CT and MRI. (3) It employs a hierarchical arrangement of objects and encodes within this arrangement the non-linear relationship that exists among objects in their geographical layout and in their size over patient populations. (For example, the size of some objects (such as pericardium) changes much less within a patient group compared to others, and similarly the geometric relationship among objects.) Such detailed prior information is exploited within the hierarchy and in devising “optimal” hierarchies. (4) It integrates naturally the fuzzy modeling and fuzzy connectedness delineation approaches. Because of combinatorial optimization and graph-based methods, this integration results in vastly more efficient algorithms than popular methods in the literature such as those using atlases. (5) It allows rapid prototyping for different applications by modifying the hierarchy to suit the application.

[0185] The AAR software system described herein is preferably used for optimized radiation therapy planning in patients with malignancies treated with PBRT or other advanced adaptive RT methods. Those skilled in the art will appreciate that PBRT allows for ultra-precise delivery of treatment due to the physical characteristics of the proton beam and eliminates exit dose received by normal critical structures. However, small changes (on the order of mm's) in anatomy can result in under-dosing of tumor and overdosing of adjacent normal structures. Achieving optimal total dose delivery to the tumor while sparing normal tissue requires accounting for changes in patient anatomy during the course of treatment.

[0186] For example, in an exemplary embodiment, the AAR methodology may be used for radiation therapy planning by using the AAR methodology to automate the process of outlining objects and critical organs in the body region that is to receive radiation without the need for recontouring. For this purpose, the AAR methodology is applied to pretreatment images (CT, MRI, etc.) to identify objects and organs in the image and to segment such objects and organs as appropriate. The contours of the objects and organs are then provided as input to the radiation plan. Each time the patient visits for radiation therapy, the process is repeated to ascertain if there has been movement of objects or organs in the body region that is to receive radiation. The system creates a file of the object and organ contours in 3D space and provides such inputs into the radiation treatment planning software. In an alternative embodiment, the system may further show the distribution of a radiation dose simulated on the image having the contour lines for the images simultaneously illustrated so that the therapist may determine, before the dose is applied, whether the radiation dose is likely to impact tissue or an organ outside of the treatment area.

[0187] In an exemplary embodiment, the invention provides a computerized method of providing radiation therapy planning in patients receiving radiation treatment by: [0188] (a) building a fuzzy anatomy model of the body region of interest from an existing set of patient images for the body region; [0189] (b) obtaining a pretreatment image of a particular patient body region of interest that is to receive radiation therapy; [0190] (c) using automatic anatomy recognition (AAR) to recognize and delineate objects in the particular patient body region of interest; and [0191] (d) providing contours of delineated objects as input to radiation treatment planning software.
Steps (b)-(d) may be repeated prior to respective patient visits for radiation therapy in order to assess changes between visits.

[0192] FIG. 30 illustrates a system 50 that provides radiation therapy planning in patients receiving radiation treatment in accordance with such a method. As illustrated in FIG. 30, the system 50 nominally includes a database 60 that stores image data from patient image sets, an imaging device 70 that obtains a pretreatment image of a particular patient's body region of interest that is to receive radiation therapy, radiation treatment planning software 80 that controls a device 90 for applying radiation to the region of interest, a memory 100 storing computer instructions for implementing the method, and a processor 110 that processes the computer instructions to build a fuzzy anatomy model of the body region of interest from existing patient image sets, to recognize and to delineate objects in the particular patient's body region of interest using automatic anatomy recognition (AAR), and to provide contours of delineated objects as input to the radiation treatment planning software.

[0193] Those skilled in the art will appreciate that automating the radiation therapy planning process in this fashion saves time and improves accuracy so as to significantly reduce dosage toxicity and to potentially save a patient's organs and other tissues from the adverse side effects of radiation therapy. Those skilled in the art will further appreciate that this approach may be used with any kind of high energy radiation beam therapies.

[0194] It will be appreciated that all of 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 methods 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.

[0195] 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. For example, the techniques of the invention are not limited to CT images but may also be used on other medical imaging modalities such as Mill, ultrasound, PET/CT, etc. Also, the adaptive radiation therapy methods described herein include any such methods known in the art including intensity modulated ratiotherapy (IMRT) and proton beam radiation therapy (PBRT). Accordingly, any such modifications are intended to be included within the scope of this invention as defined by the following exemplary claims.