Position tracking of moving liver lesion based on real-time registration between 2D ultrasound and 3D preoperative images.

Position tracking of moving liver lesion based on real-time registration between 2D ultrasound and 3D preoperative images Chijun Weon, Woo Hyun Nam, and Duhgoon Lee Department of Electrical Engineering, KAIST, Daejeon 305-701, Republic of Korea

Jae Young Lee Department of Radiology, Seoul National University Hospital, Seoul 110-744, Republic of Korea

Jong Beom Raa) Department of Electrical Engineering, KAIST, Daejeon 305-701, Republic of Korea

(Received 7 July 2014; revised 24 November 2014; accepted for publication 30 November 2014; published 24 December 2014) Purpose: Registration between 2D ultrasound (US) and 3D preoperative magnetic resonance (MR) (or computed tomography, CT) images has been studied recently for US-guided intervention. However, the existing techniques have some limits, either in the registration speed or the performance. The purpose of this work is to develop a real-time and fully automatic registration system between two intermodal images of the liver, and subsequently an indirect lesion positioning/tracking algorithm based on the registration result, for image-guided interventions. Methods: The proposed position tracking system consists of three stages. In the preoperative stage, the authors acquire several 3D preoperative MR (or CT) images at different respiratory phases. Based on the transformations obtained from nonrigid registration of the acquired 3D images, they then generate a 4D preoperative image along the respiratory phase. In the intraoperative preparatory stage, they properly attach a 3D US transducer to the patient’s body and fix its pose using a holding mechanism. They then acquire a couple of respiratory-controlled 3D US images. Via the rigid registration of these US images to the 3D preoperative images in the 4D image, the pose information of the fixed-pose 3D US transducer is determined with respect to the preoperative image coordinates. As feature(s) to use for the rigid registration, they may choose either internal liver vessels or the inferior vena cava. Since the latter is especially useful in patients with a diffuse liver disease, the authors newly propose using it. In the intraoperative real-time stage, they acquire 2D US images in real-time from the fixed-pose transducer. For each US image, they select candidates for its corresponding 2D preoperative slice from the 4D preoperative MR (or CT) image, based on the predetermined pose information of the transducer. The correct corresponding image is then found among those candidates via real-time 2D registration based on a gradient-based similarity measure. Finally, if needed, they obtain the position information of the liver lesion using the 3D preoperative image to which the registered 2D preoperative slice belongs. Results: The proposed method was applied to 23 clinical datasets and quantitative evaluations were conducted. With the exception of one clinical dataset that included US images of extremely low quality, 22 datasets of various liver status were successfully applied in the evaluation. Experimental results showed that the registration error between the anatomical features of US and preoperative MR images is less than 3 mm on average. The lesion tracking error was also found to be less than 5 mm at maximum. Conclusions: A new system has been proposed for real-time registration between 2D US and successive multiple 3D preoperative MR/CT images of the liver and was applied for indirect lesion tracking for image-guided intervention. The system is fully automatic and robust even with images that had low quality due to patient status. Through visual examinations and quantitative evaluations, it was verified that the proposed system can provide high lesion tracking accuracy as well as high registration accuracy, at performance levels which were acceptable for various clinical applications. C 2015 American Association of Physicists in Medicine. [http://dx.doi.org/10.1118/1.4903945] Key words: ultrasound, 4D preoperative image, image-guided intervention, real-time registration, lesion tracking 1. INTRODUCTION Real-time imaging of a target organ is essential for monitoring and tracking a lesion during image-guided intervention, especially if the lesion is moving, as in the liver. An ultra335

Med. Phys. 42 (1), January 2015

sound (US) scanner is considered particularly useful for intervention of the liver because it can provide real-time and noninvasive imaging. Hence, 2D US scanners are widely used as a guidance tool in intervention procedures such as percutaneous needle biopsy, radio frequency ablation (RFA),

0094-2405/2015/42(1)/335/13/$30.00

© 2015 Am. Assoc. Phys. Med.

335

336

Weon et al.: Position tracking of moving liver lesion based on real-time registration

high intensity focused ultrasound therapy (HIFU), and so on. However, US images have relatively lower image quality than magnetic resonance (MR) or computed tomography (CT) images, which occasionally makes it hard to clearly identify a target lesion in a US image. As a solution to overcome this limitation, registering a US image to a high quality preoperative MR (or CT) image in real-time has become important. The result of this registration can be helpful for radiologists or surgeons who plan the treatment of a target lesion by localizing it, because the registered high quality preoperative image can provide much clearer lesion information than the US image alone. Due to the reasons above, registration with both sufficient accuracy and high speed is needed. Various voxel-based or feature-based algorithms have been presented for performing rigid registration between US and preoperative images of the liver. Since the image characteristics of US and other preoperative images of the liver are considerably different from each other, voxel-based algorithms focus on correlating two images by converting the US image and/or preoperative image into a different type of image(s) and then defining a similarity measure based on the correlation.1,2 A feature-based rigid registration algorithm has also been proposed using the vessel centerline and liver surface as features.3 In this algorithm, however, a manual segmentation step is required for the extraction of features in both US and preoperative images. Since the liver undergoes respiratory motion with local deformation, several algorithms have tried nonrigid registration in order to compensate for the local deformation. An ICP-based nonrigid registration algorithm used an additional power Doppler (PD) US image for the registration between B-mode US and CT images.4 It utilized vessel information for the nonrigid transformation modeled with multilevel B-splines. A nonrigid registration algorithm using a PD US image was also developed based on combined landmark-intensity information.5 In that algorithm, however, a few corresponding landmark pairs in both PD US and CT images need to be manually selected in advance. Instead of using PD US images, B-mode US images were directly used in a couple of nonrigid registration algorithms.2,6 In one of them, an accurate nonrigid registration scheme was developed between B-mode US and CT images of the liver, by using a similarity measure based on the intensity and the magnitude and orientation information of gradients.6 The algorithms explained above are mainly focused on 3D–3D registration between US and other types of preoperative images. In the clinical applications of intervention, however, the simultaneous display of a realtime 2D US image with its registered preoperative image is very helpful for planning a surgical procedure and for guiding a surgical instrument. For this purpose, several commercially available fusion systems exist which are based on rigid registration, and their clinical feasibilities have been reported.7–10 However, even though those systems have shown clinical benefits for the intervention, they cannot handle the local deformation problem created by the respiratory motion of the liver. To solve this problem, another fusion system Medical Physics, Vol. 42, No. 1, January 2015

336

is tried on a nonrigid registration between the intraoperative 2D US image and a 4D MR image by using models of the respiratory motion and deformation.11 The fusion systems described above need real-time position tracking of the US transducer for initial registration and for subsequent registrations. For real-time tracking of the transducer, an optical or magnetic tracking device is usually used because they can provide highly accurate coordinates of the transducer in real-time. However, it is known that optical tracking may suffer from a “line of sight” problem, and magnetic tracking can be unstable due to interference with metallic objects nearby.12,13 In this paper, we propose a purely image-based registration system for the US-guided intervention of the liver, based on our previous works.13–15 The proposed system provides real-time registration between intraoperative 2D US and successive multiple 3D preoperative MR/CT images without using an external tracking device. In the system, we first generate a sequence of 3D preoperative images along the respiratory phase (that is, a 4D preoperative image), based on the nonrigid registration of several 3D preoperative images. In the intraoperative stage, we adopt a fixed-pose manner for a 3D US transducer for image acquisitions. A reliable rigid transformation between intraoperative US and preoperative image coordinates is then determined using a couple of 3D US images obtained in this stage. Based on the obtained rigid transformation, we suggest a simple real-time registration procedure between 2D intraoperative US and preoperative MR/CT images. This registration result can be used for the simultaneous display of US and MR/CT images or for the real-time position tracking of a liver lesion. The lesion position can not only be estimated directly based on a 2D US image which includes the lesion but also indirectly based on a different 2D US image that does not include a lesion, but which does include features sufficient for registration. A laboratory version of the system was implemented and evaluated using clinical patient datasets. This paper is organized as follows. In Sec. 2, we introduce an overview of the proposed registration system and describe several algorithms adopted in the system in detail. In Sec. 3, experimental results for clinical datasets are provided, and the following discussions are described in Sec. 4. Finally, conclusions are given in Sec. 5.

2. METHODS In this section, we describe the proposed system for realtime registration and lesion tracking. The system aims to estimate the position of a moving liver lesion based on real-time registration between intraoperative 2D US and preoperative MR/CT images. In the system, to avoid using physical tracking methods such as optical and magnetic systems, we try to rely only on real-time image registration based on image voxel information. The overall procedure for using the proposed system consists of three parts: the preoperative, intraoperative preparatory, and intraoperative real-time stages, as shown in Fig. 1.

337


337

2.A.1. 4D image generation

F. 1. Overview of the proposed real-time US-preoperative image registration system.

2.A. Preoperative stage

This stage is for obtaining a 4D preoperative image which consists of a series of successive 3D images along the respiratory phase. We may use MR or CT images as the preoperative data. A 4D preoperative image may be acquired via gated imaging.16,17 Even though gated 4D imaging reflects the accurate actual movement of the organ, it requires a substantially long data acquisition time along with an additional device for respiratory-gating. Therefore, instead, we generate a 4D image using several 3D images. First, 3D preoperative images are acquired at three different respiratory phases under breath-hold. The end-inspiration, mid-inspiration, and end-expiration may be selected as the three phases. A 4D preoperative image is then generated by using the transformation parameters, which are determined by performing nonrigid registration between each nearby pair of 3D images. The generated 4D image can be regarded as densely sampled 3D images along the time axis.

Let us assume that Ipre N , Ipre N , and Ipre N denote 1 2 3 the 3D preoperative images sequentially acquired from end-expiration, mid-inspiration, and end-inspiration. After denoising those 3D preoperative images, we manually extract anatomical features such as the liver surface, and the vessels or the inferior vena cava (IVC), only from a single image (we choose Ipre N in this paper). In addition, the 3D position 1 information of a target lesion in the single image is extracted for real-time processing in the intraoperative stage. To generate a 4D image, we first need to determine the transformations between adjacent images. Two registration procedures are performed for determining the transformations, as shown in Fig. 2. Note here that the transformation from Ipre N to Ipre N is denoted by T Nk→Nk +1, where k = 1,2. A k k +1 registration procedure uses normalized mutual information (NMI) as a similarity measure and consists of two steps; the global rigid registration based on the affine transformation and the local nonrigid registration based on the B-spline free form deformation model.18 In both registration steps, a multiresolution approach is commonly used. Interpolating the transformations along the time axis, we generate successive multiple 3D preoperative images along the respiratory phase, which is called a 4D preoperative image.19 In other words, the intermediate image Ipren between Ipre N and Ipre N can be determined on every integer point k k +1 xn as ( ) (xn − xˆ n )2 1  exp − Ipre N x Nk , (1) Ipren (xn ) = 2 k Q xˆ ∈R 2σ n

xn

where ( ) (xn − xˆ n )2 exp − , 2σ 2 xˆ n ∈Rx n xˆ n = T Nk→n x Nk = x Nk + d Nk→n x Nk = x Nk + w n,k · d Nk→Nk +1 x Nk

Q=



(0 ≤ w n,k ≤ 1).

(3)

In the equations above, xˆ n denote noninteger voxel positions in Ipren corresponding to integer voxel positions x Nk in Ipre N . In addition, Rxn denotes a set of xˆ n ’s neighboring k xn , σ is the standard deviation of the Gaussiandistribution, w n,k denotes a weight factor, and d Nk→Nk +1 x Nk denotes the

F. 2. Registration process for 4D preoperative image generation. Medical Physics, Vol. 42, No. 1, January 2015

(2)

338


displacement vector from Ipre N to Ipre N k

w n,k =

338

. We define w n,k as

k +1

n − N N1, Nk , N Nk, Nk +1 − 1

(4)

where N N1, Nk and N Nk, Nk +1 denote the number of generated images between Ipre N and Ipre N and between Ipre N and 1 k k Ipre N , respectively. k +1 To maintain the naturalness and continuity of the 4D preoperative image, we regard Ipre N as a single reference 1 image in the 4D image generation. Intermediate images between Ipre N and Ipre N can then be directly produced 1 2 by using the 4D image generating algorithm. Intermediate images at time n between Ipre N and Ipre N , however, can be 2 3 determined by redefining the transformation as xˆ n = T N1→N2 x N1 + w n,2 · d N2→N3 T N1→N2 x N1 . (5) We set the total number of interpolated images, N, in the 4D image so that the average amount of liver motion between any adjacent images can have a uniform value δ of 1 mm. Hence, the interpolated images usually do not have a uniform sampling time interval. Thereby, we should note that the 4D image does not include real respiratory time information. All 3D images in the 4D image are then denoted by Ipre1,Ipre2,...,Ipre N from the end-expiration to the endinspiration. In addition, using the information predetermined from Ipre N , we generate the 4D anatomical features and the 1 4D position information of a target lesion in a similar way for the following intraoperative registration. 2.B. Intraoperative preparatory stage

This stage corresponds to the preparation for real-time image registration. In this stage, a 3D US transducer is mounted by to a holding mechanism. A 3D region of interest (ROI), or the US transducer, is then posed and fixed so that it can include target features for the registration, such as internal liver vessels, IVC, and liver surface. We acquire a couple of 3D US images at both expiration and inspiration phases under breath-hold. Via the registration of the two US images to 3D preoperative images in the 4D image, a fixed transformation between the 3D US coordinates (or transducer coordinates) and the preoperative image coordinates is determined. Note that the transformation is consistent during the whole intraoperative procedure and is one of the major factors to avoid a need of an external tracking system for real-time registration between 2D US and 4D preoperative images. Even though the transducer is fixed during the intraoperative procedure, a 2D ROI is still flexible to change within its 3D ROI. Instead of two 3D US images, a single 3D image at an arbitrary phase may be enough to determine the transformation by registering it to the 3D preoperative image corresponding to the same respiratory phase. It was found, however, that a transformation obtained using a single image is often not reliable enough to get accurate pose information. Hence, to determine a reliable transformation, we use a couple of 3D US images at different respiratory phases of expiration and inspiMedical Physics, Vol. 42, No. 1, January 2015

F. 3. Block diagram of the IVC extraction algorithm from a 3D B-mode US image.

ration. Based on the fact that the transformation is the same regardless of the phase change due to the fixed transducer pose, we determine a transformation by simultaneously registering the two US images to preoperative images. 2.B.1. Rigid registration

For each of the two US images, we first perform a rigid registration with every 3D image in the 4D image generated in the preoperative stage. For rigid registration, we may adopt our previously proposed algorithm based on the anatomical features of vessels and the liver surface.13 If a patient has a diffuse liver disease such as fatty liver, however, the internal liver vessels are not clearly visible in an US image. We fortunately observe in this case that the IVC and the liver surface can be visible. Therefore, to increase the robustness of the registration algorithm, we adopt the IVC and the liver surface as anatomical features if vessels are hardly visible in the US liver image.15 In this subsection, we will propose a newly adopted rigid registration algorithm based on the IVC and liver surface, which consists of three steps; the feature extraction, initial registration, and refinement. 2.B.1.a. Feature extraction. For feature-based initial registration, the features need to be extracted in both preoperative and US images. Even though the features of the preoperative image were extracted manually in the preoperative stage, the features of US images must be automatically extracted because they are acquired in the intraoperative stage. A liver surface (or diaphragm) extraction algorithm had already been developed in our previous research.13 We therefore needed to develop an extraction algorithm of the IVC, which is suggested for use in this paper. The block diagram of the proposed US IVC extraction algorithm is shown in Fig. 3. For extraction, we first confine the ROI for extracting the IVC. Since the IVC is close to the liver surface, we determine the ROI by eliminating the region above the surface with a marginal distance of 30 mm and the region under the surface with a distance of 10 mm. Inside the determined ROI, we binarize IVC candidate regions

339


by applying adaptive thresholding based on a local sliding window.13,20 To select the IVC among the candidate regions, we propose an IVC metric based on our observations that the IVC is longer than noise clusters and locally straight and is thicker than other liver vessels. The proposed metric can be written as VIVC = S1 + αS2 + βS3.

(6)

Here, α and β denote weight factors, S1 and S2 represent the global lengthiness and the local linearity of the IVC, respectively, and S3 represents the thickness similarity between the extracted IVC of the preoperative image and the candidates. The terms S1 and S2 can be determined by fitting a line to a candidate region of IVC via a principal component analysis (PCA). We define S1 as w (7) S1 = − , l where l and w denote the length and width of the candidate, which are determined along the first component vector w1 and the second component vector w2 obtained from the PCA, respectively, as shown in Fig. 4. Term S2 is defined as the average distance between the local center points that are obtained from the distance transform21 and the fitted line, as shown in Fig. 4, and can be written as  1 S2 = − xk,lcp − xk,proj , (8) NUS_lcp k ∈N US_lcp

where xk,lcp denote local center points, xk,proj denote their projections on the fitted line, and NUS_lcp is the number of local center points of the candidate in the US image. Meanwhile, thickness similarity S3 is defined as   1 1 S3 = − pk,pre − pk,US , (9) N N pre_lcp US_lcp k ∈Npre_lcp k ∈NUS_lcp where pk,pre and pk,US denote the thicknesses obtained by using the distance transform, and Npre_lcp is the number of local center points of the segmented IVC in the preoperative image. The IVC can then be determined as the one which

provides the maximum value of VIVC, among the candidates extracted from the US image. 2.B.1.b. Initial registration. Using the extracted features of the IVC and liver surface, we can perform the initial registration between the preoperative image and the US image. A block diagram of the initial-registration algorithm is shown as Fig. 5. We have determined a fitted line (or a center line) for the presegmented IVC in the preoperative image. In the previous feature extraction step, we also determined a fitted line (or a center line) for the IVC in a US image. Aligning the two lines using an Euler angle rotating matrix, we can determine the IVC alignment matrix TIVCpre→IVCUS from the preoperative image domain to the US image domain.22 Once the IVC center lines are aligned, only two parameters, or the translation parameter along the IVC center line of the US image and the rotation parameter around the line, are needed for the final transformation, rather than the conventional six-rigid transform parameters. To determine the two parameters after the alignment, we propose to use geometric distance maps, as shown in Fig. 6. To generate a geometric distance map between the center line of the IVC and the liver surface, the y axis of the center line and angle φ perpendicular to the y axis are sampled with the sampling intervals of ∆ y and ∆φ, respectively. The sampled points and angles can then be written as yUS = i · ∆ y, i = 0,...,N y,US, ypre = i · ∆ y, i = 0,...,N y,pre, φUS = φpre = j · ∆φ, j = 0,...,Nφ .

(10)

Note here that the number of total center line sampling points of the preoperative image, N y,pre, is always larger than that of US image, N y,US, because the preoperative image usually captures a larger portion of the IVC. Meanwhile, the number of total angle sampling directions, Nφ , is set to the same in both images. We set ∆ y to 0.77 mm, so as to be the same as the voxel width of the US image of lower resolution, and ∆φ to 1◦. Using the sampling coordinates (i, j), we generate the US image distance map d US and the preoperative image distance map d pre, respectively.

F. 4. PCA based measures for the IVC metric on (a) IVC, (b) a liver vessel, and (c) nonvessel cluster. Medical Physics, Vol. 42, No. 1, January 2015

339

340


340

ˆ can then be The translation and rotation parameters, yˆ and φ, obtained as yˆ = iˆs · ∆ y,

φˆ = jˆs · ∆φ.

(13)

Finally, by combining the obtained translation and rotation parameters with the predetermined vector alignment matrix, an initial transformation can be described with six-rigid transformation parameters as ( ) xˆ = Tinitial xpre = T y,ˆ φˆ TIVCpre→IVCUS xpre , (14) where xpre and xˆ denote a point in a preoperative image and its corresponding point in the 3D US coordinates, respectively, and Tinitial denotes the initial transformation matrix. 2.B.1.c. Registration refinement. Once the initial registration is done, we refine the result by adopting an image-voxel based registration algorithm, which we had previously proposed.7,23 The algorithm employs a gradient descent optimization to minimize a measure based on the intensity and gradient orientation of voxels. To define the measure, we use a 3D joint histogram, h(i US,i pre,o(∆θ)), where i US and i pre denote the intensity values of the US and preoperative images, respectively, and o(∆θ) denotes the orientation coincidence of the two images that can be represented as

F. 5. Block diagram of the initial-registration algorithm.

To determine the translation and rotation parameters using the distance maps, we first coincide their respective origins. We also define a normalized SAD measure (NSAD) as NSAD(i s, js) =

1 Noverlap

 i

d US(i −i s, j − js) − d pre(i, j) ,

j

i s = 0,...,N y,pre − N y,US, js = 0,...,Nφ ,

(11)

where (i s , j s ) denote shift amounts between the two maps along the i- and j-axes, respectively. We then align the two maps by determining (i s , j s ) which minimizes NSAD. Namely, iˆs, jˆs = arg min NSAD(i s, js).

i s, js

(12)

1 + cos(2∆θ) . (15) 2 Here, ∆θ is the gradient orientation difference at each voxel between two images. The measure F is then defined as F(IUS, Ipre) = 1 −W IUS,Ipre · E IUS,Ipre , (16) o(∆θ) =

where  W IUS,Ipre =

x ∈( R US∩R pre) 2 · o(∆θ x )

, NRUS + NRpre E IUS, Ipre = H IUS,Ipre,O − M IUS,Ipre .

(18)

Here, H(IUS, Ipre,O) denotes the 3D joint entropy calculated from the 3D joint histogram, and M(IUS, Ipre) represents the

F. 6. Generation of two geometric distance maps between IVC and liver surface in the IVC-aligned US and preoperative images. Medical Physics, Vol. 42, No. 1, January 2015

(17)

341


mutual information of pixel intensities of US and preoperative images. In addition, RUS and Rpre denote interested regions such as edge regions of vessels and the liver surface in the two images; and NRUS and NRpre denote the numbers of voxels in RUS and Rpre, respectively. We should note that the image characteristics of the two images are different especially at the liver surface. The liver surface (or diaphragm) provides high intensity values in the US image due to its strong reflection. The liver surface in the preoperative image, however, shows a definite boundary of an edge type and provides high gradient magnitude.6 Therefore, when the measure in Eq. (16) is calculated in the liver surface of the preoperative image, the gradient magnitude is used rather than the intensity value.

2.B.2. Transformation parameter determination

To determine a fixed transformation between the 3D US coordinates and the preoperative image coordinates, rigid

e Tn→m

registrations are performed from the preoperative 4D image to two US images. Using the rigid registration scheme described above, we first obtain rigid transformations Tn→m from each 3D preoperative image to the two US images, where m (= 1,2) and n (= 1,...,N) denote the indices of US and preoperative images, respectively. We assume that if the respiration phases of a 3D US image and a 3D preoperative image become similar, a rigid transformation between the two images will have a small registration error. This assumption is considered true, because a nonrigid deformation due to a respiratory phase difference between those images cannot be compensated by the rigid registration. Based on this assumption, we determine the fixed transformation by selecting a rigid transformation among a number of transformation candidates Tn→m , which minimizes the sum of registration errors for two 3D US images IUS1 and IUS2. Namely, Ttransducer = arg min {e(Tn→m )},

where δUS and δ denote the average amount of liver motion between IUS1 and IUS2 and between adjacent 3D preoperative images, respectively. We may determine δUS by directly registering two 3D US images. Since this registration is hard to achieve due to small FOVs, however, we use two rigid transformations from Ipren′ to each 3D US image for an arbitrary n ′. Namely,  1 δUS = Tn′→1 xk,pren′ − Tn′→2 xk,pren′ , (22) Npren′ k ∈N pre n ′

where xk,pren′ denotes the kth point in Ipren′, and Npren′ denotes the number of voxels of the liver region in Ipren′. Note here that δUS is almost the same for any n ′. Using the inverse matrix of Ttransducer obtained in Eq. (19), we convert the coordinates of the prestored positions of target Medical Physics, Vol. 42, No. 1, January 2015

(19)

Tn→m

where

1  if m = 1   2 C IUS1,Ipre n ;Tn→1 +C IUS2,Ipre n+n d;Tn→1 , =  1  C IUS ,Ipre ;Tn→2 +C IUS ,Ipre ;Tn→2 , if m = 2. 2 1 n n−n d 2

Here, C represents the registration error of two images for a given transformation. If liver vessels are well visible in the US images, we may use the error measure defined in our previous rigid registration algorithm13 based on the internal liver vessels and liver surface as C. Otherwise, the measure F based on the IVC and liver surface in Eq. (16) can be used as C. In addition, nd denotes the difference of respiratory phase indices of 3D preoperative images corresponding to IUS1 and IUS2, which can be written as   δUS nd = + 0.5 , (21) δ

341

(20)

lesion in all 3D preoperative images of different respiratory phases to the US image coordinates so that they can be used for real-time positioning.14 2.C. Intraoperative real-time stage

In this real-time stage, we acquire real-time 2D US images using the fixed 3D US transducer. Even though the 3D US transducer is fixed with a known transformation matrix, Ttransducer, we still have the freedom to select a desired realtime 2D US image plane within the 3D FOV, by changing the scanning angle. Referring to the scanning angle value of a current 2D US image, we select the corresponding 2D preoperative slice at every respiratory phase in the 4D preoperative image transformed to the US image coordinates. Among those images, we can determine a matched 2D image and its corresponding phase in real-time. The details are as follows. 2.C.1. Best-matched 2D preoperative slice selection

Since transformation T2D, which transforms the 3D US coordinates to 2D US coordinates, is intrinsically known, the preoperative image corresponding to a 2D US image can be obtained as xˆ 2D = T3D→2D xpre = T2D Ttransducer xpre , (23) where xpre denotes a point in the preoperative image coordinates, and xˆ 2D denotes the corresponding point in the 2D

342


342

F. 7. Selection of best-matched 2D preoperative slice.

US image coordinates. Based on Eq. (23), for a real-time 2D US image, we can extract 2D preoperative slice candidates from the 3D preoperative images covering all M respiratory phases, respectively, as shown in Fig. 7. We then determine the best-matched candidate by comparing the real-time US image to M 2D preoperative slices. As a similarity measure for the matching process above, we adopt a simple gradient-based one6 similar to Eq. (17) by considering both the computation time and multimodal image characteristics. Namely,   x2D ∈(R′US R′pre) 1 + cos 2∆θ x2D S2D = , (24) NR′US + NR′pre ′ ′ where RUS and Rpre denote interesting regions such as the edge regions of vessels and the liver surface in the 2D US images and preoperative slices, respectively; NR′US and NR′pre ′ ′ and Rpre , respectively; and denote the number of pixels in RUS ∆θ x2D is the difference between gradient orientation angles at the corresponding points in the two images. Finding the maximum out of M different values of S2D, we can select a correct 2D preoperative slice which corresponds to the real-time 2D US image. It should be noted that since this procedure is based only on the selection of a 2D slice among the limited number of candidates, it can be performed in realtime.

2.C.2. Lesion position estimation

If a real-time 2D US image plane including a lesion has sufficient features, we can directly determine the bestmatched preoperative 2D slice in real-time and thereby track Medical Physics, Vol. 42, No. 1, January 2015

the lesion position along the respiratory motion. However, if none of the real-time 2D US image planes which include a lesion have features that can be used for registering to preoperative images, we then have to choose a 2D US image plane which does not include a lesion, but which does have better features for registration. The position of the target lesion can then be determined using the registered preoperative image and its displacement from the lesion,14 which has been predetermined in the intraoperative preparatory stage described in Sec. 2.B.2. We can then simultaneously display the real-time 2D US image, its matched 2D preoperative slice, and the 2D preoperative slice, including a lesion with a predetermined position shift. 3. EXPERIMENTS 3.A. Data acquisition

For experiments, we acquired 23 datasets from patient volunteers. We selected volunteers with a relatively low body mass index to avoid difficulty in US imaging due to thick subcutaneous fat. The study protocol was approved by the Institutional Review Board of Seoul National University Hospital, and written informed consents were obtained from each volunteer. Ten patients among the volunteers had cysts in their livers and thirteen patients had hepatocellular carcinoma(s) (HCC). Near the delay phase, three contrast-enhanced preoperative 3D MR images were acquired at different respiratory phases by using a Siemens MAGNETOM Trio 3T scanner. The typical dimension of 3D MR images is 384×312 × 60 with a voxel size of 0.99 × 0.99 × 3.0 mm3. Meanwhile, the intraoperative US images were acquired by using

343


343

F. 9. US image of D23 with serious shadowing artifact and blur.

For the rigid registration mentioned above, we first applied our previously proposed algorithm based on internal liver vessels and the liver surface.13 Among 23 datasets we then obtained a successful registration in 13 datasets, all of which were free from diffuse liver disease such as a fatty liver. We subsequently applied the proposed algorithm for the remaining 10 datasets, and the registration was successful for 9 datasets. One dataset, D23, included an unusually serious shadowing artifact and blur, as shown in Fig. 9. Table II shows the quantitative rigid registration error in the first column, in terms of the average distance between six corresponding fiducial points on anatomical features. The average error of 3D registration is 2.80 mm which is considered acceptable.4 We also examined the registration performance between 2D US images acquired in real-time and their corresponding 2D preoperative MR images extracted from the 4D preoperative image. For reliable registration, we selected a scanning angle that enabled us to acquire 2D US images which included sufficient features with no lesion(s), rather than having insufficient features with lesion(s). Figure 10 shows the 2D USimages at exhale and inhale phases selected from two 2D US sequences D6 and D21, respectively, and their registered preoperative MR images. To permit simple visual assessment of the registration accuracy, we indicate the corresponding locations with several pairs of arrows in both images. To demonstrate how US and

F. 8. US image acquisition with an X6-1 3D transducer attached to a holding mechanism.

a Philips IU-22 ultrasound scanner equipped with an X6-1 3D transducer based on a 2D array sensor. The US transducer was positioned on the intercostal space so as to minimize any contact problem during the free-breathing acquisition. The typical dimension of 3D US images is 512 × 378 × 222 with a voxel size of 0.43 × 0.39 × 0.63 mm3. As a mechanism for holding the US transducer in the intraoperative stage, we used a five-axis arm, as shown in Fig. 8. 3.B. Experimental results

In the preoperative MR images acquired for generating a 4D image, the average amount of liver motion between endexpiration and end-inspiration was examined, as in the first column in Table I. To maintain the average liver motion between any adjacent 3D images with about 1 mm in the 4D image, the number of generated 3D images was set to the values in the second column of the table. To obtain Ttransducer, we acquired two 3D US images at the respiratory phases close to the endexpiration and end-inspiration, respectively, and performed the rigid registration with the 4D preoperative image.

T I. Liver motion between end-expiration and end-inspiration and the number of generated 3D MR images.

D1 D2 D3 D4 D5 D6 D7 D8 D9 D10 D11 D12

Average liver motion (mm)

Number of generated 3D images

22.56 44.06 47.41 6.22 17.47 26.39 50.32 24.04 17.38 33.71 19.41 24.98

23 45 48 7 18 27 51 25 18 34 20 25

Medical Physics, Vol. 42, No. 1, January 2015

D13 D14 D15 D16 D17 D18 D19 D20 D21 D22 D23

Average liver motion (mm)

Number of generated 3D images

45.64 47.49 27.01 27.15 23.36 20.44 7.14 17.84 52.61 21.41 47.08

46 48 28 28 24 21 8 18 53 22 48

344


T II. Average fiducial registration error (FRE) for 3D registration. Avg. 3D registration error (± STD) (mm) D1 D2 D3 D4 D5 D6 D7 D8 D9 D10 D11 D12

2.02 (±0.91) 3.64 (±1.41) 2.57 (±1.04) 3.22 (±1.88) 3.11 (±0.95) 2.64 (±0.66) 3.99 (±1.57) 3.45 (±1.37) 2.52 (±1.55) 2.23 (±1.53) 3.87 (±2.97) 2.16 (±1.22)

Avg. 3D registration error (± STD) (mm) D13 D14 D15 D16 D17 D18 D19 D20 D21 D22 D23 Average

3.44 (±1.56) 1.14 (±0.36) 2.92 (±2.01) 3.63 (±2.50) 3.08 (±1.33) 2.81 (±1.54) 1.67 (±0.95) 1.78 (±0.89) 2.62 (±1.22) 3.07 (±2.30) — 2.80 (±1.44)

preoperative MR images are registered, Fig. 10 also includes the curves of the values of similarity measure, according to the respiratory phase of the 2D preoperative MR image

candidates, for two US images with different phases in each dataset. Here, the maximum values of similarity measure correspond to the registered preoperative MR images, respectively. Table III shows the quantitative 2D registration error for each clinical dataset. The error is determined as the average of the distances between 16 corresponding fiducial points at four selected respiratory phases, which are represented as end-expiration, shallow-expiration, shallowinspiration, and end-inspiration. The overall average error of 2D registration is found to be 2.42 mm. To validate the accuracy of the proposed indirect lesion position estimation approach, Fig. 11 shows the registered 2D US and preoperative MR images with no lesion(s) in D5 and the 2D US and preoperative MR images including a lesion, which were consequently estimated by determining the difference in viewing angle, as was described in Sec. 2.C.2. To determine the performance of the registration, several pairs of arrows were marked and examined in both images. The results in the figure demonstrate good registration accuracy, and the target lesion appears well aligned.

F. 10. 2D registration results for D6 and D21 and their similarity measure graphs. Medical Physics, Vol. 42, No. 1, January 2015

344

345


T III. Average FRE for real-time 2D registration. Avg. 2D registration error (±STD) (mm) D1 D2 D3 D4 D5 D6 D7 D8 D9 D10 D11 D12

2.31 (±1.21) 3.02 (±0.75) 2.24 (±0.77) 3.12 (±1.24) 2.87 (±1.13) 3.07 (±0.94) 3.74 (±1.64) 3.65 (±1.43) 2.77 (±0.99) 1.71 (±0.37) 2.45 (±0.76) 1.44 (±0.48)

Avg. 2D registration error (±STD) (mm) D13 D14 D15 D16 D17 D18 D19 D20 D21 D22 D23 Average

1.57 (±0.67) 1.72 (±0.64) 2.84 (±1.01) 2.44 (±1.00) 2.74 (±0.91) 1.99 (±0.44) 0.76 (±0.27) 1.24 (±0.37) 2.49 (±1.36) 3.14 (±1.31) — 2.42 (±0.87)

For quantitative evaluation, we used 17 datasets including lesion(s) which was visible in both 3D US and 3D preoperative MR images. The lesions were manually segmented in one of the two 3D US images and delineated. Based on the registration results given in Table III, which were produced using the 2D US images with no lesion(s), we then measured the 3D distances between the centers of lesions which had been manually segmented in the 3D US image and the lesion centers which had been estimated in each 3D preoperative MR image. Table IV shows the 3D distances or the estimation errors of the lesion positions for the 17 datasets. As shown in the table, the average estimation error is 3.27 mm. This error, which is less than 5 mm, is considered acceptable for clinical applications.24

4. DISCUSSIONS We proposed a system for real-time position tracking of a moving liver lesion. The algorithm relies on image-based

registration, which means it does not require any external sensing device or any inconvenient manual interaction in the intraoperative stage. For real-time processing, prior to the intraoperative real-time stage, we generate a 4D preoperative image and also determine the fixed probe position. Thereby, we can register an incoming 2D US image to the corresponding 3D preoperative image of the same phase in real-time, by selecting it out of the generated 4D image. Through qualitative and quantitative assessments of the registration results, we found that the proposed system can provide registration performance good enough for imageguidance lesion tracking. We also found that the computation time for the registration of a 2D US image and the estimation of a target lesion position in the intraoperative stage is about 60 ms in a standard PC with a 3.2 GHz CPU, which is considered fast enough for real-time operation of the system. Meanwhile, the computation time in the intraoperative preparatory stage is about 3 min, which can be acceptable for the operation preparation. To increase the robustness of the proposed system, we introduced a new rigid registration algorithm based on the IVC and liver surface, and an indirect tracking system of the lesion position. By adopting the proposed rigid registration algorithm in addition to our previously proposed rigid registration algorithm, 22 datasets out of 23 acquired datasets were found to be applicable using the proposed system, which can be considered robust enough for various clinical cases. Furthermore, the proposed indirect position tracking system enables the estimation of the position of a moving liver lesion, even though a target US image including a lesion does not have enough anatomical features. Therefore, it may considerably increase the treatment performance of image-guided therapy systems such as the HIFU system, Cyberknife, and so on. In this study, to estimate the US transducer pose, we first applied our previous algorithm based on internal liver vessels and then applied the proposed registration algorithm to the

F. 11. (a) 2D registration result and (b) the determined images including a lesion, for D5. Medical Physics, Vol. 42, No. 1, January 2015

345

346


346

T IV. Estimated positioning error for every selected lesion.

D2 D4 D5 D6 D8 D9 D10 D11 D12

Lesion type

Estimated lesion position error (mm)

Cyst Cyst Cyst Cyst Cyst HCC HCC HCC HCC

2.39 3.31 2.71 3.04 4.61 4.12 3.63 3.59 2.88

remaining unsuccessful datasets. We may change the order in which algorithms are applied. Applying the proposed rigid registration algorithm first, we obtained successful registration results for 17 datasets. The remaining six datasets were unsuccessful, because three datasets had an invisible IVC due to a shadowing artifact or blur in the US images and the other three datasets did not include the IVC in their US image ROI. Among them, five datasets were successful because the internal liver vessels were visible. Consequently, we can note that only one dataset was unsuccessful, regardless of the order the algorithms were applied. The proposed system acquires three 3D preoperative MR images in a breath-hold manner, to generate consecutive 3D images (or a 4D image). To guarantee the accuracy of the real-time 2D registration or lesion positioning, the 4D image should include a free-breathing range of the intraoperative stage of a patient. Therefore, 3D preoperative images should be acquired with sufficient respiratory differences. Organ deformation due to the pressure of a US probe is not considered in the proposed system. In general, to acquire high quality US images, radiologists or surgeons press a US probe to a patient’s body, which may cause deformation of the liver near the surface. If a target lesion is located at that deformed liver area, the lesion position may not be accurately determined, because the liver is usually not deformed in the preoperative MR image. Methods to compensate this deformation need to be studied in further work. 5. CONCLUSIONS We proposed a real-time image-based system to provide accurate registration between 2D US images of the liver and 3D images in a preoperative 4D liver image of a different modality. Based on the real-time registration, we then suggested a system to track the position of a liver lesion. The registration procedure consists of preoperative, intraoperative preparatory, and intraoperative real-time stages. In the first stage, we generate a 4D image based on liver motion. The 4D image consists of consecutive 3D preoperative images having a uniform amount of liver motion. In the second stage, we fix a 3D US transducer to a mechanical holder so that its FOV can include a lesion and determine the transformation Medical Physics, Vol. 42, No. 1, January 2015

Lesion type D13 D14 D15 D16 D17 D18 D21 D22 Average

Estimated lesion position error (mm)

HCC HCC HCC HCC Cyst HCC HCC Cyst

4.07 1.86 3.01 3.37 2.74 3.51 3.43 3.25 3.27 (±0.65)

between the corresponding 3D US image coordinates and the preoperative image coordinates via rigid image registration. To increase the robustness of the transformation determination, we introduced a liver feature, the IVC, to the registration. In the final stage, a 2D preoperative slice which corresponds to an incoming 2D US image is determined in real-time, by examining the similarity measures of the candidates for different respiratory phases. Using the transformation information corresponding to the determined 2D preoperative slice, we also suggest an indirect scheme for real-time tracking of a target lesion. This indirect approach enables a liver lesion to be tracked even when a target US image which includes the lesion is inadequate for registration, due to the lack of anatomical features. Based on experiments which determined the method produce sufficient registration and tracking accuracies, the proposed system is considered very prospective for real-time image-guided therapy procedures.

ACKNOWLEDGMENT The authors would like to thank Samsung Advanced Institute of Technology (SAIT) for the financial support of this work.

a)Author

to whom correspondence should be addressed. Electronic mail: [email protected] 1G. P. Penney, J. M. Blackall, M. S. Hamady, T. Sabharwal, A. Adam, and D. J. Hawkes, “Registration of freehand 3D ultrasound and magnetic resonance liver images,” Med. Image Anal. 8, 81–94 (2004). 2W. Wein, S. Brunke, A. Khamene, M. R. Callstrom, and N. Navab, “Automatic CT-ultrasound registration for diagnostic imaging and image-guided intervention,” Med. Image Anal. 12, 577–585 (2008). 3G. P. Penney, J. M. Blackall, D. Hayashi, T. Sabharwal, A. Adam, and D. J. Hawkes, “Overview of an ultrasound to CT or MR registration system for use in thermal ablation of liver metastases,” in Proceedings of Medical Image Understanding Analysis (The University of Birmingham, Birmingham, UK, 2001). 4T. Lange, N. Papenberg, S. Heldmann, J. Modersitzki, B. Fischer, H. Lamecher, and P. M. Schlag, “3D ultrasound-CT registration of the liver using combined landmark-intensity information,” Int. J. Comput. Assist. Radiol. Surg. 4, 79–88 (2009). 5L. Crocetti, R. Lencioni, S. DeBeni, T. See, C. Pina, and C. Bartolozzi, “Targeting liver lesions for radiofrequency ablation: An experimental feasibility study using a CT-US fusion imaging system,” Invest. Radiol. 43, 33–39 (2008).

347


6D. Lee, W. H. Nam, J. Y. Lee, and J. B. Ra, “Non-rigid registration between

3D ultrasound and CT images of the liver based on intensity and gradient information,” Phys. Med. Biol. 56, 117–137 (2011). 7W. Birkfellner, F. Watzinger, F. Wanschitz, R. Ewers, and H. Bergmann, “Calibration of tracking systems in a surgical environment,” IEEE Trans. Med. Imaging 17, 737–742 (1998). 8T. M. Peters, “Image-guidance for surgical procedures,” Phys. Med. Biol. 51, 505–540 (2006). 9M. Schneider and C. Stevens, “Development and testing of a new magnetic tracking device for image guidance,” Proc. SPIE 6509, 17–22 (2007). 10A. Cifor, L. Risser, M. P. Heinrich, D. Chung, and J. A. Schnabel, “Rigid registration of untracked freehand 2D ultrasound sweeps to 3D CT of liver tumours,” Abdom. Imaging 8198, 155–164 (2013). 11J. M. Blackall, G. P. Penney, A. P. King, and D. J. Hawkes, “Alignment of sparse freehand 3-D ultrasound with preoperative images of the liver using models of respiratory motion and deformation,” IEEE Trans. Med. Imaging 24, 1405–1416 (2005). 12N. D. Glossop, “Advantages of optical compared with electromagnetic tracking,” J. Bone Joint Surg. Am. 91, 23–28 (2009). 13W. H. Nam, D.-G. Kang, D. Lee, and J. B. Ra, “Anatomical registration between 3D intra-operative ultrasound and pre-operative CT images of the liver based on robust edge matching,” Phys. Med. Biol. 57, 69–91 (2012). 14C. Weon, W. H. Nam, D. Lee, Y. Hwang, J.-B. Kim, W.-C. Bang, and J. B. Ra, “Position estimation of moving liver lesion based on registration between 2D ultrasound and 4D MR images,” in Proceedings of the IEEE International Conference on Image Processing (IEEE, Orlando, FL, 2012), pp. 1677–1680. 15C. Weon, W. H. Nam, Y. Hwang, J.-B. Kim, W.-C. Bang, and J. B. Ra, “Robust feature based pre-registration of 3D MR image to 3D B-mode

Medical Physics, Vol. 42, No. 1, January 2015

347

ultrasound image of the liver,” in Proceedings of the International Society for Magnetic Resonance in Medicine (ISMRM, Salt Lake City, UT, 2013), p. 1842. 16C. Wachinger, M. Yigitsoy, E.-J. Rijkhorst, and N. Navab, “Manifold learning for image-based breathing gating in ultrasound and MRI,” Med. Image Anal. 16, 806–818 (2012). 17Y. Tong, J. K. Udupa, K. C. Ciesielski, J. M. Mcdonough, A. Mong, and R. M. Campbell, “Graph-based retrospective 4D image construction from free-breathing MRI slice acquisitions,” Proc. SPIE 9038, 90380I (2014). 18D. Rueckert, L. Sonoda, C. Hayes, D. Hill, M. Leach, and D. Hawkes, “Nonrigid registration using free-form deformations: Application to breast MR images,” IEEE Trans. Med. Imaging 18, 712–721 (1999). 19W. H. Nam, I. J. Ahn, K. M. Kim, B. I. Kim, and J. B. Ra, “Motioncompensated PET image reconstruction with respiratory-matched attenuation correction using two low-dose inhale and exhale CT images,” Phys. Med. Biol. 58, 7355–7374 (2013). 20R. Fisher, S. Perkins, A. Walker, and E. Wolfart, Hypermedia Image Processing Reference, 1st ed. (Wiley & Sons, Inc., New York, NY, 1996). 21I. Ragnemalm, “The Euclidean distance transform in arbitrary dimensions,” Pattern Recognit. Lett. 14, 883–888 (1993). 22E. Lomonosov, D. Chetverikov, and A. Ekart, “Pre-registration of arbitrarily oriented 3D surfaces using a genetic algorithm,” Pattern Recognit. Lett. 27, 1201–1208 (2006). 23Y. S. Kim, J. H. Lee, and J. B. Ra, “Multi-sensor image registration based on intensity and edge orientation information,” Pattern Recognit. 41, 3356–3365 (2008). 24G. P. Penney, “Applications in image guided interventions,” in IEEE ISBI Tutorial Notes: Tutorial on Biomedical Image Registration (IEEE, Rotterdam, Netherlands, 2010).

Deformable registration of preoperative MR, pre-resection ultrasound, and post-resection ultrasound images of neurosurgery.

Continuous roadmapping in liver TACE procedures using 2D-3D catheter-based registration.

Development of fast patient position verification software using 2D-3D image registration and its clinical experience.

Ultrasound 2D Strain Estimator Based on Image Registration for Ultrasound Elastography.

The agreement between 3D, standard 2D and triplane 2D speckle tracking: effects of image quality and 3D volume rate.

A survey among Brazilian thoracic surgeons about the use of preoperative 2D and 3D images.

3D-2D registration for surgical guidance: effect of projection view angles on registration accuracy.

A Registration Method Based on Contour Point Cloud for 3D Whole-Body PET and CT Images.

2D registration.

Acceleration method of 3D medical images registration based on compute unified device architecture.

Semiautomatic registration of 3D transabdominal ultrasound images for patient repositioning during postprostatectomy radiotherapy.

3D reconstruction of a carotid bifurcation from 2D transversal ultrasound images.

Registration of 2D C-Arm and 3D CT Images for a C-Arm Image-Assisted Navigation System for Spinal Surgery.

3D Registration Accuracy in 4 Major Anatomic Sites Using On-Board 2D Kilovoltage Images for 6D Patient Setup.

Texture-based classification of different single liver lesion based on SPAIR T2W MRI images.

2D Registration of the Left Atrium: Fast versus Consistent.

2D Model-to-Image Registration for Quantitative Dietary Assessment.

Random Walk Based Segmentation for the Prostate on 3D Transrectal Ultrasound Images.

Aortic root 3D parametric morphological model from 2D-echo images.

3D joint speaker position and orientation tracking with particle filters.

Consequences of Intermodality Registration Errors for Intramodality 3D Ultrasound IGRT.

Feature point based 3D tracking of multiple fish from multi-view images.

Location detection and tracking of moving targets by a 2D IR-UWB radar system.

Calibration for position tracking of swept motor 3-D ultrasound.