A patient-specific geometrical model for coronary arteries: introducing Yağmur

In order to perform realistic numerical simulations of the blood flow problem or the fluid-structure interaction problem between the blood and the arterial walls, the geometrical representation of the vascular network has to be as detailed as possible. Within this framework, a precise patient-specific reconstruction of the arterial tree is compulsory and in this blog post we present a possible approach to this problem.

Yağmur is a Master student at Istanbul Technical University and my supervisor for the Summer of HPC program at the Turkish Center for High Performance Computing. Moreover, she is the kind owner of the heart I am presenting for my analysis thus she will be our patient for today :-).

Data assimilation: CT scan

CT stands for Computed Tomography and is an X-rays based medical imaging procedure. This data assimilation method produces a set of slices of the area under analysis along three axis and then digitally combines them in order to generate a 3D image of the inside of the object. The cross-sectional outputs of the CT scan can be used both for clinical purposes and to reproduce a 3D tomographic representation of the structures as in figure 1.

Inverse problem in imaging

Starting from the information contained in the images produced by a CT scan, we aim at retrieving a full description of the geometry of Yağmur’s Right Coronary Artery (RCA from now on). This is an example of a so-called inverse problem, that is we seek the best model (in our case, the geometry) that approximates a phenomenon or a situation (Yağmur’s RCA) given some known observations (the images arising from the CT scan). In general, this kind of problems is very difficult to solve and its numerical counterpart is usually ill-posed: this results in unstable numerical schemes whose solution is very sensitive to perturbations of the data. Therefore an alternative strategy is required.

An evolutionary strategy for image segmentation

In order to reconstruct the 3D geometrical model of the RCA, first we have to identify its contours within the sets of images provided by the CT scan. These files are in DICOM format (Digital Imaging and COmmunications in Medicine) which is a standard in medical image analysis: DICOM files are grayscale images within which the intensity in the voxels – meaning volumetric pixels – is associated with the corresponding Hounsfield Unit, that is a quantitative information of the attenuation coefficient of the material.

To bypass the issues related to the consistent solution of an inverse problem, a classical approach in the literature relies on the parameterization of the contours of the objects to be identified by means of a level-set curve. Thus the original imaging problem consisting of the identification of the most suitable geometry to describe the RCA is replaced by an algorithm that analyzes the evolution of its contours driven by a partial differential equation well known in the literature [1].

The Vascular Modeling Toolkit

VMTK, that is the Vascular Modeling Toolkit, is an open-source utility that provides several features in the field of medical image segmentation and data processing of blood vessels information. This toolkit is based on VTK and ITK libraries and can be easily installed for Windows, Mac OSX and Ubuntu users starting from the precompiled binary packages available on the official website. Since I am a Fedora user, I had to compile it from source: the overall installation may take some time if you do not have the VTK libraries already installed on your system but the procedure is straightforward even if you are not an expert user of the Linux terminal.

Balloon, fast marching algorithm and colliding fronts

VMTK provides several algorithms for image segmentation. The Balloon segmentation is an extension of the level-set algorithm to the 3D case and relies on the idea of inserting an empty balloon inside the object to be segmented and starting to inflate it until it achieves all the boundary walls. However, VMTK offers a smarter option that assists the segmentation when the target is an arterial vessel: the so-called colliding fronts feature allows to restrict the segmentation to the vessels where two wave fronts travelling in opposite directions meet, neglecting smaller branches (see figure 2, courtesy of Luca Antiga [2]).

Thus my first task relied on the identification of the RCA location within the DICOM images. Studying some anatomy of the heart, I figured out that the RCA can be easily recognized within a CT scan or a Magnetic Resonance Imaging output when observing an horizontal section of the chest from above. As we can observe in figure 3, the image is not symmetrical and we can use this information to orientate ourselves: the heart is located in center-left part of the scan, so the left lung is the smallest of the two black regions we are observing. Scrolling among the different sections acquired by the scan we seek the image where the big light gray circle in the center is deforming towards an arch shape: this is the starting point of the RCA and the location to put the seed for segmentation in.

A critical step in the segmentation of the RCA is the positioning of the seeds for the colliding fronts sources. As a matter of fact, if the seed is too close to the aorta – whose dimension is considerably higher than the one of coronary arteries – the segmentation algorithm will tend to segment the bigger vessels too. A simple way to avoid this drawback relies on putting the seed already inside the coronary artery even if this results in the loss of the top section of the vessel. The resulting segmentation is presented in figure 4.

Surface smoothing and post-processing

The outcome of the segmentation procedure may present irregular or noisy regions thus a smoothing procedure is needed. However the user has to be particularly careful while performing this operation in order to increase the regularity of the surface preserving as many details as possible. Eventually to transform the surface into a volumetric representation of a network of arteries, the ending parts of the vessels have to be cut and open as displayed in figure 5.

Mesh generation

The final step we have to perform before being able to simulate our blood flow relies on the generation of a computational mesh that allows us to approximate the problem by means of a Finite Element or Finite Volume based tool. VMTK helps us in this task too, providing support to export the mesh in several different formats. For example, we recall the .msh format compatible with the commercial software FLUENT and the .lifev one used by the LifeV library developed by Ecole Polytechnique Fédérale de Lausanne, Politecnico di Milano, INRIA Rocquencourt and Emory University.

During the mesh generation procedure, it is important to precisely describe the regions near the walls both for fluid and fluid-structure interaction problems: in particular, boundary layers develop in these areas since no-slip conditions are imposed on the walls and the computation of wall shear stresses is subject to a detail description of the structural deformation of the arteries. Here is the final mesh exported in a FLUENT compatible format and visualized in ParaView.

References

[1] T.F. Chan and L.A. Vese. Active contours without edges. IEEE Transactions on Image Processing, 10(2):266–277, 2001.

[2] L. Antiga. Patient-Specific Modeling of Geometry and Blood Flow in Large Arteries. Ph.D. thesis in BioEngineering at Politecnico di Milano, 2002.