Breast Tumor Localization by Prone to Supine Landmark Driven Registration for Surgical Planning

Breast cancer is the most common cancer in women worldwide. Screening programs and imaging improvements have increased the detection of clinically occult non-palpable lesions requiring preoperative localization. Wire guided localization (WGL) is the current standard of care for the excision of non-palpable carcinomas during breast conserving surgery. Due to the current limitations of intraoperative tumor localization approaches, the integration of multimodal imaging information may be especially relevant in surgical planning. This research proposes a novel method for performing preoperative image-to-surgical surface data alignment to determine the position of the tumor at the time of surgery and aid preoperative planning. First, the volume of the breast in the surgical position is reconstructed and a set of surface correspondences is defined. Then, the preoperative (prone) and intraoperative (supine) volumes are co-registered using landmark driven non-rigid registration methods. We compared the performances of diffeomorphic and Bspline based registration methods. Finally, our method was validated using clinical data from 67 patients considering as target registration error (TRE) the distance between the estimated tumor position and the reference surgical position. The proposed method achieved a TRE of 16.21 ± 8.18 mm and it could potentially assist the surgery planning and guidance of breast cancer treatment in the clinical practice.


FIGURE 1.
Workflow of the proposed method for the registration of the preoperative image to the intraoperative surface to localize the tumor in the surgical position. The proposed scheme combines the information from the preoperative volume obtained from a prone MRI (Left) and the intraoperative surface obtained with an optical scanner (Right). Using the nipple as fiducial marker, first artificial surface fiducials are generated and then used in the subsequent landmark driven non-rigid registration step. The proposed workflow estimates the large volumetric transformation and localizes the tumor in the surgical position.
female breast cancer cases, representing almost one in four of all cancers in women. This disease is the most frequently diagnosed cancer in the vast majority of countries and is also one of the leading causes of cancer death [1]. Screening programs and imaging improvements have increased the detection of clinically occult non-palpable breast lesions that require preoperative localization [2], [3]. In women with non-palpable breast cancer, several randomized trials have shown that breast conservative surgery (BCS) is the treatment of choice [4]. The main challenge in the resection of nonpalpable tumors is to obtain appropriate surgical margins by minimizing the resection of healthy breast tissue with good aesthetic results. Currently, wire-guided localization (WGL) is the most widely used method for locating non-palpable breast lesions [5]. The limitations of WGL include technical complications such as wire transection and migration, patient discomfort, and poor cosmetic outcome [6], [7], [8]. The female breast undergoes different deformations due to external forces or changes in position during routine medical imaging. Magnetic Resonance Imaging (MRI) acquisitions in the prone position provide an accurate delineation of the size and extent of the tumor and offer the highest sensitivity for intraductal extension typical of invasive cancers. However, while prone breast MRI is well suited for diagnosis, it is less suitable for guiding BCS where the patient typically lies in a supine position, with the arm extended. Therefore, breast lesions undergo significant displacement between the preoperative image and the intraoperative positions. Methods using biomechanical models for the prone to supine pose transformation have been suggested for use in image-guided breast surgery [9], [11], [12], [13], [14]. However, the accuracy of alignment of biomechanical models alone is limited due to the lack of knowledge of the boundary constraints and loading conditions, which can contribute to errors in the location of lesions during clinical procedures. It is also important to note that the validation of these methods has been limited to a few cases.
Several groups have investigated the acquisition of an additional preoperative MRI acquired in the supine position [10], [15], [16], [17], [19] despite its disadvantages: longer protocol, breathing artifacts, limited diagnostic information, low contrast, cost, and the requirement, in any case, of a registration step to the final surgical position.
Prone to supine MRI registration based on pure intensitybased algorithms is likely to fail because the initial overlap is too small to guide the registration in the correct direction. For this reason, hybrid registration approaches based on biomechanical models and intensity-based image registration techniques have been widely used to map motion between diagnostic prone images and supine MRI [10], [15], [16], [17]. Preoperative supine MRI rigidly registered to an intraoperative optical scan to adjust the MRI to the exact breast position in the operative room has also been used for surgical guidance [19]. Although MRI-based supine approaches can benefit from comparable imaging sensitivity or a initial better alignment to the surgical presentation, is an addition to clinical practice. Furthermore,recent studies have shown that there is a significant deformation (including nonrigid) between surgical position and supine magnetic resonance positioning that need to be accounted for through image guidance [18]. Other groups proposed the use of the intraoperative surface acquired with an optical scan and to deform the preoperative MRI derived volume to match the surgical surface using finite element method (FEM) models with boundary conditions relying on skin-attached or anatomical markers [20], [21]. The patient surface extracted from a preoperative computerized tomography (CT) image in supine position has been proposed as a reference image in the surgical position to guide surface-based registration algorithms in conjunction with FEM-based biomechanical models [22], [23].
This work presents a workflow for performing a prone image-to-surgical physical data alignment strategy to determine the correspondence between the tumor identified in the preoperative MRI and the final position of the tumor in the surgical position. The developed method is based only on the prone MRI volume and the surface of the patient in the surgical orientation. The tumor localization error has been evaluated using retrospective cases with preoperative prone MRI and supine surface and tumor position extracted from a CT image. The position of the tumor has been estimated using a workflow based on non-rigid registration algorithms driven by artificial landmarks. This approach has been tested in 67 clinical cases, which is one of the largest study cohorts ever reported in the literature for prone to supine breast image registration.

II. METHOD
This work proposes a novel method for performing a volumeto-surface registration for tumor localization in breast surgery planning. The proposed workflow is explained in Figure 1. Given the preoperative volume obtained from a prone T2 SPAIR (SPectral Attenuated Inversion Recovery) MRI ( Fig. 1 Left) and the intraoperative surface obtained with an optical scanner ( Fig. 1 Right), the aim is to estimate the large volumetric transformation and localize the tumor in the surgical position. The proposed method involves a preprocessing step to obtain an initial rigid alignment between the volume and the surface and to generate artificial surface fiducials for the subsequent landmark driven non-rigid registration.

A. PREPROCESSING
The volume of the breast was obtained from the segmentation of the preoperative MRI. However, the preoperative MR scan was a T2 SPAIR MRI obtained with a fat suppression technique. This sequence is used in the breast diagnostic imaging protocol for its capability to enhance suspicious  lesions however, it is not the most accurate to differentiate breast contours especially the interior boundary of the breast. Therefore, we adopted a semi-automatic segmentation strategy. First, an automatic segmentation step (including thresholding and morphological filters) was applied, obtaining a breast volume mask with a precise outer boundary. Then, to correct the delineation of the inner boundary, we manually segmented the internal boundary of the breast in 5 slices of the volume and we interpolated them to obtain an internal volume mask. The latter was subtracted from the automatic segmentation to obtain the final preoperative volume mask, as shown in Figure 2c. The segmented preoperative volume and the intraoperative surface are first translationally aligned in the cranial-caudal direction assuming that the nipples in both positions lie approximately at the same level. Furthermore, to VOLUME 10, 2022 account for different tilts of the thorax, the intermammary space is aligned using a point based rigid registration algorithm. In this procedure, it is assumed that the intermammary cleft does not vary significantly between the two positions because the overall movement of the breast is generally in the posterior-lateral direction. After the initial alignment, the information from the optical scanned surface and the internal delineation of the breast volume in the prone position are combined to obtain an artificial volume mask of the breast in the surgical position, as shown in Figure 3d. The prone breast volume mask and the computed surgical volume mask are forced to have the same total volume.

B. ARTIFICIAL SURFACE FIDUCIALS GENERATION
To perform the landmark driven non-rigid registration between the prone and artificial supine volumes, it is necessary to define a set of surface correspondences to guide the registration process.
The nipples are the only known markers available and can be easily identified in the prone MRI and retrieved from the texture map of the optical scanned intraoperative surface. Taking the nipple position as reference, we have designed different strategies to define the set of corresponding fiducials points: • Nipple Areola Fiducials Set: a set of radial points that arise from the nipple and cover the nipple areola, as shown in Figure 4a.
• Breast Axes Fiducials Set: a set of points representing the horizontal and vertical axes of the breast, as shown in Figure 4b.
• Laplacian Surface Fiducials Set: a set of dense points that cover the entire surface. Their correspondence is based on a Laplacian transformation guided by the Breast Axes Fiducials Set [24] as shown in Figure 4c.
In order to define the Breast Axes Fiducials Set it was necessary to find corresponding boundary points of the breast volume in the medial, lateral, lower, and upper direction, as shown in Figure 5. Given the prone breast volume and the artificial supine volume masks, the main idea is to find a subvolume corresponding to the pendulous portion of the prone breast and the corresponding portion in the supine position. The portion of the prone volume relating to the pendulous breast is obtained by considering the part of the volume anterior to the coronal plane passing through the center of mass of the intermammary space. The corresponding supine breast volume mask is obtained by first removing the volume below an oblique plane passing through the medial point identified for the prone breast and the extreme lower point of the supine surface. Subsequently, the final match of volume with respect to the pendulous breast is achieved subtracting the required portion from the lateral part of the breast following the posterior to anterior direction (coronal plane).
The horizontal breast axes boundary points are found by intersecting the breast surface, corresponding to the pendulous breast portion, with an axial plane passing through the nipple. While for defining the vertical breast axes boundary points we considered the intersection between the breast surface and a vertical plane passing through the normal at the breast surface nipple point. To define the intermediate points between the upper-lower and medial-lateral boundary points, we first fixed the total number of points in the two directions and then we used geodesic distance to find equidistant points on the surface along those directions. Furthermore, we define an Internal Boundary Fiducials Set, which consists of a uniform set of points located on the internal surface of the prone breast corresponding to the boundary surface of the chest wall (Figure 4d).

C. NON-RIGID REGISTRATION
The main idea underlying this registration step is to smoothly propagate the transformation associated with the surface data onto all the positions in the breast domain in order to be evaluated at any desired set of positions.
Free-form deformation models have been successfully employed in non-rigid registration of medical images to handle situations in which organs have large, nonlinear deformations.
We adopted a landmark based version of the large deformation diffeomorphic metric mapping (LDDMM) registration [25]. This framework is specifically designed to deal with a large amount of deformation without distorting the topology of the object.
The LDDMM framework computes a transformation ϕ 1 between the prone image (moving) I P and the supine image (fixed) I S such that I S = I P • ϕ −1 1 . The computed transformation ϕ is the end point of a flow of velocity fields given by the ordinary differential equation v t (ϕ t ) = ϕ t , where ϕ 0 is the identity transformation. We have adopted a symmetric approach [26] where the optimal transformation is calculated by integrating the vector field that is found by minimizing the following equation.
where σ I represents the image noise variance, Sim(A, B) the similarity measure between the images and the norms are taken in the spaces of square integrable functions L 2 and allowed fields V. Then, by enforcing a certain smoothness on V, ϕ t is guaranteed to lie in the space of diffeomorphisms. In practice, smoothness is promoted by an operator L = (α + id) γ with the Laplacian and α and γ controlling the level and properties of smoothness respectively. The fourth term of (1) penalizes the divergence of the velocity field and induces better Jacobians. The term E S in (1) represents a symmetric spring energy formulation [27] that penalizes the distance between the landmarks in the prone and supine images and is calculated as following.
where S is the number of landmarks, w are weighting factors and x i P and x i S are the landmark positions in the prone and supine images. We adopted as similarity measure the sum of squared differences (SSD), and the registration algorithm was optimized by Polak-Ribière conjugate gradient conjugate gradient methods.
The performance of LDDMM was compared with a landmark driven Bspline registration [28], [29]. The optimization criterion is the one proposed by [28] and is composed of three different terms, a data term based on a similarity measure, a regularization term, and a spring term. Given a deformation function g, the cost function can be written as follows.
where Sim (A, B) is the similarity measure between the images, α is a parameter to weight the smoothness, S is the number of landmarks, w are weighting factors and x i P and x i S are the landmark positions in the prone and supine images.
The similarity measure was also SSD and the registration algorithm was optimized by a gradient descent method, embedded in a multi-resolution scheme.
Once the non-rigid registration is performed, the original preoperative breast volume mask is transformed using the deformation field to obtain the final deformed breast mask in the supine position and the estimated position of the tumor.

A. DATA
We collected 72 retrospective cases of breast cancer from the Gregorio Marañón General University Hospital (HGUGM) in Madrid after ethical committee approval. Each case included a preoperative T2 SPAIR MRI and a Subtraction MRI in the prone position, acquired several weeks before breast conserving surgery, and a CT image in the supine position used to simulate the intraoperative surface and in which the tumor in the supine position can be located for evaluation purposes. It is important to note that the acquisition of a CT image does not belong to the standard preoperative protocol, but is often required for cancer staging purposes. Two cases were excluded from the collection of cases, as the breasts in the MRI were in contact with the MR coil array. As a result, the deformation of the breast is caused not only by gravity, but also by the coil's holding force. One case was excluded because there was no clearly identifiable cancer lesion in the CT image. Other two cases were excluded because the breast in the CT images was cropped in the axial plane and could not be used to extract the surface in the supine position. Therefore, a total of 67 image sets (both MRI and CT) were used in this study. Table 1

B. DATA PREPARATION
The breast was segmented in both CT and MRI using 3D Slicer [30]. The segmentation of the breast in the CT image was only used to extract the surface used to simulate the patient's surface in the surgical position. The CT image was also used to segment the tumor in the supine position to validate the methodology. The volume of the breast in the prone position was obtained from the T2 SPAIR MRI in a semi-automatic manner as described in Section II-A. The tumor in the prone preoperative position was segmented in the Subtraction MRI. The position of the nipples in both the CT and MRI was marked manually in 3D Slicer to be used as reference real correspondence between the prone and supine surfaces. The prone and supine surfaces were rigidly registered aligning the nipples and the intermammary space (see Section II-A). After a visual inspection of the rigid registration process for 4 cases was necessary to correct the adjustment in the z direction to achieve a better alignment of the total breast. Furthermore, in 2 cases with large breasts, a correction in x-y plane was also necessary because the intermammary cleft was not matched in the supine position due to the movement of the breast in the medial direction.

C. IMPLEMENTATION AND EVALUATION
After the reconstruction of the supine volume masks and the generation of the corresponding surface fiducials, a landmark driven non-rigid registration was performed to estimate the prone to supine volumetric displacement field and to localize the tumor in the surgical position. We compared two non-rigid registration methods, one based on LDDMM and one based on a multiresolution Bspline. The landmark driven LDDMM registration algorithm was implemented in MATLAB, the maximum number of iterations were 300 and the average running time was less than 2 minutes for each case. A low resolution grid (8 mm) was considered to ensure robust performance and computational efficiency. The multiresolution landmark driven Bspline registration was implemented in Python 2.6 and a compiled C library to accelerate multidimensional Bspline computation. The maximum numbers of iterations in each stage were 150 and the average running time was more than 30 minutes for each case. None of the algorithms has been optimized to reduce computation time. The evaluation of the registration algorithm has been carried out using as target registration error (TRE) the Euclidean distance between the estimated position of the tumor centroid and the reference position of the tumor centroid in the supine CT image. The cutaneous projection of the estimated tumor position, defined as the closest point between the center of mass of the tumor and the skin, was also calculated and compared with the cutaneous projection of the center of mass of the supine tumor from the CT image. The distance between the two projected points is a key metric to evaluate the effectiveness of the method in supporting the decision making in the preoperative planning and it has been previously used to assess this type of algorithms [19], [20]. In addition, the Dice coefficient of the artificial supine mask and the deformed prone mask was computed to evaluate the global deformation.

D. PARAMETER ADJUSTMENT
To adjust the registration parameters and obtain the most satisfactory results we have considered a set of 5 cases with different characteristics. Multiple experiments have been conducted varying the weights of the landmark and regularization terms, the density of the deformation field and other parameters involved in the registration algorithms. After tuning the registration parameters, we used a randomly sampled subset of cases (training set) to find the set of fiducials that better guided the landmark driven registration process. We considered a training set of 45 cases with the following distribution of breast sizes: 22 % small, 38 % medium and 40 % large. The remaining 22 cases were used as an independent test set to evaluate the method and had the following distribution of breast size: 32 % small, 41 % medium and 27% large. Table 2 shows the average tumor distance obtained with different sets of fiducials for all the training set and for each breast size subset.   training set are considered, and specifically for M size breast group, there is a significant difference as proved by Wilcoxon left-tailed signed rank test (p ≤ 0.05). Additionally, the  workflow to obtain this fiducials set is simpler and it does not require the computation of the laplacian deformation of the surface like the Laplacian Surface fiducials set and therefore it ensures faster running time of the registration process. It is important to note that this fiducials configuration provides good performance for medium and large breast volumes for which the problem is more challenging due to large deformations. Table 2 also reports as None Fiducials Set the average tumor distance obtained when the registration is not landmark driven. The results are significantly worse as confirmed by paired Wilcoxon left-tailed signed rank test (p 0.05 ) and justify the use of the landmark driven approach.

E. RESULTS
The results of the performance of LDDMM registration are shown Table 3 in comparison with those obtained with the multi-resolution Bspline based registration. When all 67 patients were considered, the average tumor distance was 16.21 ± 8.18 mm for LDDMM and 18.24 ± 9.69 mm for Bspline while the Dice score was 0.93 ± 0.02 for LDDMM and 0.82 ± 0.06 for Bspline. The average distance between the cutaneous projections of the actual and estimated tumor was 13.86 ± 7.59 mm for LDDMM and 14.83 ± 8.12 mm for Bspline. Visual results of 3 cases from the training set are shown in Figure 7, from smaller (Fig. 7a) to larger tumor displacement (Fig. 7c) between prone to supine positions. Figure 8 shows the results of three cases of the test set with different breast volumes. Figure 9 shows visual results with a 3D model of the supine surface and the projections on the skin of the supine tumor and the estimated tumor.

IV. CONCLUSION AND DISCUSSION
The goal of this research is to develop a landmark driven non-rigid registration method to assist surgeons in planning, and possibly guiding, BCS by estimating tumor position. The method was developed to require only the surface of the patient in the surgical position and a standard of care diagnostic preoperative prone MRI in which the tumor is clearly visible. After rigid alignment, we generated an artificial supine volume based on the delineation of the chest wall boundary and then defined a set of artificial correspondences between the prone and the supine surfaces to guide the volumeto-volume landmark driven registration. We extended the LDDMM registration approach to include a landmark based spring term and compared the results to a multiresolution non-rigid Bspline registration that also includes such a term. Our experiments demonstrate that the LDDMM method is more accurate, robust and faster than the Bspline implementation used and it ensures more stability of the retrieved deformation field. The proposed method was evaluated in 67 clinical cases with an average localization error of 16.21 ± 8.18 mm. The results for the test set (13.88 ± 7.30 mm) reported in Table 3 are better than those reported for the training set. This can be attributed to the fact that as a result of the random split in this subset there is a smaller proportion of large breasts (24 % in the test set vs 40 % in the training test). As reported in Table 2 the method is more effective in the range of breast volumes including small and medium size breasts. No significant variation in the localization error was found according to tumor depth and size.
This preoperative localization approach could be useful to guide the surgeon towards the location of the tumor without the limitations of WGL. The benefit of use could be of interest even in combination with other guidance mechanisms such as the wire itself, radar guidance, radiofrequency identification, or radioactive seeds [31], [32]. On the other hand, the target of resection always includes a 5-10 mm margin of healthy tissue surrounding the tumor in order to allow for the inaccuracies inherent to the different localization techniques and to ensure complete removal of cancerous tissue.
In addition, the error associated with the skin projection of the lesion is smaller and may be the main interesting information for the procedure. The use of a visualization system in the operating room with a 3D model of the breast surface and the skin projection of the tumor could complement the WGL information and help determine the most appropriate incision site to ensure a better cosmetic outcome or to endorse the localization achieved by other less invasive localization techniques.
The tumor target error results obtained are very competitive considering that they are obtained without including biomechanical models in the registration workflow nor the requirement of acquiring supine MRI or other preoperative supine imaging. Some prior hybrid methods using supine image information achieve better results. However, they have been tested in a very reduced datasets (less than 5 cases) without ensuring generalization of the proposed technique [10], [15], [16], [17]. The use of surgical surface information combined with biomechanical models produced similar results (average landmarks distance of 15.3 mm) [22] or worse results (average cutaneous distance of 20.8 mm) [20]. Recent research using a hybrid approach and combining biomechanical models and supine surface information has obtained very competitive results (average tumor distance of 8.05 mm), although the dataset is still limited (25 cases) and with a smaller percentage of medium and large breast volume cases [23].
The validation of the current proposal has been developed using the surface of the patient from a supine CT image for validation purposes. The proposed workflow has been designed to use an optical scanner for the acquisition of the breast surface in the operating room, altering minimally the surgical procedure. Further research is guaranteed in a prospective study during surgery or biopsy using the optical scanner. Preliminary experiments of this setup have confirmed the feasibility of acquisition and processing of the surface from the optical scan of the supine position, including the proper identification of the nipple from the texture map. Unlike previous proposals in the literature that involve supine MRI, our proposal does not require any change in the standard of care imaging protocol.
Future research will consider improving tumor localization results by incorporating deep learning algorithms to accelerate and learn the plausible deformation fields underlying the transformation from prone to supine breast pose.
In conclusion, we have designed and developed a tool that using only the surface of the patient in the surgical position and a preoperative MRI achieves competitive results in terms of localization error and could potentially be used in the clinical practice to assist the planning and guidance of breast cancer treatment. Additionally, this study cohort represents the largest reported for prone to supine registration, ensuring that different ages and sizes are represented.  (2011-2012 and 2013-2014). He is currently an Associate Professor at Universidad Complutense de Madrid, Madrid. He is also a Researcher with the Instituto de Investigación Sanitaria Gregorio Marañón. He has authored or coauthored more than 50 articles in indexed journals and books on the multidisciplinary treatment of breast cancer pathology.
ANDRÉS SANTOS (Senior Member, IEEE) received the M.Eng. and Ph.D. degrees from the Universidad Politécnica de Madrid. He is currently a Professor and the Director of the Biomedical Image Technologies Laboratory, Universidad Politécnica de Madrid. He is an International Research Fellow at SRI International, Menlo Park, CA, USA. He has been the Co-Director of the master's program on biomedical technology and instrumentation. He has authored or coauthored over 100 articles in indexed journals on medical imaging. He has coordinated more than 50 national and international research projects and technology transfer contracts with private companies. His main research interest includes biomedical image analysis and acquisition. She additionally completed two different master's programs, such as Biomedical Engineering and Medical Physics from the University of Patras, Greece; and Bioengineering from UNED, Spain. She is currently a Full Professor with the Biomedical Image Technologies Laboratory, Universidad Politécnica de Madrid. She has authored or coauthored over 100 publications in indexed journals and conferences. Her main research interests include biomedical image analysis, especially cardiac and pulmonary imaging, image-guided therapy, microscopy image analysis, registration, and motion estimation and compensation. VOLUME 10, 2022