Centerline detection of (cardiac) vessels in CT images Martin Korevaar

1. Centerline detection of (cardiac) vessels in CT imagesMartin Korevaar Supervisors: Shengjun Wang Han van Triest Yan Kang Bart ter Haar Romenij

2. Overview • Introduction • Method • Feature space • Minimal Cost Path (MCP) search • Results • Conclusion • Discussion

3. Introduction • Coronary artery and heart diseases are one of the main causes of death in Western world • Therefore improvement of diagnosis, prevention and treatment is needed • Diagnosis is improved by techniques like CTA and MRA • CAD is needed to analyze huge amount of data

4. Introduction • Centerline of vessel is interesting feature for CAD • Needed e.g. for CPR and lumen size measurements • Should be • Independent of vessel segmentation • Robust with respect to image degradations.

5. CPR • Visualization method Images from: Armin Kanitsar in IEEE Viz. 2002 Okt

6. CPR • Visualization method • Maps 3D path on 2D image Images adapted from: Armin Kanitsar in IEEE Viz. 2002 Okt

7. CPR • Visualization method • Maps 3D path on 2D image • Artery can be investigated from just 1 image

8. CPR • Visualization method • Maps 3D path on 2D image • Artery can be investigated from just 1 image • Centerline needs to be correct Images from: Armin Kanitsar in IEEE Viz. 2002 Okt

9. Algorithm Create feature space Feature space with: • Center of vessel lower value then peripheryof vessel • Vessel lower value then background Find Minimum Cost Path with Dijkstra’s algorithm

10. Feature space • Filter based on eigenvalues of Hessian

11. Dxx*G(x,y,) G(x,y,) Dyy*G(x,y,) Dxy*G(x,y,) Hessian (1) • Matrix with all second order derivatives2D Hessian: • Derivative of image:Convolve derivative of Gaussian with image • Gaussian:

12. 1v1 2v2 Hessian (2)Eigenvectors / -values • 1< 2(<3) • measure of curvature • i relates to vi • v2 points to max. curvature • v1 points to min. curvature, perpendicular to v1 • Rate of change of intensity is curvature of an image

13. Hessian (2)Eigenvectors / -values 3v3 Eigenvalue of Hessian of a pixel gives information about the local structure (e.g. tube) 2v2 1v1

14. Tubular structure filter Frangi Distinguishes blob-like structures, cannot distinguish between plate and line-like structures Distinguishes between line and plate-like structures Filters noise.

15. Tubular structure filters Frangi Wink et al. α = β = 0.5 and c = 0.25 Max[Greyvalue] Olabarriagaet al. α = 1, β = 0.1 and c > 100 => better discrimination center / periphery vessel Chapman et al. α = 0.5, β = ∞ and c > 0.25 Max[Greyvalue] => drops RB-term => better discrimination vessel / background

16. Tubular structure filters HessDiff Better discrimination background / vessel

17. Tubular structure filters • Hessian is calculated at multiscale (2.6    18.6) • Scale with highest response is scale of a voxel • Highest response is in center vessel • Hessian is scaled with 2

18. Algorithm Create feature space Feature space with: • Center of vessel lower value then peripheryof vessel • Vessel lower value then background Find Minimum Cost Path with Dijkstra’s algorithm

19. Invert pixel values(1/pixelvalue) • Response highest at center of the vessel. • Minimal cost path needs lowest • Pixel values are inverted 1 / I

20. Minimal cost path • Select start and end point • Find minimal cost path in between • Dijkstra’s algorithm to find that minimal cost path

21. Dijkstra 2a 4 5 5 2 End 7 6 3 7 7 6 1 Begin 7 5 6 10 20 3

22. Dijkstra 2b Find min neighbours (green) Add it to investigated nodes (red). Remember predecessor . 4 5 5 2 (1) [6] . . 7 6 3 7 7 6 1 7 5 (1) [7] . 6 8 20 3 (1) [8]

23. Dijkstra 2c Add neighbours (green) Remember predecessor . . 4 4 (2) [11] 4 5 5 2 (1) [6] 2 (1) . 7 6 3 3 7 7 7 6 1 7 5 5 (1) [7] . 6 8 20 3 (1) [8]

24. Dijkstra 2d Find min neighbours (green) Add it to investigated nodes (red). Remember predecessor . . 4 (2) [11] 4 5 5 2 (1) [6] 2 (1) . 7 6 3 3 7 7 7 1 6 7 5 (1) [7] . 6 8 20 3 (1) [8]

25. Dijkstra 2e • Add neighbours (green) • Remember predecessor Update predecessor and cost if node already in neighbours . . . 4 (5) [10] 4 5 5 2 (1) [6] 2 (1) . 7 (5) [13] 7 6 3 3 7 7 7 1 6 7 5 5 (1) 5 (1) [7] . 6 8 20 3 (1) [8]

26. Dijkstra 2f Find min neighbours (green) Add it to investigated nodes (red). Remember predecessor . . . 4 (5) [10] 4 5 5 2 (1) [6] 2 (1) . 7 (5) [13] 7 5 3 3 7 7 7 1 6 7 5 5 (1) 5 (1) [7] . 6 8 20 3 (1) [8]

27. Dijkstra 2g Add neighbours (green) Remember predecessor . . . 4 (5) [10] 4 5 5 2 (1) [6] 2 (1) . 7 (5) [13] 7 6 3 3 7 7 . 7 1 6 7 5 (1) [7] 5 5 (1) . 6 (3) [28] 8 20 3 (1) [8]

28. Dijkstra 2h Find min neighbours (green) Add it to investigated nodes (red). Remember predecessor . . . 4 (5) [10] 4 5 5 2 (1) [6] 2 (1) . 7 (5) [13] 6 3 3 7 7 . 7 1 6 7 5 (1) [7] 5 5 (1) . 6 (3) [28] 8 20 3 3 (1) [8]

29. Dijkstra 2i . Goal! . . 4 (5) [10] 4 (5) 5 5 2 (1) [6] 2 (1) . 7 7 (5) [13] 6 3 3 7 7 . 7 1 6 7 5 (1) [7] 5 5 (1) . 6 (4) [17] 8 20 3 (1) [8] 3

30. Dijkstra 2i . Backtrack: 7 => 5 => 1 . . 4 (5) [10] 4 (5) 5 5 2 (1) [6] 2 (1) . . 7 7 (5) [13] 6 3 3 7 7 . 7 1 6 7 5 5 (1) [7] 5 (1) . 6 (4) [17] 8 20 3 3 (1) [8]

31. Minimal cost path Defined cost: V(i) is the voxel value a is weight factor i is ith voxel of the path

32. Minimal cost path Defined cost: V(i) is the voxel value a is weight factor i is ith voxel of the path Higher values of a => Relative difference between center and surrounding increases => Will follow minimum better instead of shortest path

33. Algorithm Calculate response to (Frangi’s) filter Get Eigenvalues of the Hessian Matrix Invert pixel values (1/pixelvalue) Find minimum cost path with dijkstra’s algorithm

34. Experiments • Original method on different datasets • On worst performing dataset • a = 1 • a = 5 • Different filters • Frangi with different parameters • Wink • Olabarriaga • Chapman • HessDiff

35. Result (CPR) 3 Datasets (1) LAD (2) RCx (3) LAD

36. Different filters and cost functions (CPR) Proximal part Wink a=1 Wink a=5 Chapman a=1 Chapman a=5

37. Different filters and cost functions (CPR) Proximal part Olabarriagaa=5 Olabarriaga a=1

38. Olabarriaga (a=5) Wink (a=5) Different filters: Wink and Olabarriaga Sagittal slice at the stenosis Frangi’s filter with Wink’s constants Frangi’s filter with Olabarriaga‘s constants

39. Different filters: Wink and Olabarriaga Sagittal slice heart wall

40. Different filters and cost functions (CPR) Proximal part HessDiff a=1

41. Different filters: HessDiff Low response at stenosis Lot of false positives Strong false positives at the heart wall Response CT

42. Different filters and cost functions Distal part a=1 a= 5 a=1 a= 5 a=1 a= 5 a=1 a= 5 Chapman Wink Olabarriaga HessDiff

43. Different filters Doesn’t follow vessel at heart wall Wink Grey Olabarriaga Chapman HessDiff

44. Multiscale Wink • Vessel response at low scale • Heart wall response at high scale • Heart wall response is stronger σ= 2.6 σ= 6 σ= 10 σ= 16 σ= All Grey value

45. Different cost function a=1 and a=5 a=1 • α = 5 vs. α =1 • Investigates nodes in smaller area • Less computations • More able to follow local minima • Less able to pass local maxima (stenosis)

46. Discussion • Used scales were high (2.6    18.6) • High responses of the heart wall => bad centerline extraction • HessDiff and Olabarriaga track the centerline badly: • Low response at stenosis. • HessDiff lot of false positive response • Wink and Chapman track the centerline excellent.

47. Conclusion • The centerline is tracked in most cases (more or less) accurate • Wink and Chapman are best filters. • They can even coop with a stenosis. • It returns to the center even if it gets outside the vessel (robust) • Different cost functions yield different results: • High power more precise in details and faster. • Low power more robust and slower.

48. Further research • Smaller scales might improve results • Use Wink’s constants

49. Questions?