System and method for 3D vessel segmentation with minimal cuts
View Patent ↗A method for segmenting tubular structures in digital medical images includes extracting a subregion from a 3-dimensional (3D) digital medical image volume containing a vessel of interest, identifying potential vessel centerpoints for each voxel in the subregion by attaching to each voxel a tip of a 3D cone that is oriented in the direction of the voxel's image gradient and having each voxel within the cone vote for those voxels most likely to belong to a vessel centerline, selecting candidates for a second vote image that are both popular according to a first vote image, as well as being consistently voted upon by a radius image, reconfiguring the subregion as a graph where each voxel is represented by a node that is connected to 26 nearest neighbors by n-link edges, and applying a min-cut algorithm to segment the vessel within the subregion.
1. A computer implemented method for segmenting tubular structures in digital medical images, said method executed by a computer comprising the steps of:
extracting a subregion from a 3-dimensional (3D) digital medical image volume containing a vessel of interest;
for each voxel in the subregion, identifying potential vessel centerpoints by attaching to each voxel a tip of a 3D cone that is oriented in the direction of the voxel's image gradient and having each voxel within the 3D cone vote for those voxels most likely to belong to a vessel centerline;
selecting candidates for a second vote image that are both popular according to a first vote image, as well as being consistently voted upon by a radius image;
reconfiguring the subregion as a graph wherein each voxel is represented by a node that is connected to 26 nearest neighbors by n-link edges; and
applying a min-cut algorithm to segment the vessel within the subregion.
2. The method of claim 1 , wherein the segmented vessel includes branches.
3. The method of claim 1 , further comprising detecting intersections of the segmented vessel with the subregion's boundaries, wherein intersection voxels serve as starting voxels for extracting a new subregion overlapping with the current subregion, wherein said new subregion serves to track the vessel's centerline.
4. The method of claim 1 , wherein said subregion is a cubic subregion of size
x
i
=
[
-
d
2
,
d
2
]
for each axis i centered on a user supplied starting voxel wherein d is predetermined based on a vessel category.
5. The method of claim 1 , wherein a voxel {right arrow over (v)} within said 3D cone is defined with respect to said cone tip voxel {right arrow over (x)} as:
v
→
=
Θ
→
(
∇
I
(
x
→
)
∇
I
(
x
→
)
)
[
r
cos
(
θ
)
sin
(
α
)
,
r
sin
(
θ
)
sin
(
α
)
,
r
cos
(
θ
)
]
T
,
wherein
∇
I
(
x
→
)
∇
I
(
x
→
)
is the image gradient,
Θ
→
(
∇
I
(
x
→
)
∇
I
(
x
→
)
)
is a rotation matrix aligning the central axis of the cone with the image gradient direction, α:[−α max ,α max ], is an angle {right arrow over (v)} forms with the image gradient vector,
r
:
[
0
,
h
cos
(
α
max
)
]
,
is a length of a vector represented by {right arrow over (v)}, and θ: [0,2π] is a sweep out parameter.
6. The method of claim 5 , wherein voting comprises incrementing a first vote image as V({right arrow over (x)}+{right arrow over (v)})=V({right arrow over (x)}+{right arrow over (v)})+K({right arrow over (v)})log(|∇I({right arrow over (x)})|) for each voxel {right arrow over (v)} in said cone, wherein
(
-
d
2
,
-
d
2
,
-
d
2
)
≤
(
x
→
+
v
→
)
≤
(
d
2
,
d
2
,
d
2
)
,
incrementing a radius image as R({right arrow over (x)}+{right arrow over (v)})=R({right arrow over (x)}+{right arrow over (v)})+K({right arrow over (v)})|{right arrow over (v)}−{right arrow over (x)}|, wherein K({right arrow over (v)}) is a kernel that weights a point {right arrow over (v)} in the cone by
K
(
v
→
(
r
,
α
)
)
=
exp
(
-
α
2
/
2
σ
2
)
r
cos
(
α
)
.
7. The method of claim 6 , further comprising incrementing a weight image Σ({right arrow over (x)}) for each voxel vote according to Σ({right arrow over (x)}+{right arrow over (v)})=+Σ({right arrow over (x)}+{right arrow over (v)})+K({right arrow over (v)}), and after voting, normalizing said radius image by dividing by said weight image.
8. The method of claim 7 , wherein selecting candidates for a second vote image comprises finding, for each voxel {right arrow over (x)} in the subregion, a voxel, {right arrow over (v)} max within the cone whose tip is attached to voxel {right arrow over (x)} that has the maximum number of votes, and forming said second voting image V′({right arrow over (x)}) by voting for those {right arrow over (v)} max voxels whose distance to said cone tip {right arrow over (x)} is approximately equal to the normalized radius image value for that point, according to V′({right arrow over (x)}+{right arrow over (v)} max )=V′({right arrow over (x)}+{right arrow over (v)} max )+K ({right arrow over (v)} max )log(|∇I({right arrow over (x)})|).
9. The method of claim 8 , further comprising normalizing said second voting image by dividing by said weight image.
10. The method of claim 1 , wherein reconfiguring the subregion as a graph further comprises including two additional nodes, {s,t} representing a source and a sink respectively, linking said source and sink to the voxel nodes pεP forming t-link edges;
t-link weights in the graph are determined by the second vote image, and the p-link weights are determined by the subregion's intensity and gradient.
11. The method of claim 10 , wherein said min-cut algorithm minimizes an energy of the form
E
(
f
)
=
∑
p
∈
P
D
p
(
f
p
)
+
∑
p
,
q
∈
N
V
p
,
q
(
f
p
,
f
q
)
,
wherein f=[0,1] is a label distinguishing a vessel node from a non-vessel node, V p,q is a weight for an n-link connecting voxel nodes p, q,
V
p
,
q
,
∈
N
=
C
1
exp
(
(
(
I
(
p
)
-
I
(
q
)
)
2
+
cos
-
1
(
∇
I
(
p
)
∇
I
(
p
)
·
∇
I
(
q
)
∇
I
(
q
)
)
)
2
σ
2
)
,
D p is a t-link weight defined as D p =C 2 A(p) for voxel nodes p connected to the source node, and D p =C 2 (1−A(p)) for voxel nodes p connected to the sink node, σ 2 , C 1 and C 2 are predetermined constants, and A(p) is a normalized second vote image evaluated for node p.
12. A computer-implemented method for segmenting tubular structures in digital medical images, said method executed by a computer comprising the steps of:
extracting a subregion of a digital medical image volume containing a vessel of interest;
having each voxel in the subregion vote on a vessel center wherein votes are cast in a conical region whose tip is on each said voxel and whose central axis is oriented in a gradient direction, wherein the votes are weighted proportional to the intensity of a gradient at a cone voxel and how close the cone voxel is to the central axis of the cone to create a first vote image, and are weighted proportional to the distance between each voxel in the subregion and the cone voxel to create a radius image;
searching said first vote image for a most popular centerpoint candidate within the conical region for each subregion voxel;
voting for said centerpoint candidate if the distance to the centerpoint candidate from each said subregion voxel is in agreement with the value for a corresponding location in the radius image to form a second vote image;
constructing a graph from the second vote image; and
searching said graph for an optimal path from a user supplied endpoint to any point on a vessel surface to recover said vessel of interest.
13. The method of claim 12 , wherein extracting a subregion comprises:
segmenting a left ventricle from the image volume using an initial point found by searching for bright circular regions at an orientation commonly assumed by the left ventricle;
determining a second point in an aorta based on spatial and intensity priors relative to said left ventricle segmentation; and
extracting an aorta mask, T({right arrow over (x)}), by finding a connected component of voxels with intensities crossing a threshold, wherein the threshold is computed using said second point inside said aorta.
14. The method of claim 13 , further comprising incrementing a weight image for each voxel vote by a weighting factor associated with said cone defined by
K
(
v
→
(
r
,
α
)
)
=
exp
(
-
α
2
/
2
σ
2
)
r
cos
(
α
)
,
wherein voxel {right arrow over (v)} is within said cone, α:[−α max ,α max ] is an angle {right arrow over (v)} forms with the image gradient vector,
r
:
[
0
,
h
cos
(
α
ma
x
)
]
is a length of a vector represented by {right arrow over (v)}, and σ 2 is a predetermined constant, and normalizing said second voting image by dividing by said weight image to form an axial symmetry image.
15. The method of claim 14 , wherein vertices V={B,H,U} of said graph are defined as B={{right arrow over (x)}:∇T({right arrow over (x)})>ε 2 } wherein T({right arrow over (x)}) is the aorta mask, H={{right arrow over (x)}: A({right arrow over (x)})>ε 3 }, the axial symmetry voxels over a threshold; and U is a user supplied endpoint, edges E ij are created between any two vertices i, j which are less than a predetermined distance apart, and the edges are weighted by line integral between two vertices through a multiplicative inverse of the axial symmetry image,
w
(
E
ij
)
=
∫
C
ij
1
A
′
(
s
)
ⅆ
s
,
wherein C ij is a line segment connecting vertices v i and v j .
16. The method of claim 12 , wherein said graph is searched using Dijkstra's algorithm.
17. A program storage device readable by a computer, tangibly embodying a program of instructions executable by the computer to perform the method steps for segmenting tubular structures in digital medical images, said method comprising the steps of:
extracting a subregion from a 3-dimensional (3D) digital medical volume containing a vessel of interest;
for each voxel in the subregion, identifying potential vessel centerpoints by attaching to each voxel a tip of a 3D cone that is oriented in the direction of the voxel's image gradient and having each voxel within the 3D cone vote for those voxels most likely to belong to a vessel centerline;
selecting candidates for a second vote image that are both popular according to a first vote image, as well as being consistently voted upon by a radius image;
reconfiguring the subregion as a graph wherein each voxel is represented by a node that is connected to 26 nearest neighbors by n-link edges; and
applying a min-cut algorithm to segment the vessel within the subregion.
18. The computer readable program storage device of claim 17 , wherein the segmented vessel includes branches.
19. The computer readable program storage device of claim 17 , the method further comprising detecting intersections of the segmented vessel with the subregion's boundaries, wherein intersection voxels serve as starting voxels for extracting a new subregion overlapping with the current subregion, wherein said new subregion serves to track the vessel's centerline.
20. The computer readable program storage device of claim 17 , wherein said subregion is a cubic subregion of size
x
i
=
[
-
d
2
,
d
2
]
for each axis i centered on a user supplied starting voxel wherein d is predetermined based on a vessel category.
21. The computer readable program storage device of claim 17 , wherein a voxel {right arrow over (v)} within said 3D cone is defined with respect to said cone tip voxel {right arrow over (x)} as:
v
→
=
Θ
→
(
∇
I
(
x
→
)
∇
I
(
x
→
)
)
[
r
cos
(
θ
)
sin
(
α
)
,
r
sin
(
θ
)
sin
(
α
)
,
r
cos
(
θ
)
]
T
,
wherein
∇
I
(
x
→
)
∇
I
(
x
→
)
is the image gradient,
Θ
→
(
∇
I
(
x
→
)
∇
I
(
x
→
)
)
is a rotation matrix aligning the central axis of the cone with the image gradient direction, α:[−α max ,α max ], is an angle {right arrow over (v)} forms with the image gradient vector,
r
:
[
0
,
h
cos
(
α
ma
x
)
]
,
is a length of a vector represented by {right arrow over (v)}, and θ: [0,2π] is a sweep out parameter.
22. The computer readable program storage device of claim 21 , wherein voting comprises incrementing a first vote image as V({right arrow over (x)}+{right arrow over (v)})=V({right arrow over (x)}+{right arrow over (v)})+K({right arrow over (v)})log(|∇I({right arrow over (x)})|) for each voxel {right arrow over (v)} in said cone, wherein
(
-
d
2
,
-
d
2
,
-
d
2
)
≤
(
x
→
+
v
→
)
≤
(
d
2
,
d
2
,
d
2
)
,
incrementing a radius image as R({right arrow over (x)}+{right arrow over (v)})=R({right arrow over (x)}+{right arrow over (v)})+K({right arrow over (v)})|{right arrow over (v)}−{right arrow over (x)}|, wherein K({right arrow over (v)}) is a kernel that weights a point {right arrow over (v)} in the cone by
K
(
v
→
(
r
,
α
)
)
=
exp
(
-
α
2
/
2
σ
2
)
r
cos
(
α
)
.
23. The computer readable program storage device of claim 22 , the method further comprising incrementing a weight image Σ({right arrow over (x)}) for each voxel vote according to Σ({right arrow over (x)}+{right arrow over (v)})=Σ({right arrow over (x)}+{right arrow over (v)})+K({right arrow over (v)}), and after voting, normalizing said radius image by dividing by said weight image.
24. The computer readable program storage device of claim 23 , wherein selecting candidates for a second vote image comprises finding, for each voxel {right arrow over (x)} in the subregion, a voxel, {right arrow over (v)} max within the cone whose tip is attached to voxel {right arrow over (x)} that has the maximum number of votes, and forming said second voting image V′({right arrow over (x)}) by voting for those {right arrow over (v)} max voxels whose distance to said cone tip {right arrow over (x)} is approximately equal to the normalized radius image value for that point, according to V′({right arrow over (x)}+{right arrow over (v)} max )=V′({right arrow over (x)}+{right arrow over (v)} max )+K({right arrow over (v)} max )log(|∇I({right arrow over (x)})|).
25. The computer readable program storage device of claim 24 , the method further comprising normalizing said second voting image by dividing by said weight image.
26. The computer readable program storage device of claim 17 , wherein reconfiguring the subregion as a graph further comprises including two additional nodes, {s,t} representing a source and a sink respectively, linking said source and sink to the voxel nodes pεP forming t-link edges;
t-link weights in the graph are determined by the second vote image, and the n-link weights are determined by the subregion's intensity and gradient.
27. The computer readable program storage device of claim 26 , wherein said min-cut algorithm minimizes an energy of the form
E
(
f
)
=
∑
p
∈
P
D
p
(
f
p
)
+
∑
p
,
q
∈
N
V
p
,
q
(
f
p
,
f
q
)
,
wherein f=[0,1] is a label distinguishing a vessel node from a non-vessel node, V p,q is a weight for an n-link connecting voxel nodes p, q,
V
p
,
q
,
∈
N
=
C
1
exp
(
(
(
I
(
p
)
-
I
(
q
)
)
2
+
cos
-
1
(
∇
I
(
p
)
∇
I
(
p
)
·
∇
I
(
q
)
∇
I
(
q
)
)
)
2
σ
2
)
,
D p is a t-link weight defined as D p =C 2 A(p) for voxel nodes p connected to the source node, and D p =C 2 (1−A(p)) for voxel nodes p connected to the sink node, σ 2 , C 1 and C 2 are predetermined constants, and A(p) is a normalized second vote image evaluated for node p.