Vascular Structure Identification in Intraoperative 3D Contrast-Enhanced Ultrasound Data

In this paper, a method of vascular structure identification in intraoperative 3D Contrast-Enhanced Ultrasound (CEUS) data is presented. Ultrasound imaging is commonly used in brain tumor surgery to investigate in real time the current status of cerebral structures. The use of an ultrasound contrast agent enables to highlight tumor tissue, but also surrounding blood vessels. However, these structures can be used as landmarks to estimate and correct the brain shift. This work proposes an alternative method for extracting small vascular segments close to the tumor as landmark. The patient image dataset involved in brain tumor operations includes preoperative contrast T1MR (cT1MR) data and 3D intraoperative contrast enhanced ultrasound data acquired before (3D-iCEUSstart) and after (3D-iCEUSend) tumor resection. Based on rigid registration techniques, a preselected vascular segment in cT1MR is searched in 3D-iCEUSstart and 3D-iCEUSend data. The method was validated by using three similarity measures (Normalized Gradient Field, Normalized Mutual Information and Normalized Cross Correlation). Tests were performed on data obtained from ten patients overcoming a brain tumor operation and it succeeded in nine cases. Despite the small size of the vascular structures, the artifacts in the ultrasound images and the brain tissue deformations, blood vessels were successfully identified.


Introduction
Intraoperative ultrasound imaging is nowadays commonly used in neurosurgery during brain tumor operations [1]. At the beginning of the intervention, the ultrasound images show the surgeon the intraoperative state of the tumor (Figure 1b) [2]. The tumor size or position can be possibly slightly different at the operation time point from the preoperative state depicted in the preoperative MR data (Figure 1a). During the operation ultrasound imaging is a valuable tool to detect the residuals of tumor with the goal to optimize the tumor removal (Figure 1c,d) [3]. However, the interpretation of the ultrasound can be complex [4,5]. The acquisition of the images through the skull opening, called craniotomy, requires skill and experience. Firstly, the sweep of the ultrasound probe is limited by the small opening. The presence of air between the probe and the brain surface blocks the propagation of ultrasound burst. Secondly, the scanned brain area in the images is limited. The orientation of the image and the interpretation of the information it shows can be complex, mostly for images acquired at the end of the operation. Figure 1. Patient image data acquired during tumor resection. The brain tumor (white arrow) is represented in the preoperative cT1MR data (a) and in the 3D-iCEUS start (b) acquired at the beginning of the operation. After tumor removing, the resection cavity (indicated by two white arrows in (c) and (d) is well visible in the B-mode ultrasound image (c). In the 3D-iCEUS end acquired at the end of the operation (d); the borders of the cavity can be easily interpreted as a blood vessel.
A common method to assist the neurosurgeon in the visual analysis of the intraoperative ultrasound images is to visualize them overlapped to the preoperative MR data in the navigation system because the preoperative image data offers an overview of the entire head [2,6,7]. This technique has limitations for ultrasound images acquired at the end of the operation [5]. The brain structures are so largely deformed, mainly because of tissue resection, so that the information included in the visualization planes does not correspond at all anymore [8].
The evaluation and correction of the brain shift has been extensively studied. Methods are based on similarity measures based on image intensities [9][10][11][12], anatomical landmarks [13][14][15] or biomedical deformation models [16,17]. An interesting previous work proposed to use the vascular structures surrounding the tumor as landmarks in order to estimate the brain deformations [14,18]. In this work, the extraction of the blood vessels was performed in preoperative MR angiographic (MRA) data and in intraoperative Doppler ultrasound images. MRA examination is in routine not acquired in case of brain tumor. Doppler ultrasound can be performed intraoperatively but has limitations to depict the diameters with accuracy and to provide 3D visualization. Moreover, information about brain tissue itself is missing. Several medical studies were conducted to evaluate contrast-enhanced ultrasound imaging for the visualization of brain tumors and the detection of residuals of tumor [19][20][21]. This imaging modality enables to enhance vascularized structures, like lesions and tumors, but also blood vessels [22,23].
In this paper, an alternative method of vascular structures identification in contrast-enhanced ultrasound data is presented with a possible application in the brain shift estimation. Instead of considering the preoperative MRA and intraoperative Doppler ultrasound as previously mentioned, this work introduces the possibility of using cT1MR and 3D-iCEUS data which are routinely involved in brain tumor operations. The method has three main steps i.e., the selection of a vascular segment in cT1MR, the identification of this pattern firstly in 3D-iCEUS start , then in the 3D-iCEUS end .
Vesselness-based segmentation methods and rigid registrations techniques were employed, ensuring a computing time compatible with an intraoperative use. Another contribution of this work is the evaluation of the method on CEUS images acquired at the end of the operation, when tissue deformations are large.

Patient Image Dataset
Tumor operations were guided by using a neuronavigation system (SonoNavigator, Localite, Sankt Augustin, Germany), including an AplioXG ultrasound device (Toshiba Medical Systems Europe, Zoetermeer, Netherlands) with 2D ultrasound transducers. At the beginning of the intervention, preoperative 3D contrast T1 Magnetic Resonance (cT1MR) data were registered with the patient's head, based on anatomical landmarks. The procedure was then improved using a surface-based registration technique. During the operation and in addition to the acquisition of intraoperative B-mode ultrasound volumes (noted 3D-iUS), two 3D intraoperative contrast enhanced ultrasound (3D-iCEUS) data were acquired after the injection of 4.8 milliliter of an ultrasound contrast agent (SonoVue, Bracco s.p.a., Milan, Italy) at a rate of 3.0 mL/min using a syringe pump (ACIST VueJect, Bracco s.p.a, Milano, Italy) and using the contrast harmonic imaging (CHI) method. The size and density of the micro-bubbles in the SonoVue contrast agent was 1.9 ± 0.1 micromillimeter and 3.4 ± 0.5 micro-bubbles/mL as recently reported in [24]. The first volume, denoted 3D-iCEUS start , was obtained transdurally immediately after the craniotomy when the tumor was still entire, and the second one, denoted 3D-iCEUS end , was acquired after removing the tumor at the end of the surgery. The surgeon scanned the cerebral region of interest (ROI) with the ultrasound probe (large linear array transducer with a field of view of 38 mm, a range of frequencies of 4.8 MHz to 11.0 MHz and an average frequency of 8 MHz was selected for the acquisitions). The transducer was tracked. The neuronavigation system reconstructs a 3D volume from the 2D-iCEUS slices and overlaid it onto the cT1MR data. The pixel size in the original 2D ultrasound images is 0.422 mm × 0.422 mm and the voxel size of the reconstructed 3D-iCEUS data is 1 × 1 × 1 mm 3 . The data involved in the surgery process are shown in Figure 1.

Vascular Structure Segmentation
A common solution to distinguish the blood vessels in the 3D-iCEUS data is their extraction by using segmentation methods. However, the segmentation of vascular structures in ultrasound image data is generally a complex task because of the speckle and the blood vessel diameters less than 2 millimeters (Figure 1b,d). Model-based techniques, which include a priori knowledge of the object shape or its data representation, are required to guide the segmentation process and improve the success of algorithms [25,26]. In this work, the extraction of blood vessels is performed in two steps: the computation of a Hessian based vesselness measure for enhancing the vascular structures, followed by the Otsu thresholding algorithm [27]. Moreover, the segmentation is performed in limited target regions of interest in the 3D-iCEUS data, as it will be explained in Section 2.3. The vesselness measure, which is used, was introduced by Sato [28,29] to describe tubular structures within an image using the Hessian matrix. The Hessian matrix H of an image I is defined by: The second order derivative of the image I is computed from convolution with a Gaussian kernel of standard deviation σ.
where G σ is a Gaussian function with a standard deviation σ.
Moreover, it has been shown that the variation of sigma values σ enables to describe vascular structures with different radii sizes. Therefore, v(σ) is computed at each voxel position in the image for different values of the standard deviations σ min ≤ σ ≤ σ max and the maximum response kept v I is described by: Finally the obtained image of vesselness responses is thresholded using the Otsu method [27] in order to extract the vascular structures. Figure 2 presents an illustration of blood vessel segmentation process from cT1MR, 3D-iCEUS start and 3D-iCEUS end . An additional step of 3D representation is added here.

Vascular Structure Identification
The proposed method consists in identifying the vascular structures in the image data involved in the patient treatment. During the planning stage, some vascular structures close to the tumor are selected by the neurosurgeon. These patterns are automatically recognized in the 3D-iCEUS data by using image registration techniques. The identification method of vascular structures is depicted in Figure 3 and is described as follows:

•
Step 1: Selecting a vascular segment pattern in cT1MR Interaction with the application in the operating room has to be limited because of sterilization constraints and restricted time. Due to this, during the operation planning, the user delineates interactively a region of interest including a blood vessel near to the tumor in the cT1MR data. The vascular structure is segmented using the method described in Section 2.2. It performs well because the blood vessels are enhanced in the cT1MR data due to the contrast agent, and any other anatomical structure represented with similar intensities is included in the region of interest.
The segmented blood vessel represents the pattern (white frame, Step 1, Figure 3) that is searched for in the 3D-iCEUS start and 3D-iCEUS end data.

•
Step 2: Blood vessel identification in 3D-iCEUS start The pattern is firstly searched in the 3D-iCEUS start data acquired before resection. In order to reduce the computing time, the search space (large white frame, Step 2, Figure 3) is smaller than the entire 3D-iCEUS start , but large enough to take the tissue deformations into account. It has twice the volume of the region of interest defined in the cT1MR data, and is centered on the same image position. The enhanced structures in this region are segmented and then a rigid registration method is used to find the sample in the 3D-iCEUS start data (yellow frame, Step 2, Figure 3) which corresponds best to the pattern. A rigid transformation is sufficient here, since the goal is the identification of the position of the blood vessel in the 3D-iCEUS start data, which looks like the vascular segment pattern selected in the cT1MR data. The blood vessel detected in the 3D-iCEUS start then becomes the new pattern, which has to be identified in the 3D-iCEUS end data after resection.

•
Step 3: Blood vessel identification in the 3D-iCEUS end A similar method as described in the step 2 enables to find the position of the target blood vessel (red frame, Step 3, Figure 4) within the search space (yellow frame, Step 3, Figure 3) defined in the 3D-iCEUS end data acquired after tumor resection.
Finally, the preselected vessel in the cT1MR is overlapped on the 3D-iCEUS start and 3D-iCEUS end images for visualization purpose. Figure 4 illustrates the workflow of the entire method and illustrations of the blood vessel identification on patient data is shown in Section 3.

Blood vessel identification Identification
Tissue deformation Step 2 Step 3 Pattern After tumor resection Before tumor resection The day before the surgery It should be noted that the method is so far not able to recognize an incorrect identification of the selected blood vessel. In this case, the user has to reposition interactively the area including the vascular segment corresponding to the pattern. However, the approach is to minimize the interaction between surgeon and computer in operating room because of the sterilization requirement.

Validation
One key point in image registration is the similarity measure. The Normalized Cross Correlation (NCC) and the Normalized Mutual Information (NMI) are typical similarity measures commonly used [32,33]. The NCC computes pixel-wise normalized cross-correlation between two images and usually aims at registering images of the same modality. The NMI describes how well one image can predict the other one and is, therefore, suitable for multi-modality registration. Robust and easier to interpret, the Normalized Gradient Field (NGF) measures a normalized distance between the gradients of the images to be registered [34,35]. Some works proposed calculating the similarity measures from binary data, arguing that the registration process is, in this way, accelerated [36][37][38][39]. These three similarity measures were used for comparison in the registration algorithm. In order to validate the method, a neurosurgeon visually checked that the vascular segments, automatically identified in the 3D-iCEUS start and 3D-iCEUS end , corresponded to the blood vessel selected in the cT1MR data. Moreover, the results of the algorithm were quantitatively compared with a registration, interactively performed, using the Dice Similarity Index (DSI) and the Hausdorff distance. The DSI represents the percentage of overlapping of two objects. Its value is 1 if they overlap perfectly and 0 if they do not intersect at all. The Hausdorff distance represents the maximal distance between two objects. The computing time of the registration process was also measured for each similarity measure. The NCC, NMI and NCC measures were alternately computed by using the following equations: where x i and y i are respectively the set of points of images X (fixed image) and Y (moving image). In addition,x andȳ are the mean of X and Y; and σ x and σ y are their standard deviation. H(X) and H(Y) are the entropies of random variables X and Y, and H(X, Y) represents their joint entropy.
n(I, x) := where n(I, x) is a regularized normalized gradient field of a given image I(x). Then, d c (Y, X) and d d (Y, X) represent the distance between the gradients of images X and Y. Besides, the Dice Similarity Index (DSI) and the Hausdorff distance are described as follows: where X is the fixed image and Y the moving image.

Results
The collection of data has been performed at the University Hospital, Department of Neurosurgery, University of Leipzig, Germany, in the context of a previous research project funded by the German Research Society (Deutsche Forschungsgemein-schaft) and accepted by the ethic commission of the University of Leipzig. The implementation was done with an Intel Celeron, 1.5 GHz and 2 GB of memory using MeVisLab tool. The manual registration was validated by the neurosurgeon of this institution. The identification algorithm was tested on ten datasets of patients who overcame a brain tumor operation. A simple vascular segment near to the tumor was selected in the preoperative cT1MR data. Table 1 includes the DSI, Hausdorff distances and the processing time values computed between the segmented blood vessels in the cT1MR and 3D-iCEUS start data and between the segmented blood vessels in the 3D-iCEUS start and 3D-iCEUS end data, both normalized with the scores obtained by the expert registration. The values followed by a star indicate that the algorithm failed to find the correct vascular segment. If it occurred before resection, the identification in the 3D-iCEUS end data was not performed. It is clearly observed that the proposed approach was able to find correctly the targeted vascular structure with at least one similarity measure used.
Besides, Table 2 includes the mean DSI values and the mean Hausdorff distances, averaged on the successful cases, while Table 3 presents the mean computing time. The outcomes show that the NGF achieves the highest mean rate of registration than the NMI and NCC. Nevertheless, the last ones performs the registration in a shortest time.
In Figure 5, the results obtained step by step for five cases are illustrated. First, the ROI is defined in the cT1MR by the surgeon (column (a)) in the planning stage and it encloses the selected vascular structure so as to reduce the space of segmentation. Second, the vessel is extracted by using the process described in Section 2.2 (column (b)). Third, the vascular structure extracted from the cT1MR data is overlaid on the 3D-iCEUS start (column (c)). It should be observed that the extracted structure is not spatially aligned with its corresponding in the 3D-iCEUS start . Fourth, after the segmentation of the vascular structure in 3D-iCEUS start , column (d) presents the registration results, in which, vascular structures are aligned. At this step, the blood vessel found in 3D-iCEUS start becomes the new pattern and it is superimposed on the 3D-iCEUS end image in column (e). In most of cases, there is no matching between the new pattern with its corresponding vessel in the 3D-iCEUS end because of large tissue deformations at this stage of the operation. Finally, the registration of vascular structures from the 3D-iCEUS star to 3D-iCEUS end is carried out (column (f)).  Figure 5. Illustration of vascular structures identification for five cases: (a) preoperative cT1MR; (b) vascular structure selection in cT1MR; (c) overlapping of vascular structure segmented in cT1MR on CEUS start before registration; (d) overlapping of vascular structure segmented in cT1MR on CEUS start after registration; (e) overlapping of vascular structure segmented in CEUS start on CEUS end before registration; (f) overlapping of vascular structure segmented in CEUS start on CEUS end after registration.  Table 3. Mean processing time in s during the registration computed for the successful cases (patient 10 excluded).
Similarity Measure cT1MR − 3D-iCEUS start 3D-iCEUS start − 3D-iCEUS end An illustration of using vascular structures from cT1MR and 3D-iCEUS start for correcting the brain shift is depicted in Figure 6. First, the blood vessels are extracted from the both modalities and by overlapping them it is clearly observed that they are not aligned. Second, a registration process is performed to align the segmented structures. By registering the targets, a matrix denoted T which describes the spatial transformation necessary for aligning structures is obtained. Then, this transformation matrix is applied on the cT1MR images so as to achieve its matching with the 3D-iCEUS start data. Finally, the superimposing of these image modalities is presented for the input before registration and for the output after registration. Figure 6. Application of vascular structure identification in brain shift correction between cT1MR and 3D-iCEUS start .

Visual Validation
The visual checking for identification a selected blood vessel was performed without problem in the 3D-iCEUS start . The ultrasound images were easy to interpret because the anatomical structures are, in general, similarly represented in the 3D-iCEUS start and in the preoperative cT1MR data. Brain tissues were already moved at this stage of the operation (after the craniotomy) by comparing the images. The corresponding planes in the 3D-iCEUS start and preoperative cT1MR data were coarsely only translated and the visual comparison of anatomical structures in these planes were still possible. After tumor resection, the comparison of the 3D-iCEUS start and the 3D-iCEUS end was much more complex for the expert. Firstly, the loss of tumor in the 3D-iCEUS end deprives of a valuable reference structure for the comparison of the images. Secondly, the image quality was by trend lower in the 3D-iCEUS end . The acquisition was performed with the resection cavity filled with liquid in order to conduct the ultrasound waves. The liquid outflow outside the cavity caused artifacts in the images. Thirdly, brain tissues deformed largely between the beginning of the operation and the stage after tumor resection. Corresponding anatomical structures were elastically deformed. They were represented in different slices and the mental reconstruction of the information was difficult. Therefore, the postoperative MR image were important for assisting the visual evaluation.

Quantitative Validation
Additionally, the Dice similarity index and the Hausdorff distance were used for quantitatively assess the results obtained from the experiments. It was observed that the NGF similarity provided the highest mean DSI values and the lowest mean Hausdorff distances (Table 2). In contrast, it was found as the slowest (Table 2). Compared to the NMI and NCC, the particularity of the NGF measure is its robustness in cases of curved vascular segments. A T-test was performed to compare the DSI values obtained with the different similarity measures and it showed that the NGF and NCC are statistically different (p-value: 0.0304). As well as, the p-values obtained by comparing the NGF and NMI, then the NMI and NCC were 0.0611 and 0.1747, respectively. In general, the computing time is lower than 15 s on average. The calculation of the NCC and NMI involves the image intensities and was more quickly performed (computing time lower than 10 s on average, Table 3) while the NGF requires the previous computation of the image gradients (computing time lower than 13 s on average, Table 3). However, these mean computing time values are still acceptable for use in the operating theater. The size of the ROI plays an important role on the value of the computing time. It should be noticed that the results from the Table 1 show that the method will always perform a registration even if the correct vascular structure was not found. Despite the tissues deformation which tends to decrease the DSI value, a higher value than 0.60 was found acceptable compared to the neurosurgeon validation. However, the DSI close to 1 indicate the perfect matching obtained.

Image Quality
The quality of the ultrasound images limits the performance of the registration algorithm. Mainly, it is affected by factors such as: the size of the craniotomy, the location and depth of the tumor in the brain that may complicate the acquisition process. Also, a precise sweeping of the 2D US probe has to be performed within a restricted time window, when the contrast agent enhances optimally the structures. The time window for an optimal acquisition is about 30 s [21]. Nevertheless, a technical solution to improve the image quality is the reduction of the size of voxel in the 3D ultrasound volume, for example 0.5 mm × 0.5 mm × 0.5 mm. The use a 3D ultrasound transducer will enable to perform the acquisition faster to overcome the washed in and out of contrast agent [40].

Misidentification of Vascular Structures
The algorithm succeeded better if a vessel with a specific shape (e.g., curve) is selected. A single straight vascular segment can be in certain cases, confused with other vascular segments or with the resection cavity edges in the 3D-iCEUS end data if they are close. Possible approaches to solve this problem would be to segment no-vascular structures in the images in order remove them in the registration process. For example, the tumor can be extracted using a model-based segmentation technique [41]. The resection cavity can be easier identified in the B-mode ultrasound images. Another interesting possibility to test in the future could be the use of the Scale-Invariant Feature Transform (SIFT) proposed by David Lowe [42,43] for registering images as applied in [44,45]. However, in spite of the difficulties mentioned above, the proposed methodology was able to find the correct vascular structures via at least one similarity measure used and, furthermore it was capable to carry out successfully the registration.

Vascular Structures Segmentation
The traditional vesselness based on Hessian matrix was used in this study for vascular structure segmentation. The future intent could be, for instance, the use of alternative methods based on Gabor filter and multiobjective optimization applied in [46] on automatic segmentation of coronary arteries.

Conclusions
Despite the small size of the vascular structures surrounding the brain tumors, the low signal to noise ratio in ultrasound image and the brain tissue deformation, it was possible to correctly identify vascular segments in 3D-iCEUS patient's data. Moreover, we showed that the NGF is more robust, especially in cases of vascular segments with specific shapes. However, the NCC and NMI have a lower computation time than the former. But it is important to notice that the computation time achieved by the algorithm by using these three similarity measures is still compatible with intraoperative use. Besides, an application of this work is the use of the vascular structures as landmarks for the estimation and correction of brain shift. The future works should be the improving of the vascular segmentation in 3D-iCEUS due to the low signal noise ratio in this modality in general. On the other hand, a solution has to be found in cases where the blood vessel is not visible in the intraoperative ultrasound data. For instance, the combination of information included in the B-mode and contrast-enhanced ultrasound data should make the method more robust.