Region-growing algorithm

Information

  • Patent Grant
  • 9595111
  • Patent Number
    9,595,111
  • Date Filed
    Tuesday, May 12, 2015
    9 years ago
  • Date Issued
    Tuesday, March 14, 2017
    7 years ago
Abstract
A technique for automatically generating a virtual model of a branched structure using as an input a plurality of images taken of the branched structure. The technique employs an algorithm that avoids inaccuracies associated with sub-optimal threshold settings by “patching” holes or leaks created due to the inherent inconsistencies with imaging technology. By “patching” the holes, the algorithm may continue to run using a more sensitive threshold value than was previously possible.
Description
BACKGROUND OF THE INVENTION

Navigating the airways of the lungs has always presented challenges to physicians attempting to diagnose and treat lesions transluminally. As such, numerous navigational aids and imaging tools have been developed and/or utilized to provide physicians a “map” of the lungs.


One such tool is a CT scanner. CT scanners use X-ray technology to take multiple scans or “slices” of the lungs. These scans each represent a cross-section of the lungs and can be viewed individually or assembled, via computer programs, to form a three-dimensional CT model. However, CT scans, like most images using X-ray technology, are somewhat cloudy and translucent in nature and difficult to view. As such, computer graphics techniques are employed to interpret the information provided by the CT model and “grow” a virtual model of the airways which mimics what might be seen by a bronchoscope navigating the airways. An example of this process is shown and described in U.S. patent application Ser. No. 11/939,537, entitled Adaptive Navigation Technique For Navigating A Catheter Through A Body Channel Or Cavity, the entirety of which is incorporated by reference herein.


This graphical technique is sometimes referred to as “region growing or 3D map generation,” and presents its own challenges. For example, region growing typically involves a processing of the CT data by analyzing each two-dimensional pixel, or, more pertinently, three-dimensional voxel, for brightness or “density” and assigning the voxel a value that indicates either tissue or air based on whether the density meets a certain threshold value. CT scans are grayscale images composed of a plurality of pixels (2D) or voxels (3D—if the scans have been assembled into a volume), each pixel or voxel varying in intensity from white (most intense) to black (least intense). Each intensity level between white and black appears as a shade of gray. By designating the various shades of gray from the CT scans either “tissue” or “air” the resulting image of the lungs becomes much more clear. However, if the voxels are designated incorrectly, the model of the lungs becomes inaccurate. Incorrect voxel designation results from setting the threshold level at an incorrect value, which is an inherent problem when attempting to assign discreet values (air or tissue) to voxels which are actually various shades of gray.


A presently-used technique for optimal threshold setting is beginning with a conservative threshold and performing a region-growing iteration. A conservative threshold is one that is not likely to result in leakage, which occurs when tissue is designated as air and creates a virtual image of the airways that looks like air (color) is spilling out of the airways. However, even with a conservative threshold, inaccuracies in the CT scans can result in “holes” after the segmentation process. These holes result in false branches.


Moreover, a conservative threshold results in airways that end prematurely. Therefore, after a conservative iteration is performed, resulting in a stunted branched structure, the threshold is incrementally increased and a second iteration is performed. These steps are repeated until leakage occurs. Thus, the maximum threshold that does not result in leakage is determined and used. This approach, however, naturally results in the least-dense portion of the CT image dictating the threshold level. There are other problems that arise from this method as well.


For example, during a full inhalation, the airways stretch and thin in order to accommodate the additional air volume. The thinning of the airways results in reduced tissue imaging density, and leakage thus arises even at lower threshold values.


Another example is that each time the threshold is increased, the algorithm runs from the initial seed point. Hence, if the threshold is increased ten times before leakage arises, each of the voxels analyzed in the initial iteration, is analyzed nine more times. This algorithm is thus inherently taxing on processing resources.


As such, a need is identified for a region-growing algorithm that is able to identify localized weaknesses in image data and “repair” them such that more distal branches of a bronchial tree can be segmented and “grown”.


SUMMARY OF THE INVENTION

The present invention provides a technique and algorithm for controlling leakage while using a high threshold during a region growing procedure. Rather than beginning with a conservative (low) threshold, a high threshold is used but growing is analyzed and controlled.


One aspect of the present invention uses an algorithm that defines iteration as a single voxel layer of growth, rather than growth until leakage is detected. Each of the iterations is analyzed for a predicted increase (e.g., parabolic increase) in the number of voxels. If the voxels are not growing according to an expected rate, the unexpected increase is identified as leakage.


Rather than ending the region-growing process when leakage is detected, the point of leakage is isolated or “patched” so that further iterations on the current path do not continue from that point, but the rest of the model is allowed continued growth.


Another aspect of the present invention provides a post-iteration step which, in effect, repairs or salvages branches encountering leakage. After all of the iterations are performed, branches that had leakage are retrieved from a stored location and regrown with a reduced threshold. These branches are then added to the model.


The advantages of the algorithm of the present invention are numerous. One advantage is a more meaningful increase in the number of branches as a result of the utilization of the highest possible threshold.


Another advantage of the present invention is that the appearance of “holes” in the central segmented area of the model is greatly reduced.


Yet another advantage of the present invention is a reduced number of false branches. Because each iteration is analyzed, growth is controlled and false branches are prevented from growing.


Another advantage of the present invention is that the branches grown have an increased length, due to an optimized threshold as well as secondary growing efforts.


Still another advantage of the present invention is reduced importance on initial threshold value selection. Through experimentation, it has been found that this algorithm is self-correcting, yielding virtually identical results regardless of the initial threshold used, so long as the thresholds are relatively high.





BRIEF DESCRIPTION OF THE DRAWINGS


FIG. 1 is a flow-chart of a preferred embodiment of an algorithm of the present invention.





DETAILED DESCRIPTION OF THE INVENTION

The technique 10 of the present invention is charted in the flowchart presented as FIG. 1. The technique 10 begins at 20 by selecting a starting point for the segmentation of the CT data. For example, selecting a point in the trachea is a logical starting point as it is the largest, and most proximal airway, and easily recognizable on a CT scan. Preferably, the starting point is a single voxel inside the trachea that is centered and meets a high threshold value. In other words, a voxel which is clearly air is selected. Any point from which further adjacent voxels are analyzed is hereinafter referred to as a “seed point”. The starting point is the first seed point.


At 30 the propagation process is initiated by designating adjacent voxels around the starting point. This designation is known as segmentation, and it indicates that the new voxels met the threshold level for air. Because the starting point is preferably selected in the trachea, and it is not desired to grow the airway back toward the mouth of the patient, only voxels in a distal direction of the starting point are segmented. As such, a “wave” is generated that travels distally into the airways. It is understood that the branches of the lungs fan out in all directions. As such, “distal direction” as used herein is interpreted as getting further away from the starting point along a path that remains in the airway. In other words, some “distal” points in the airways may be relatively close to the starting point if one were to cut through tissue to make a straight line between the points.


It is understood that the starting point, being a single voxel, is surrounded by 26 adjacent voxels, 17 of which are extending in a direction of desired growth (in a direction not yet segmented and not in a reverse direction, such as toward the mouth).


At 40, it is determined whether any new voxels were segmented. In other words, if no voxels in that “wave” met the threshold level for air, there is no growth and the tree is complete.


If the answer at 40 is “yes”, the algorithm continues to 50, where the segmented voxels are counted and analyzed to determine whether there is leakage. The new voxels may be more numerous than the previous iteration, but the growth should be controlled. In other words, if the increase in the number of voxels is expected (more specifically, a parabolic increase has been observed in normal growth patterns without leakage) then it is determined that there is no leakage and the process returns to step 30 to begin a new iteration. For example, as indicated above, the starting point consisted of one voxel. Assuming the starting point was surrounded by air, the next “wave” of voxels would be seventeen new seed points. However, because many of these are adjacent to each other, the next successive wave would not give rise to seventeen new seed points for each of the previous seventeen seed points. Additionally, the voxels behind each seed point that have been already analyzed, are not segmented again. If the airway being segmented were perfectly cylindrical, as soon as the seed points reached the walls of the airway, the “wave-front” would be a convex sheet, a single voxel in thickness, and would remain constant. Hence, the mathematical model for growing is somewhat parabolic, except when new branches are introduced, and considering that the airways are narrowing in the distal direction. Leakage, however, results in an abrupt increase in the number of segmented voxels.


If at 50 the analysis results in an unexpected or abnormally high increase in segmented voxels, it is determined that leakage exists and the process moves to step 60, which identifies and records the segmented voxels from the previous iteration and labels them as accurate.


Leakage determination is derived from two important conclusions: (1) It is expected that the front size has an inverse (not necessarily linear) dependence on the iteration number, e.g. [front size]˜(1/iteration number). (2) Bifurcations and changes in airway shape may result in somewhat linear growth in front size.


At 70, upon the detection of leakage, an analysis is done in order to select the most recent “good” iteration that does not contain leakage. This decision is based on the same principals used to satisfy the mathematical model in compliance with the natural structure of the bronchial tree.


The voxels segmented up to the last good iteration are assumed to be committed to the segmented voxel list, while voxels that belong to each iteration above the “good” one are analyzed in order to separate the voxels that led to the leakage. In order to make this analysis, recently segmented voxels are labeled to connected objects (part of the branch). Each object is then analyzed in order to detect the object that caused leakage.


The coordinates of voxels that belong to inaccurate objects are stored separately and treated differently thereafter. They are locked in the binary map to prevent their participation in the segmentation process. Voxels, belonging to accurate branches are returned to the segmentation and the process returns to 30.


At 80 the objects identified as leaking are removed from further segmentation and stored in a queue for post processing. After 80 the process returns to 30 for the next iteration.


If at 50 the answer was “no” for leakage, the process returns to 30 for the next iteration. It should be noted that the flow chart 10, though presented in series for clarification, is actually a parallel operation. In other words, each new voxel is a seed point and the flow chart is performed on each next iteration therefrom simultaneously. Hence, viewing the growth of the bronchial tree real time, one would see a near-instant tree appear, depending of course on the power of the processor running the algorithm.


If at 40 there are no new voxels detected, the algorithm proceeds to 90, which is a decision step asking whether any leakage objects were identified. If the answer is “yes” step 100 is executed, which retrieves the last object from the storage queue (step 80).


Next at 110, the threshold is reduced and the algorithm is performed on only the selected leakage object. Because a reduced threshold is being used, the likelihood of leakage is reduced.


If at 90 it is determined that there are no leakage objects, either because there were none or they have all been reprocessed, the process is completed at 120.


Although the invention has been described in terms of particular embodiments and applications, one of ordinary skill in the art, in light of this teaching, can generate additional embodiments and modifications without departing from the spirit of or exceeding the scope of the claimed invention. Accordingly, it is to be understood that the drawings and descriptions herein are proffered by way of example to facilitate comprehension of the invention and should not be construed to limit the scope thereof.

Claims
  • 1. A method for performing a region growing process via a region growing algorithm, the method comprising: comparing voxels adjacent to a first seed voxel with a predetermined threshold level for air;identifying the adjacent voxels representing interior cavities or lumens of an anatomical structure, as seed voxels, when each of the adjacent voxels meets the predetermined threshold level;iterating the identification of seed voxels a plurality of times with the predetermined threshold level;analyzing each of the plurality of iterations for an increase in a number of voxels;determining whether the increase in the number of voxels exceeds a predefined rate; andisolating one or more voxels of the plurality of voxels when it is determined that the increase in the number of voxels exceeds the predefined rate without terminating the region growing process.
  • 2. The method of claim 1, wherein each of the plurality of iterations includes segmenting voxels.
  • 3. The method of claim 2, wherein the segmenting of the voxels is commenced at the first seed voxel.
  • 4. The method of claim 3, wherein the first seed voxel is a point in a trachea.
  • 5. The method of claim 4, wherein the first seed voxel in the trachea meets a high threshold value.
  • 6. The method of claim 1, further comprising: identifying leakage, wherein leakage is identified when it is determined that the increase in the number of voxels exceeds the predefined rate.
  • 7. The method of claim 6, wherein the number of the segmented voxels is counted and analyzed to determine leakage levels.
  • 8. The method of claim 7, wherein the segmented voxels are identified and recorded.
  • 9. The method of claim 8, wherein the segmented voxels from each of the plurality of iterations are used to create a voxel list.
  • 10. The method of claim 9, wherein the voxel list is used to create a virtual model of the anatomical structure.
  • 11. A method for creating a virtual model of an anatomical structure, the method comprising a region growing algorithm for: comparing voxels adjacent to a first seed voxel with a predetermined threshold level for air;identifying the adjacent voxels representing interior cavities or lumens of an anatomical structure, as seed voxels, when each of the adjacent voxels meets the predetermined threshold level;iterating the identification of seed voxels a plurality of times with the predetermined threshold level;analyzing each of the plurality of iterations to determine a rate of change in a number of voxels;determining whether an increase in the number of voxels exceeds at a predefined rate; andisolating one or more voxels of the plurality of voxels when it is determined that the increase in the number of voxels exceeds the predefined rate without terminating the region growing process.
  • 12. The method of claim 11, further comprising segmenting voxels.
  • 13. The method of claim 12, further comprising commencing the segmenting of the voxels at first seed voxel.
  • 14. The method of claim 13, wherein the first seed voxel is a point in a trachea.
  • 15. The method of claim 14, wherein the first seed voxel in the trachea meets a high threshold value.
  • 16. The method of claim 15, further comprising identifying leakage when it is determined that the increase in the number of voxels exceeds the predefined rate.
  • 17. The method of claim 16, further comprising counting and analyzing the number of the segmented voxels to determine leakage levels.
  • 18. The method of claim 17, further comprising identifying and recording the number of the segmented voxels.
  • 19. The method of claim 18, further comprising creating a voxel list from the segmented voxels.
  • 20. The method of claim 19, further comprising using the voxel list to create the virtual model of the anatomical structure.
CROSS-REFERENCE TO RELATED APPLICATIONS

This present application is a continuation application of U.S. application Ser. No. 14/492,561 filed Sep. 22, 2014, which is a continuation application of U.S. application Ser. No. 13/867,217 filed Apr. 22, 2013, now U.S. Pat. No. 8,842,898, which is a continuation application of U.S. application Ser. No. 13/019,261 filed Feb. 1, 2011, now U.S. Pat. No. 8,428,328, which claims priority to U.S. Provisional Application Ser. No. 61/300,423 filed Feb. 1, 2010, which are hereby incorporated herein by reference in their entirety.

US Referenced Citations (154)
Number Name Date Kind
4346384 Raab Aug 1982 A
4791934 Brunnett Dec 1988 A
5383454 Bucholz Jan 1995 A
5386828 Owens et al. Feb 1995 A
5402801 Taylor Apr 1995 A
5417210 Funda et al. May 1995 A
5425367 Shapiro et al. Jun 1995 A
5445166 Taylor Aug 1995 A
5472441 Edwards et al. Dec 1995 A
5480422 Ben-Haim Jan 1996 A
5558091 Acker et al. Sep 1996 A
5572999 Funda et al. Nov 1996 A
5600330 Blood Feb 1997 A
5630431 Taylor May 1997 A
5668844 Webber Sep 1997 A
5695500 Taylor et al. Dec 1997 A
5704361 Seward et al. Jan 1998 A
5729129 Acker Mar 1998 A
5730129 Darrow et al. Mar 1998 A
5738096 Ben-Haim Apr 1998 A
5740802 Nafis et al. Apr 1998 A
5749362 Funda et al. May 1998 A
5752513 Acker et al. May 1998 A
5772594 Barrick Jun 1998 A
5797849 Vesely et al. Aug 1998 A
5800352 Ferre et al. Sep 1998 A
5810007 Holupka et al. Sep 1998 A
5829444 Ferre et al. Nov 1998 A
5840025 Ben-Haim Nov 1998 A
5871445 Bucholz Feb 1999 A
5873822 Ferre et al. Feb 1999 A
5891034 Bucholz Apr 1999 A
5902239 Buurman May 1999 A
5913820 Bladen et al. Jun 1999 A
5919188 Shearon et al. Jul 1999 A
5928248 Acker Jul 1999 A
5944023 Johnson et al. Aug 1999 A
5976127 Lax Nov 1999 A
6016439 Acker Jan 2000 A
6019724 Gronningsaeter et al. Feb 2000 A
6032675 Rubinsky Mar 2000 A
6073043 Schneider Jun 2000 A
6077257 Edwards et al. Jun 2000 A
6115626 Whayne et al. Sep 2000 A
6147480 Osadchy et al. Nov 2000 A
6149592 Yanof et al. Nov 2000 A
6161032 Acker Dec 2000 A
6188355 Gilboa Feb 2001 B1
6201387 Govari Mar 2001 B1
6203493 Ben-Haim Mar 2001 B1
6211666 Acker Apr 2001 B1
6226543 Gilboa et al. May 2001 B1
6233476 Strommer et al. May 2001 B1
6236875 Bucholz et al. May 2001 B1
6289235 Webber et al. Sep 2001 B1
6314310 Ben-Haim et al. Nov 2001 B1
6331116 Kaufman et al. Dec 2001 B1
6332089 Acker et al. Dec 2001 B1
6335617 Osadchy et al. Jan 2002 B1
6366799 Acker et al. Apr 2002 B1
6373240 Govari Apr 2002 B1
6405072 Cosman Jun 2002 B1
6427314 Acker Aug 2002 B1
6453190 Acker et al. Sep 2002 B1
6484118 Govari Nov 2002 B1
6556696 Summers et al. Apr 2003 B1
6558333 Gilboa et al. May 2003 B2
6574498 Gilboa Jun 2003 B1
6580938 Acker Jun 2003 B1
6591129 Ben-Haim et al. Jul 2003 B1
6593884 Gilboa et al. Jul 2003 B1
6615155 Gilboa Sep 2003 B2
6618612 Acker et al. Sep 2003 B1
6650927 Keidar Nov 2003 B1
6690816 Aylward et al. Feb 2004 B2
6690963 Ben-Haim et al. Feb 2004 B2
6702780 Gilboa et al. Mar 2004 B1
6706041 Costantino Mar 2004 B1
6711429 Gilboa et al. Mar 2004 B1
6788967 Ben-Haim et al. Sep 2004 B2
6833814 Gilboa et al. Dec 2004 B2
6947788 Gilboa et al. Sep 2005 B2
6995729 Govari et al. Feb 2006 B2
6996430 Gilboa et al. Feb 2006 B1
7176936 Sauer et al. Feb 2007 B2
7197354 Sobe Mar 2007 B2
7233820 Gilboa Jun 2007 B2
7236567 Sandkamp et al. Jun 2007 B2
7286868 Govari Oct 2007 B2
7301332 Govari et al. Nov 2007 B2
7321228 Govari Jan 2008 B2
7324915 Altmann et al. Jan 2008 B2
7343195 Strommer et al. Mar 2008 B2
7353125 Nieminen et al. Apr 2008 B2
7366562 Dukesherer et al. Apr 2008 B2
7370656 Gleich et al. May 2008 B2
7373271 Schneider May 2008 B1
7386339 Strommer et al. Jun 2008 B2
7397364 Govari Jul 2008 B2
7517318 Altmann et al. Apr 2009 B2
7822461 Geiger Oct 2010 B2
7831076 Altmann et al. Nov 2010 B2
8165385 Reeves Apr 2012 B2
8218847 Averbuch et al. Jul 2012 B2
8218848 Lenglet Jul 2012 B2
8842898 Averbuch et al. Sep 2014 B2
9042625 Averbuch et al. May 2015 B2
20010031919 Strommer et al. Oct 2001 A1
20020022837 Mazzocchi et al. Feb 2002 A1
20020065461 Cosman May 2002 A1
20020082498 Wendt et al. Jun 2002 A1
20020193686 Gilboa Dec 2002 A1
20030074011 Gilboa et al. Apr 2003 A1
20030099390 Zeng et al. May 2003 A1
20030144658 Schwartz et al. Jul 2003 A1
20030216639 Gilboa et al. Nov 2003 A1
20040006268 Gilboa et al. Jan 2004 A1
20040019350 O'Brien et al. Jan 2004 A1
20040097804 Sobe May 2004 A1
20040138548 Strommer et al. Jul 2004 A1
20040215181 Christopherson et al. Oct 2004 A1
20040249267 Gilboa Dec 2004 A1
20040254454 Kockro Dec 2004 A1
20050033149 Strommer et al. Feb 2005 A1
20050090818 Pike et al. Apr 2005 A1
20050107688 Strommer May 2005 A1
20050182295 Soper et al. Aug 2005 A1
20050197566 Strommer et al. Sep 2005 A1
20060058647 Strommer et al. Mar 2006 A1
20060064006 Strommer et al. Mar 2006 A1
20070167738 Timinger et al. Jul 2007 A1
20070167743 Honda et al. Jul 2007 A1
20070167806 Wood et al. Jul 2007 A1
20070232896 Gilboa et al. Oct 2007 A1
20070287901 Strommer et al. Dec 2007 A1
20080008367 Franaszek et al. Jan 2008 A1
20080033452 Vetter et al. Feb 2008 A1
20080086051 Voegele Apr 2008 A1
20080097187 Gielen et al. Apr 2008 A1
20080118135 Averbuch May 2008 A1
20080132909 Jascob et al. Jun 2008 A1
20080132911 Sobe Jun 2008 A1
20080139915 Dolan et al. Jun 2008 A1
20080147000 Seibel et al. Jun 2008 A1
20080157755 Kruger et al. Jul 2008 A1
20080161682 Kendrick et al. Jul 2008 A1
20080162074 Schneider Jul 2008 A1
20080183071 Strommer et al. Jul 2008 A1
20080183073 Higgins Jul 2008 A1
20080188749 Rasche et al. Aug 2008 A1
20080247622 Aylward et al. Oct 2008 A1
20090182224 Shmarak et al. Jul 2009 A1
20100016658 Zou et al. Jan 2010 A1
20100063410 Avila Mar 2010 A1
Foreign Referenced Citations (4)
Number Date Country
2007-312837 Dec 2007 JP
9605768 Feb 1996 WO
9625882 Aug 1996 WO
2008125910 Oct 2008 WO
Non-Patent Literature Citations (15)
Entry
WIPO, U.S. International Search Authority, International Search Report mailed Sep. 24, 2009 in International Patent Application No. PCT/IL2009/000569.
WIPO, U.S. International Search Authority, International Search Report mailed Nov. 13, 2008 in International Patent Application No. PCT/IB07/04567.
WIPO, U.S. International Search Authority, International Search Report mailed Oct. 26, 2009 in International Patent Application No. PCT/IB2009/005609.
WIPO, International Search Report mailed Mar. 31, 2000 in International Patent Application No. PCT/US99/26826.
WIPO, International Search Report mailed Aug. 14, 2001 in International Patent Application No. PCT/IL01/00224.
Palagyi et al., “Quantitative analysis of pulmonary airway tree structures”, Computers in Biology and Medicine 36 (2006), pp. 974-996.
WIPO, U.S. International Search Authority, International Search Report in International Patent Application No. PCT/IB11/00243.
Japanese Search Report dated Mar. 18, 2014 for JP Appln. No. 2012-550530.
Australian Patent Examination Report No. 1 for AU 2014-201021 dated May 22, 2015.
Japanese Notice of Allowance for JP 2014-229546 dated Apr. 12, 2016.
Huang et al., “RGVis: Region Growing Based Techniques for Volumn Visualisation”, Proceedings of the 11th Pacific Conference on Computer Graphics and Applications, 2003, 9 pages.
Pinho et al., “Robust Region Growing Based Intrathoracic Airway Tree Segmentation”, Proc. of Second International Workshop on Pulmonary Image Analysis, 2009, 11 pages.
Sherbondy et al., “Fast volume segmentation with simultaneous visualization using programmable graphics hardware”, In Visualization, 2003, VIS 2003, IEEE (pp. 171-176), 6 pages.
Australian Patent Examination Report in connection with application No. AU2015238846, dated Jul. 19, 2016, 4 pages.
Japanese Office Action for Japanese Appln. No. JP 2014-229546 dated Oct. 15, 2015.
Related Publications (1)
Number Date Country
20150243042 A1 Aug 2015 US
Provisional Applications (1)
Number Date Country
61300423 Feb 2010 US
Continuations (3)
Number Date Country
Parent 14492561 Sep 2014 US
Child 14710113 US
Parent 13867217 Apr 2013 US
Child 14492561 US
Parent 13019261 Feb 2011 US
Child 13867217 US