The present invention relates to image segmentation, and more particularly, to segmentation of tubular structures in images.
Tubular structures can appear as strokes or stroke-like structures in 2-dimensional images. As used herein, the term “stroke” refers to variable width curves in 2-dimensional images. Segmentation of stroke-like structures, such as blood vessels, is a fundamental problem in medical imaging, and is an important component of clinical applications involving diagnosis (e.g., stenosis, aneurysm, etc.), surgical planning, anatomical modeling and simulation, and treatment verification. Segmentation of stroke-like structures is a problem that also arises in other contexts including industrial applications and aerial/satellite image analysis.
Using manual segmentation techniques to segment stroke-like structures, it is possible to obtain highly accurate results. However, such manual techniques typically require too much tedious labor to be practical in clinical applications. Accordingly, various fully automatic segmentation methods have been developed for segmenting vessels and other stroke-like structures. However, due to poor contrast, noise, and clutter that is common to medical images, it is often difficult for fully automatic segmentation methods to yield robust results. Furthermore, one may be interested in extracting only a subset, for example a specific path through a branching network of stroke-like structures. Therefore, there is a salient need for an interactive segmentation method that is mostly automatic, but accepts input from an operator to guide the segmentation in a particular direction, quickly correct for errant segmentations, and add branches to an existing segmentation result. Computational efficiency is crucial to such a semi-automatic segmentation method, so that the operator will not have to wait for segmentation results while interacting with data. Typical conventional automatic segmentation techniques have runtimes that are too slow for use in such an interactive segmentation method.
The present invention provides a method for segmenting tubular or stroke-like structures in images. Embodiments of the present invention are directed to segmenting stroke-like structures in images using pearling. Pearling is the generation an ordered series of pearls, which are variable-radius 2D disks, as a discrete representation of the stroke geometry. Pearling is robust to fluctuations in image intensities (due to noise, etc.) as the forces acting on a pearl are integrated over the region inside the pearl. Pearling is computationally efficient and well suited to user interactivity. Such interactivity can allow operator guidance of the segmentation in a particular direction, as well as operator correction of errant segmentation results.
In one embodiment of the present invention, user inputs identifying a first region on the image inside of a tubular (stroke-like) structure and a second region of the image outside of the tubular structure are received. Based on the pixel intensities in the first and second regions, probability densities for inside and outside of the structure can be estimated. An ordered series of pearls are generated along the structure. Pearls are 2D disks, each having a center location and a radius determined based on local pixel intensities in the image. The probability densities are used to iteratively estimate the center point and the radius for each pearl. A continuous model of the structure is generated by interpolating the center locations and radii of the ordered series of pearls. The center locations and the radii of the ordered series of pearls can be interpolated using an iterative subdivision interpolation, which at each step introduces a new pearl between each pair of consecutive pearls in the ordered series of pearls. The continuous model is output as a segmentation result and can be used in interactive segmentation.
These and other advantages of the invention will be apparent to those of ordinary skill in the art by reference to the following detailed description and the accompanying drawings.
The present invention is directed to a method for segmenting tubular or stroke-like structures in 2D images. As used herein, the term “stroke” refers to variable width curves in 2-dimensional images. The terms tubular structure, stroke structure, and stoke-like structure are used herein interchangeably. Examples of such structures include, but are not limited to, blood vessels, bones, roads, rivers, electrical wirings, and brush-strokes.
Embodiments of the present invention are described herein to give a visual understanding of the segmentation method. A digital image is often composed of digital representations of one or more objects (or shapes). The digital representation of an object is often described herein in terms of identifying and manipulating the objects. Such manipulations are virtual manipulations accomplished in the memory or other circuitry/hardware of a computer system. Accordingly, is to be understood that embodiments of the present invention may be performed within a computer system using data stored within the computer system.
Embodiments of the present invention are directed to segmenting strokes (or stroke-like structures) in 2D images using pearling. Pearling extracts a high level parametric representation of a stroke in a 2D image. This representation is called a string.
where Ii is the intensity of the ith pixel in the region, m is the mean of intensities of the n pixels in the region, and h is the bandwidth of the estimator. A Gaussian kernel,
can be used. Performing this operation on the user input first and second regions results in first and second probability densities, {circumflex over (p)}in(I) and {circumflex over (p)}out(I), for the first region (inside the stroke structure) and the second region (outside the stroke structure), respectively.
At step 204, pearls are iteratively generated along the stroke-like structure in the 2D image. As defined herein, a pearl is a 2D disk having a radius and a center point. The radius of each pearl can be different, and the radius and center point of each pearl is determined based on local intensities at that pearl. An initial pearl is generated based on an initial seed c0. The initial seed c0 can be received as a user input or can be automatically determined. The initial seed c0 is a point in the stroke structure at which the initial pearl is generated. The initially seed c0 can be chosen to lie close to the centerline of the stroke structure and possibly at one end of the stroke. Based on the initial seed c0, an iterative method is used to generate an ordered series of pearls, one at a time. During this iterative method, at each step, the center ci and radius ri of the current pearl is calculated so as to maximize ri subject to the image data, given the constraint that the distance di between ci and the center ci−1 of the previous pearl is bound by functions dmin(ri−1,ri) and dmax(ri−1,ri). The linear functions dmin,max (ri−1,ri)=ari−1+bri can be used, with bounds dmin (ri−1,ri)=ri/2 and dmax (ri−1,ri)=ri−1+ri (i.e., the current pearl is tangent to the previous pearl). Since di can fluctuate, it is possible to capture rapid changes in the thickness of the stroke. After several iterations, the pearl converges to a local minimum. In addition to the initial seed point, it is possible that a direction be input by a user to indicate the direction from the initial seed point in which the stroke structure runs.
Referring to
ø(x) is a classification function that classifies a pixel l(x) as being inside the stroke or outside the stroke. The function f(ci,ri) sums the vectors (ci−x), where x is the vector coordinate of the current pixel, across the entire area of pixels for pearl Pi, using only those pixels such that ø(x)=1, i.e., pixels determined to be outside the stroke structure. Each of these vectors is weighted by its distance from ci, such that pixels nearer ci have a stronger influence on the result, as reflected in the
component of Equation (2). Intuitively, Equation (3) states that each point inside the ith pearl but outside the stroke imparts a force on the pearl that pushes the pearl away from the stroke boundary. When the forces are all balanced on the pearl, the pearl is typically centered in the stroke center. The
term Equation (2) is a normalization factor calculated based on the case where a pearl is cut into two equal halves by a linear boundary separating good (inside the stroke structure) and bad (outside the stroke structure) pixels. By assuming this ideal case, the offset is calculated that is needed to move the pearl half the minimum distance that would move it entirely in the stroke structure. Determination of whether a pixel lies inside the stroke structure or outside of the stroke structure is achieved based on the first and second intensity probability densities {circumflex over (p)}in (I) and {circumflex over (p)}out (I) determined using non-parametric density estimation. A pixel having an intensity I is classified as outside of the stroke structure if {circumflex over (p)}out (I)>{circumflex over (p)}in (I); otherwise it is classified as inside the stroke structure.
At step 304, the radius ri of the current pearl is estimated using g(ci,ri). The function g(ci,ri) is used to adjust ri to better fit the stroke. For robustness, pearls are designed to have a percentage p of their pixels inside the stroke structure and the rest of its pixels outside of the stroke structure. In a possible implementation, g(ci,ri) is negative when less than p percentage of the pearl's pixels are inside the stroke structure, and positive when less than 1−p percentage of the pearl's pixels are outside of the stroke structure. The portion of the enlarged pearls that lies outside of the stroke structure can be calculated by sampling the pearl's interior. Alternatively, it is possible the pearl's boundary can be sampled. The result of g(ci,ri) is used to scale ri to better fit the stroke structure. Accordingly, when the less than p percentage of the pearl's pixels are inside the stroke structure, ri is reduced, and when less than 1−p percentage of the pearl's pixels are outside of the stroke structure, ri is increased. g(ci,ri) can be expressed as:
At step 306, a constraint is applied on the distance between the current pearl and the previous pearl. The distance di between ci and the center ci−1 of the previous pearl is constrained by the functions dmin(ri−1,ri) and dmax(ri−1,ri). The linear functions dmin,max(ri−1,ri)=ari−1+bri can be used, with bounds dmin(ri−1,ri)=ri/2 and dmax(r−1,ri)=ri−1+ri. Using these bounds, at the maximum distance, the current pearl is tangent to the previous pearl.
At step 308, it is determined whether the center location ci and the radius ri of the current pearl have converged. The estimation of ci and ri are interleaved for the current pearl. For a given position ci and radius ri, both parameters can be updated independently using the results given by f(ci,ri) and g(ci,ri), respectively. In both cases, the quality of the fit is measured and an adjustment is calculated and returned. It is possible the center location ci and the radius ri are determined to converge when f(ci,ri) and g(ci,ri) fall below a threshold value. If the center location ci and the radius ri of the current pearl have not converged, the method returns to step 302 and repeats until convergence of both ci and ri. The adjustments of ci and ri at each iteration are calculated while enforcing the constraint on di. When the center location ci and the radius ri of the current pearl have converged, the values of ci and ri are saved for that pearl, and the method of
The method of
In applications in which the method of
A loop can also exist in the stroke structure. An intersection check between the current pearl and all previous pearls will reveal if such a loop has occurred. In this case, the current pearl is linked to the intersected pearl and its string is ended. The intersected pearl can then become a bifurcation pearl, if it already has more than one adjacent pearl, or it could be the end of another stroke still in the growing phase. If the latter is the case, and two growing strokes collide, they are both ended and become part of the same stroke.
Returning to
Note that the formulation of the four-point interpolation corresponds to a local cubic fit. It is obtained by solving for the coefficient d=F(0) of a cubic function F(u)=au3+bu2+cu+d, given the constraints: F(−3)=q(Pa), F(−1)=q(Pb), F(1)=q(Pc), and F(3)=q(Pd). Strings for which the beginning pearl is not equal to the ending pearl are extended by two pearls so that the subdivision has enough region of support to calculate Equation (5) at the end pearls.
As described above, subdivision can be used to generate the continuous model of the series of pearls, but the present invention is not limited thereto. For example, the continuous model may be generated by fitting a piecewise polynomial model to the pearl centers and the radii of the series of pearls, or by fitting a tangent curve to consecutive pearls in the ordered series of pearls, on opposite sides of each pearl corresponding to the stroke boundary.
Since the values of ri are larger than the actual stroke structure, as described above, the values can be adjusted to better fit the stroke structure. The radii are calculated to have a percentage p of pixels inside the stroke structure and the rest outside of the stroke structure. Accordingly, all values of ri can be multiplied by p, as used in calculating g(ci,ri), to get the true radii which can be used in the continuous model W. This reduced radius defines the core of each pearl. The area between the reduced radius and the original radius for each pearl is the crust.
Returning to
For more control, a user may provide a rough trace of the centerline of a stoke-like structure in an image, rather than providing the initial seed point and a direction. As the user is manually tracing, pearls are generated at samples along the curve and converge in real time, thus forming the stroke as the user is tracing it with a stylus. In this mode, the user can use a stylus to quickly trace a rough curve through the desired stroke structure, and then add braches as desired. As the operator is still tracing the curves, the pearling segmentation method can calculate the idealized strokes (centerline and radius) and display them. Each new trace can then either add or remove a branch.
The above-described method steps of
The foregoing Detailed Description is to be understood as being in every respect illustrative and exemplary, but not restrictive, and the scope of the invention disclosed herein is not to be determined from the Detailed Description, but rather from the claims as interpreted according to the full breadth permitted by the patent laws. It is to be understood that the embodiments shown and described herein are only illustrative of the principles of the present invention and that various modifications may be implemented by those skilled in the art without departing from the scope and spirit of the invention. Those skilled in the art could implement various other feature combinations without departing from the scope and spirit of the invention.
This application claims the benefit of U.S. Provisional Application No. 60/887,888, filed Feb. 2, 2007, the disclosure of which is herein incorporated by reference.
Number | Name | Date | Kind |
---|---|---|---|
6754376 | Turek et al. | Jun 2004 | B1 |
7116810 | Miller et al. | Oct 2006 | B2 |
20020168110 | Al-Kofahi et al. | Nov 2002 | A1 |
20040114800 | Ponomarev et al. | Jun 2004 | A1 |
20050110791 | Krishnamoorthy et al. | May 2005 | A1 |
20060159326 | Rasche et al. | Jul 2006 | A1 |
20060251304 | Florin et al. | Nov 2006 | A1 |
20070216678 | Rouet et al. | Sep 2007 | A1 |
20090129671 | Hu et al. | May 2009 | A1 |
Number | Date | Country |
---|---|---|
WO 2004104939 | Dec 2004 | WO |
Number | Date | Country | |
---|---|---|---|
20080187197 A1 | Aug 2008 | US |
Number | Date | Country | |
---|---|---|---|
60887888 | Feb 2007 | US |