Email updates

Keep up to date with the latest news and content from BMC Medical Imaging and BioMed Central.

Open Access Research article

Multimodal image registration of the scoliotic torso for surgical planning

Rola Harmouche12*, Farida Cheriet12, Hubert Labelle2 and Jean Dansereau12

Author Affiliations

1 École Polytechnique de Montréal. 2500, chemin de Polytechnique, Montreal, H3T 1J4, Canada

2 Hôpital Ste-Justine 3175, Chemin de la Côte-Sainte-Catherine, Montreal, H3T 1C5, Canada

For all author emails, please log on.

BMC Medical Imaging 2013, 13:1  doi:10.1186/1471-2342-13-1


The electronic version of this article is the complete one and can be found online at: http://www.biomedcentral.com/1471-2342/13/1


Received:26 May 2012
Accepted:26 November 2012
Published:4 January 2013

© 2013 Harmouche et al.; licensee BioMed Central Ltd.

This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

Background

This paper presents a method that registers MRIs acquired in prone position, with surface topography (TP) and X-ray reconstructions acquired in standing position, in order to obtain a 3D representation of a human torso incorporating the external surface, bone structures, and soft tissues.

Methods

TP and X-ray data are registered using landmarks. Bone structures are used to register each MRI slice using an articulated model, and the soft tissue is confined to the volume delimited by the trunk and bone surfaces using a constrained thin-plate spline.

Results

The method is tested on 3 pre-surgical patients with scoliosis and shows a significant improvement, qualitatively and using the Dice similarity coefficient, in fitting the MRI into the standing patient model when compared to rigid and articulated model registration. The determinant of the Jacobian of the registration deformation shows higher variations in the deformation in areas closer to the surface of the torso.

Conclusions

The novel, resulting 3D full torso model can provide a more complete representation of patient geometry to be incorporated in surgical simulators under development that aim at predicting the effect of scoliosis surgery on the external appearance of the patient’s torso.

Background

Idiopathic scoliosis is a disease characterized by a complex three-dimensional curvature of the spine and the rib cage; these internal curvatures are externally manifest as a lateral trunk asymmetry and/or a rib hump. Such external deformations are often aesthetically undesirable for patients and can cause psychological problems, and in more severe cases, chronic back problems or pulmonary problems [1]. Treatments aim at slowing down the curvature progression or at correcting some of the undesired curvature. They include a brace in less severe cases and surgery in the form of vertebral fusion in more severe cases. Surgeons rely on their experience and intuition in order to establish the adequate instrumentation that would lead to the desirable post-operative external trunk appearance. However, the effects of the brace or surgical instrumentation on the external shape of the trunk cannot be reliably predicted prior to completion of treatment. Our group has recently developed a simulator with the aim of predicting the effect of scoliosis surgery on the spine and the torso [2]. Interesting preliminary results were obtained. However, the simulation outcome lacked generalizability across patients. Current research aims to integrate both bone and soft tissue in the 3D model of the trunk in order to verify whether surgical simulators that model soft tissue information could be useful in improving the prediction of the effects of surgery on the external appearance of the patient’s trunk. Such a model would require fusion of soft tissue information, typically obtained from MRIs in prone position, spine and rib cage data, typically obtained from X-rays in standing position, and a representation of the external surface of the trunk, obtained from an active vision system. This registration task is not trivial because of a mixture of both rigid and non rigid deformations that occurs between the acquisition of the different image modalities, in particular between the MRI and the remaining modalities, since the patient is lying down during MRI acquisition. Adding to the difficulty is the lack of correspondences between the tissue visible on the MRI and the structures visible on the remaining modalities. The aim of the present work is to register MRI data acquired in prone position with X-ray and TP data acquired in standing position, all while compensating for the postural changes that occur between the acquisitions and while preserving the rigidity of bone structures, in order to obtain a 3D representation (incorporating bone structures, trunk surface, and soft tissue information) of the torso of a patient with scoliosis.

Multimodal medical image registration has been applied to several types of images such as MRI/X-ray, MRI/CT, CT/X-ray, MRI/ultrasound, etc. Registration techniques can either be rigid, affine, semi-rigid, or deformable. Rigid or affine registration techniques apply either a rigid or an affine transformation to the entire source image being registered with a target image. Early work used rigid registration techniques in order to register MRI/CT and MRI/PET data of the brain [3]. In terms of non-rigid methods, several techiniques have been used in medical image registration. Those consist of thin-plate spline [4], free-form deformations using B-splines [5-8], elastic models [9], fluid models [10-13], and Markov random field approaches [14-16].

Most of the work on MRI/X-ray registration was not focused on the spine and consisted mainly of rigid registration methods [17-25] to the exception of [26] which use perspective transformations and [27], which use free-form deformations. A review of these techniques can be found in [28] and [29].

Rigid registration techniques cannot capture the complex deformation that occurs in the shape of the spine between the standing and prone positions in which the different image modalities are acquired. Furthermore, with the knowledge that the vertebrae themselves are rigid structures, traditional purely non-rigid registration algorithms are also not appropriate for the task at hand. For example, Skerl et al. [7] perform registration of vertebrae using MR and CT images, but do not model the rigidity of these structures. As a compromise, semi-rigid techniques have been used for the registration of spine data. For example, vertebral structures extracted from MRI data have previously been modeled as rigid bodies for registration purposes [30]. Soft tissue was registered using modified thin-plate splines that allow segmented vertebral structures to be constrained to rigid deformations. Similar work was done by Rohr [4], requiring only a few correspondence points instead of a full segmentation of the rigid structures. Huesman et al. [31] register CT and MR images of the torso using thin-plate spline transformations. Vertebral rigidity constraints were incorporated into the thin-plate spline approximation parameter. However, acquiring CT data for the entire torso is not possible in clinical settings due to radiation issues. More importantly, CT and MR images are both acquired in prone position, which implies that non-rigid deformations due to posture changes are not taken into account. Loeckx et al. [8] use B-splines in order to register PET/PET and CT angiography/CT angiography while incorporating a rigidity constraint in the cost function. Wang et al. [32] use triangular B-splines to deform 3 sagittal 2D MR images of vertebrae by placing knots at the boundaries of segmented vertebrae in order to insure their rigidity. B-spline functions have the advantage of allowing multiresolution registration. However, they require a regular grid of image data, which is not available in the image modalities that we wish to register. None of the works above register MRI and X-ray data, and, since all the images are acquired in prone position, nor do they model the spinal deformation between postures. Some works that have registered MRI to X-ray spine data have either required fiducials on cadaveric data [33] or required 9 X-ray images per patient for adequate X-ray/ MRI registration [34]; both strategies which are not possible in clinical settings. In later work, Van de Kraats et al. [35] generate CT-like data from MRI in order to perform rigid registration of vertebral bodies from X-ray images of cadavers. Markelj et al. [36] use 2 2D X-rays in order to rigidly register MRI and X-ray images by matching 3D gradients of 3D images to 3D gradients reconstructed from the 2D X-ray images. Tang et al. [37] apply rigid registration to fit a model of the vertebrae only onto Sagittal MRI slices in order to measure spinal curvature. Harmouche et al. [28] used an articulated model representation of the spine in order to register the spine extracted from MRI to that obtained from X-ray data of the same patient. The remaining tissues on the MRI were not registered. Finally, to our knowledge, no previous works address the simultaneous registration of MRI, X-ray and TP data.

This paper proposes a method that registers MRI to both TP and X-ray data in order to construct a 3D model covering the thoracic and lumbar vertebral levels of the torso of a scoliotic patient. MRI slices are acquired in prone position and X-ray and TP data are acquired in standing position. External markers are placed on the patient’s skin prior to X-ray and TP acquisition in order to register the two modalities using thin-plate splines. The MRI is then registered to the TP and X-ray data using a two step process: first, bone structures guide an articulated model representing the 3D spinal deformations that occur between the two postures and serves as an initial transformation of the MRI data. This transformation is then refined by confining the soft tissue to the volume delineated by the trunk surface and the surface of the vertebrae using weighted thin-plate splines with constraints for bone rigidity and soft tissue elasticity. The availability of the surface contours allows us to model the rigidity constraint in a simpler fashion when compared to previous works. To our knowledge, this is the first work that combines bone, soft tissue, and surface topography information for a model of the human torso.

This article is separated as follows: The Methods section describes the data acquisition process and the articulated model based registration. Qualitative and quantitative results of the proposed method are compared with both rigid and simple articulated model registration in the Results and discussion section. The key conclusions and proposed future works are presented in the Conclusions section.

Methods

In order to obtain the 3D representation of a patient incorporating bone, soft tissue and surface information, MRI, X-ray images and TP data are obtained for scoliosis patients. The 3D spine models are then extracted from the MRI and X-ray images and the 3D positions of external markers are obtained from the TP and X-ray data. The goal is to then align MRIs of the torso with reconstructed thoracic and lumbar vertebrae from biplanar X-ray images and surface topography data (TP). Thus, we are searching for the independent transformations TTPXray and TMRIXray, which would transform MRI and TP data into X-ray space, respectively.

Data acquisition and preprocessing

MRI, X-ray images and TP data are obtained for three adolescent patients with scoliosis at the Sainte-Justine Hospital in Montreal. Patient 1 is 153.4 cm tall, weighs 48.8 kg, and has a right thoracic curvature with a 50° Cobb angle. Patient 2 is 146.1 cm tall, weighs 41.7 kg, and has a left thoraco-lumbar curvature with a 57° Cobb angle. Patient 3 is 160.5 cm tall, weighs 46.2 kg, and has a right thoracic curvature with a 49° Cobb angle. The present research involving human subjects has been approved by the ethics committee of the CHU Sainte-Justine Research Center (Comité d’Éthique de la Recherche du CHU Sainte-Justine). Written consent is obtained from the participants and a parent or guardian. These 3 image modalities are obtained during the same day but at different times for each patient.

First, four 3D digitizers (Creaform inc., average resolution of 1.1 mm), each covering a different view of the patient, are used in order to obtain the surface topography of each patient’s torso (Figure 1). Each digitizer consists of a color CCD camera and a structured light projector. The patient is instructed to stand in the middle of the camera setup, legs slightly apart, arms slightly raised to the sides. This posture was found to minimize occlusions to the patient’s torso all while maximizing patient comfort (thus minimizing patient movement). Deformed patterns are then obtained from each digitizer and are used to retrieve shape and texture information. The shape information obtained from the different views is merged using the Inspeck EM software resulting in a mesh representation of the surface topography (Figure 2(b)). Prior to acquisition, several visible green adhesive markers are placed on the surface of the patient’s torso. The markers are placed at the following locations: vertebra prominens, pilonidal dimples, antero-posterior iliac spine, sternum, spinous process, and inferior angle of the shoulder blades. These markers are regularly used in our clinical setting since their localisation was found to be the most reproducible. The 3D position of the adhesive markers is then obtained by manually identifying the markers on the digitized texture information using in-house software.

thumbnailFigure 1. Inspeck system. Inspeck system used for capturing surface topography data. Four digitizers allow the capturing of the full surface of the patient torso.

thumbnailFigure 2. Medical image modalities. Multimodal image data used for the 3D patient representation developed. The spine representation in (c) is reconstructed from bi-planaer X-rays. The externally placed landmarks are displayed in blue on the surface topography and in red on the X-ray image. The vertebral landmarks used to construct the articulated model are displayed in green on the X-ray (c).

Postero-anterior and lateral radiographs are then obtained for each of the patients while they are in standing position, legs slightly appart, and arms to the front and bent upwards. An explicit calibration method [38] is then used in order to calculate the 3D position of vertebral landmarks manually identified by experts on both X-rays. The obtained landmarks allow the mapping of generic vertebral surface meshes onto the patient space; resulting in a patient-specific 3D geometric model of the spine from the X-ray data (Figure 2(a)). In this case, six of those landmarks are used per vertebra in order to generate the X-ray articulated model. These are placed on the centers of the superior and inferior plates of each of the 17 thoracic and lumbar the vertebral bodies, and below and above the left and right pedicles. The chosen landmarks are regularly used for vertebral reconstruction in our clinical setting since they are clearly visible on the X-rays. In addition, prior to image acquisition, radio-opaque markers which are clearly identifiable on the obtained X-rays are placed on the surface of the patient’s torso and affixed on top of the adhesive markers used for TP data acquisition. The 3D position of the markers is then calculated using the same methodology as for the vertebral landmarks. These external markers are used for TP/X-ray registration.

Sagittal MR images of the spine are routinely acquired in our clinical setting for scoliosis patients awaiting surgery (1.5 Tesla, TR/TE = 771/15, 704x704, 350 FOV, with a 0.5 mm by 0.5 mm in-plane resolution and 3 mm thickness, with a 3.6 mm separation between slices). In order to visualize the entire torso, T1-weighted axial MRI slices of 1 mm by 1 mm in-plane resolution, 2 mm thickness and 12 mm spacing between slices and covering the entire torso are acquired (Siemens Symphony system 1.5 Tesla, TR/TE = 650/12, 704x704) (Figure 2(c)). Since acquisition of the sagittal slices already requires a total of 30 minutes inside the MRI, the larger spacing between the axial slices allows for the entire trunk to be imaged in a reasonable amount of time for patient comfort. The time inside the MRI thus totals around 50 minutes including rest time between acquisitions. Since the sagittal and axial slices are acquired during the same scan, their location with respect to each other is known and they are thus pre-registered. The sagittal slices are thus used in order to situate the axial MRI acquisitions to the appropriate vertebral level and consequently to calculate the transformations required. This is required because the distance between the axial slices is large and the contours of the vertebrae are not all visible on these acquisitions. It must be noted that in cases where sagittal slices are not required for clinical purposes, we would be able to use the 50 minute total scan time in order to acquire the axial slices only, thus allowing for a significantly better resolution and better visibility of the vertebrae. In order to generate a 3D geometric model of the spine from MRI data, the 3D shape of the seventeen thoracic and lumbar vertebrae is manually segmented from the sagittal slices using TomoVision’s SliceOmatic software. Landmarks are then manually placed on the left and right edges of the posterior, anterior, inferior and superior ends of the vertebral body for all thoracic and lumbar vertebrae also using TomoVision’s SliceOmatic software. These landmarks are used in order to calculate the articulated model for MRI/X-ray registration.

X-ray / TP registration

The TP data is registered to the X-ray reconstruction by applying a thin-plate spline transform (TTPXray) [39]. The transform has the following form:

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M1','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M1">View MathML</a>

(1)

where f(x,y,z) is a vector valued function which maps each point (x,y,z) of the original TP data onto a point (x,y,z) in X-ray space. The ai are the coefficients of the affine transformation, wiare the weights, and r is the distance between the control point and the point to be transformed. U(r)=|r|for a 3D transformation, as was stated in the works of Bookstein [39]. The ai and wi are estimated by solving the following system of linear equations using the manually placed markers on the surface of the patient prior to acquisition of both TP and X-ray data:

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M2','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M2">View MathML</a>

(2)

where P is the vector of source points obtained from the TP data, V is the vector of target points obtained from the X-ray data, and K is the matrix containing U(r).

Transforming MRI data

Since the spacing between consecutive axial MRI slices is large, they are treated as independent 2D slices and the continuity between successive slices is not adressed. Thus, each slice is registered separately to the remaining model using TMRIXray, which a composition of an articulated-model transformation (TMRIbone) and a thin-plate spline transformation (TMRIsoft). The goal is to transform each MRI voxel into the space of the 3D X-ray model, taking into account the non-rigid deformation following the posture change subject to the following constraints: First, the spine extracted from the MR images has to be aligned with the X-ray spine model. Second, MRI data of the torso has to be contained within the TP volume, such that the contour of the torso on the MRI corresponds to the surface topography.

In order to register each MRI slice such that the spine information extracted from MRI data and that extracted from X-ray data are aligned, an articulated model previously used to align 3D models of the spine obtained from both modalities is calculated [28]. This articulated model allows us to define the spine as a combination of local intervertebral transformations which can be obtained in a number of ways. In our case, we segment the vertebrae from the sagittal MRI data in order to extract the vertebral landmarks. The sagittal slices are used since they are readily available in our clinical setting and have a better resolution than the axial slices. The local coordinate system of each of the thoracic and lumbar MRI vertebrae v is then obtained by calculating the center and the orientation of the vertebrae. Similarly, the local coordinate system of each of the thoracic and lumbar X-ray vertebrae is obtained using the X-ray vertebral landmarks. Intervertebral transformations consisting of rotations and translations are then defined between the local coordinate system of each vertebra v and its lower neighboring vertebra v−1. Although landmarks are used in order to obtain the local coordinate system of each vertebra, the landmarks extracted are different for the X-ray and MRI data and thus point correspondences cannot be established. Thus, they have not been explicitly used for registration. In fact, future work consists of obtaining the local coordinate system from automatic vertebral segmentations, which can replace the need for manually placed landmarks. The global transformation for each of the vertebrae is calculated separately using a composition of the local rigid intervertebral transformations.

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M3','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M3">View MathML</a>

(3)

where T0,1 is the global rigid transformation between the world coordinates and the first vertebra on each of the two image modalities. The overall articulated model transformation (TMRIbone)from MRI to X-ray space is then obtained using a concatenation of the global MRI and X-ray vertebral transformations:

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M4','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M4">View MathML</a>

(4)

Once the transformation model between the MRI data acquired in prone position and X-ray spine acquired in standing position is obtained, it can be used to provide a preliminary registration for the axial MRI slices. Since both the sagittal and axial MRI data are acquired at once, they are already registered to the same space and thus can be used interchangeably. Each axial MRI slice s is transformed using the vertebral transformation TvsMRIbone where vs is the vertebra that has the closest z value to the slice. This part of the transformation takes into account the changes in the general alignment of the MRI slices due to changes in posture. This is followed by a weighted thin-plate spline transformation TMRIsoft in order to approximate the non-rigid deformations on each slice separately. This transformation is illustrated in Figure 3 and is calculated in the following manner: first, the plane corresponding to each of the registered MRI slices (i.e. having the same location and angle) is obtained from the surface topography data. In order to guide the non-rigid registration between the MRI and TP planes, correspondence points are extracted from the torso contour of the MRI slice (visible in red in Figure 3 ) and from the corresponding plane on the surface topography (visible in green in Figure 3). These points are extracted at 30 degree intervals, angle 0 being a vector passing from the center through the anterior point of the vertebra. The MRI slices have already been transformed to the appropriate location and angle with respect to the vertebrae following the articulated model registration. Thus, points in the same direction on a slice can be used as correspondences. These correspondences are used given the fact that we do not have access to anatomically significant landmarks and that the purpose of the correspondence points is to drive the thin-plate spline registration by ensuring that the interior of the surface topography is fully covered by the MRI slice. A preliminary thin-plate spline transformation TMRIsoftedge is approximated by solving the same system of equations as in equation 2. The points extracted from the edges of the segmented MRI are used as source landmarks and those obtained from the corresponding plane on the surface topography are used as target landmarks. The transformation of an arbitrary point p is then a modified version of the thin-plate spline transformation TMRIsoftedge: A weight is incorporated depending on the relative distance of the point p from the edge of the torso and from the closest vertebra, given that this transformation should be equal to the identity matrix I within the vertebra. The overall transformation becomes:

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M5','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M5">View MathML</a>

(5)

where, for each point p:

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M6','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M6">View MathML</a>

(6)

thumbnailFigure 3. Ratio for weighted thin-plate spline. Calculation of the weighted thin-plate spline ratio based on the distance of the point p to be registed, the closest vertebra, and the surface of the torso. The points in green are obtained from the contours of the surface topography and the points in red are obtained from the MRI.

where Dvertebrais the distance between point p and the border of its closest vertebra obtained from the MRI data and Dsurface is the distance between point p and the border of the segmented MRI torso. Since the meshes representing the X-ray vertebrae were originally obtained from cadaveric data and registered to patient space using a few interest points, the use of the MRI vertebral meshes were favoured in order to calculate this distance. The MRI vertebrae were segmented from real patient data. The weighting applied to the thin-plate spline at each voxel is such that the further we are from the vertebra, the more weight is given to the non-rigid deformation. The proposed transformation was chosen as a simple and sufficient way to model the deformation between the vertebrae and the torso surface, as there is a lack of correspondences between the different image modalities in that area. In addition, the present goal is to simply fill the volume contained within the surface topography with approximate soft tissue information that would allow for the creation of a 3D patient model for surgical simulations.

Results and discussion

Following the desciption of the validation methods, qualitative and quantitative results for our registration method will be presented.

Validation

The proposed registration method is compared to rigid registration and simple articulated registration. Qualitative results are followed by overlap comparisons between the surface contained within the torso on the TP and the registered MRI slices using the Dice similarity coefficient [40]. In order to calculate the Dice measure, a mask of the torso data is obtained for each axial MRI slice. A corresponding planar cut through the surface topography at the same position in 3D patient space and in the same orientation as the registered MRI slice is obtained and a second mask of the interior of the surface topography is then calculated. The overlap between the two masks is then calculated using the following equation:

<a onClick="popup('http://www.biomedcentral.com/1471-2342/13/1/mathml/M7','MathML',630,470);return false;" target="_blank" href="http://www.biomedcentral.com/1471-2342/13/1/mathml/M7">View MathML</a>

(7)

Higher Dice values signify better overlap between the 2 surfaces, with perfect overlap corresponding to a Dice value of 1. The assessment of Dice values is somewhat dependent on the size of the areas being compared. That is, larger areas have higher area to perimeter ratio resulting in inherently higher Dice values. Since in the present case the different methods are comparing the same area being deformed, the Dice value is adequate for the comparison. In addition, the Dice values calculated assume that the MRI slices are in the proper orientation following registration and correspond anatomically to the surface topography slice at the same location and orientation. We make this assumption as the MRI slices have been aligned using the articulated model. We tolerate any bias obtained following this assumption given the fact that our main aim is to verify whether the soft tissues fit within the space encompassed within the surface topography. The Dice measure is able to verify this aim. The mean and variance of the Dice values is calculated for all patients, in addition to testing for statistically significant differences between the Dice values of the different methods.

In order to provide a more precise Dice measure that does not depend on the orientation of the slices, a 3D Dice analysis is preformed. This analysis reduces the bias present in the 2D Dice calculations. In order to obtain a volume from the registered MRI slices, we first build a surface mesh over the entire torso by triangulating points sampled from the MRI contours using a Delaunay triangulation. Since the MRI slices have a large distance between them, a more detailed mesh is then obtained using linear subdivision. We then create a volume with a resolution of 1 mm x 1 mm x 3 mm and having slices parallel to the x-y plane, and set all voxels inside the surface to be equal to 1, otherwise 0. This resolution is chosen since it is only slightly larger than the MRI slice thickness, and such a voxel size has been successfully used for MRI segmentation applications in the past. A similar volume with the same resolution is obtained from the surface topography mesh. We then compare every voxel of the obtained MRI and TP volumes using the Dice measure.

In order to assess the quality of the deformation field resulting from our proposed registration, and to examine the areas that display the highest amount of deformation, the determinant of the Jacobian of the deformation field is calculated both with and without the weighting applied on the thin-plate spline transformation. The Jacobian of the deformation field J is defined as the matrix of partial derivatives of the deformation field. It has been used in previous works on image registration in order to ensure a smooth and topology preserving transformation [41]. For each voxel, the transformation vectors are calculated using the methods described previously. The discrete partial derivatives making up the Jacobian are then approximated using the first order derivative of a Gaussian (sigma = 1). The determinant of the Jacobian is obtained from the partial derivatives and indicates the voxel-wise relative volume change. For example, a determinant value of 3 indicates that the original volume has decreased by a factor of 3. In order to insure volume preservation of the tissues and organs during registration, a determinant value of 1 is required. The determinant of the Jacobian for the deformation of all voxels of MRI slices are shown and areas of highest deformation are addressed. The mean and variance values for all slices are then calculated.

Qualitative results

Figure 4 shows results for 1 slice of patient 1. A clear misalignment can be seen between the vertebrae extracted from the X-ray data and those extracted from MRI data when rigid registration is used (Figure 4(a)), illustrating the difference in the shape of the spine between the two postures in which these two modalities are acquired. As a result, the axial MRI slice displayed is misaligned to the surface topography and the X-ray vertebrae. The articulated model provides a better vertebral alignment (Figure 4(b)), but does not provide an adequate fitting between the MRI soft tissue data and the surface of the trunk. The torso on the MRI slice shows better alignment with the surface topography compared to both rigid and articulated-model registration (Figure 4(c)). The proposed method is thus able to correct for some of the deformations that can be visible on the MRI due to the change in posture between the image modalities.

thumbnailFigure 4. Registration results. Results comparing rigid registration with registration using the proposed method for the MRI, X-ray and TP. The X-ray vertebrae are displayed in pink and the MRI vertebrae are displayed in cyan. A decimated mesh representing the surface topography is displayed in blue. It can be seen the results with the proposed method yield a better alignment of the MRI and X-ray vertebrae and a better fitting of the MRI within the surface topography.

Quantitative results

Dice anaylsis

Figure 5(a) presents a graph of the Dice values for all 33 axial slices of patient 1 covering the thoracic and lumbar vertebral levels. An average Dice value of 0.973±0.008 is observed for the proposed method compared to 0.806±0.064and 0.826±0.071 for the rigid and articulated registrations, respectively. A statistically significant improvement is observed when the proposed method is used compared to both the rigid and the articulated model registration (p < 0.01 in both cases). A statistically significant improvement is also observed when the articulated model registration is compared to rigid registration (p<0.01). The articulated registration performs worst on slices containing breast tissue (slices 15 onwards), where our method is better able to approximate the non-rigid deformations of the breasts. Our proposed results still contain inaccuracies in the breast area as their deformation is somewhat independent from the rest of the tissues. However, proper modeling of breast deformation was not addressed in the present work as it does not affect the accuracy of the outcome of a scoliosis surgery.

thumbnailFigure 5. Quantitative registration results. Dice values for patients 1,2 and 3 measuring the overlap between torso data obtained from axial MRI slices and corresponding surface topography axial cuts for all patient slices for rigid (red), articulated (green) and the proposed (blue) registration methods. Results range from inferior to superior axial slices.

Figures 5(b) and 5(c) show the Dice values for patients 2 and 3, respectively. A significant improvement in Dice values is also observed in these cases as well when our method is compared with rigid and articulated registrations. For the case of patient 2, the average Dice values are 0.844±0.028, 0.830±0.031 and 0.974±0.008 for rigid, articulated, and proposed registrations, respectively. A statistically significant improvement is also observed in the case of patient 2 when the proposed method is used compared to both the rigid and the articulated model registration (p<0.01 in both cases). However, an improvement was not seen when rigid and articulated model Dice values were compared.

In the case of patient 3, the average Dice values obtained are 0.878 ± 0.052, 0.871 ± 0.010, and 0.976 ± 0.010 for rigid, articulated, and proposed registrations, respectively. As in the case of patient 2, a statistically significant improvement is observed when the proposed method is used compared to both the rigid and the articulated model registrations (p < 0.01 in both cases), and no improvement is observed when rigid and articulated model registrations are compared.

The 3D Dice measure is obtained for the 3 patients using the results from the proposed method. In this case, values of 0.939, 0.938, and 0.940 were obtained for patients 1, 2, and 3, respectively. These values are only slightly lower (by an average of 3.66%) than those obtained using the 2D Dice measure. This is to be expected since, unlike the 2D Dice measure, the 3D measure compares voxels within the patient’s surface that were not used in the registration process. More importantly, these errors are inherent to the fact that an interpolation, more so a linear interpolation, was used between consecutive MRI slices due to the low MRI resolution and not to the method itself. Nonetheless, the values obtained using the 3D measure are still considered excellent in the literature, and are still considerably higher (by an average of 11.52%) than the 2D values obtained for rigid and articulated registration. We can thus conclude that, although the spacing between consecutive MRI slices is high, our method still provides a volume of soft tissue information contained within the surface topography.

One obvious caveat with the Dice similarity coefficient to validate our work is that it does not measure the anatomical correctness with which our soft tissues were deformed. However, for the desired precision, given that the MRI resolution is low in the z-direction thus leading to a less precise overall model, that a measurable difference in trunk shape changes due to treatment is in the order of a tenth of a cm, and that no anatomical correspondences are present in the space contained between the vertebrae and the patient’s trunk surface, the current results are deemed satisfactory. Our main concern is to obtain a model with which the vertebrae maintain their rigid characteristics and the space contained within the surface of the patient is filled with soft tissue information. However, we must note that the presented registration framework does not take the physical characteristics of the tissues into account when modeling the deformations.

Assessment of the transformation

The quality of the transformation is studied using the determinant of the Jacobian. Figure 6 shows the determinant of the Jacobian of the deformation field overlayed on top of the MRI data for slice 9 of patient 3, which is situated around the stomach area. Results without (a) and with (b) the rigidity constraint are displayed. The spectral mapping shows values closer to 1 as blue and values further from 1 as red. Since we did not impose any requirements on volume preservation, the fact that that most values are far from 1 in (b) is not problematic. However, it is interesting to see which areas display the most significant amounts of variation in the deformation. Both cases have higher Jacobian values towards the anterior portion of the torso (left of the images). This is because there has been an increased amount of volume expansion around the stomach following registration, which was intended to counter tissue compression that occurrs while the patient is lying down during MRI acquisition. The main difference between the values with and without the rigidity constraint is that, in the case of our proposed method, the determinant is equal to one inside the vertebra (due to the rigidity constraint imposed), and increases as voxels go radially outwards towards the surface of the patient’s skin, since a higher weight is placed on the non-linear portion of the deformation.

thumbnailFigure 6. Determinant of the Jacobian for slice 9. Determinant of the Jacobian of the deformation field resulting from registration without (a) and with (b) rigidity constraints for slice 9 of patient 3. A spectral intensity mapping is displaying |Jacobian(i)−1|for pixel i superimposed on top of the registered MRI slice. The spectral mapping shows values closer to 1 as blue and values further from 1 as red. The anterior portion of the torso is to the left of the image.

Figure 7 shows the determinant of the Jacobian of the deformation field for slice 25 of patient 3, which is situated around the breast and shoulder area. In each of these images, the red and yellow areas are located in the breast region. In this case, deformation varies much less smoothly around the breasts than in the remaining areas. These results are expected, as the breasts are likely to deform independently from the articulated model used.

thumbnailFigure 7. Determinant of the Jacobian for slice 25. Determinant of the Jacobian of the deformation field resulting from registration without (a) and with (b) rigidity constraints for slice 25 of patient 3. The intensity mapping is the same as in the previous figure.

Discussion

The results above show that the proposed method is able to register the MRI, X-ray and TP data of a human torso with satisfactory precision, and doing so while still compensating for the deformations that occur between images due to differences in posture in which these images are acquired. The residual registration errors still present in the results may be due to several factors. For example, the precision of the manual intervention required for the localization landmarks on all images has inherent limits. X-ray and MRI landmark localisation errors have been studied in previous works to be 2.1±1.5 mm[42] and 3.17±3.3 mm[28]. The MRI landmark localisation errors have been shown to significantly decrease to 1.57±1.13 mm when the centroids are compared, as is the case for our proposed registration method. In the case of the surface topography, where the resolution of the equipment is 1.1 mm, the landmark localisation error is assumed to be equivalent to the radius of the adhesive markers, thus 2.5 mm. The registration error between the surface topography and the X-ray data has been previously established at 2.7 mm[43]. Ongoing work within our group is aimed at automating the landmark extraction process. Furthermore, unlike existing methods, our proposed articulated model method does not require actual correspondence points for registration. This is due to the fact that the center and the orientation of the vertebral bodies are used in order to calculate the transformation between vertebrae. Thus, the need for landmark extraction can be eliminated and replaced with the use of higher order primitives. This has the potential to reduce registration variability and to improve precision. In terms of the correspondence points used to drive the thin-plate spline registration between the MRI and surface topography, the accuracy is not studied. These correspondences might suffer from a lack of accuracy due to the fact that deformations in the z axis between the 2 modalities that are due to gravity have not been taken into account by the articulated model. However, a lack of anatomically significant correspondences makes it difficult to measure the accuracy. Since the landmark selection is automatic, we assume that the precision is mostly dependent on the resolution of the images. In this case, the image with the lower resolution is the MRI. The MRI resolution is 1 mm by 1 mm in-plane and has 2 mm thickness leading to a maximum distance of 2.45 mm between voxels.

Further refinement of the method can be made by incorporating tissue-specific elasticity constraints, which would require tissue-specific biomechanical analysis resulting in a more complex model. This would provide for more realistic anatomical deformations. It should also be noted that tissues contained within the boundaries of the ribcage are believed to deform differently from tissues outside of this boundary. Modeling the transformation of these tissues differently has the potential of improving registration precision. In addition, a thorough biomechanical analysis would allow us to model the effect of gravity on the non-rigid deformations resulting from a posture change. The implication is that gravity would have a different effect on the various anatomical structures being registered. However, the 3D voxel deformations resulting from gravitational forces are very likely to be outside of the plane of the acquired MRI slices. Thus, incorporating these deformations would require a higher resolution in the z plane, which is currently infeasible in a clinical setting due to prohibitive acquisition times.

Conclusions

A method to register MRI, X-ray and surface topography data was proposed. This method first registered surface topography and X-ray data using thin-plate splines, and then fit the MRI data onto the model by taking into account the non-rigid deformations that are due to the postural difference between acquisitions. An articulated model was used in order to approximate the vertebral deformations in the MRI, and the remainder of the soft tissues was deformed using thin-plate splines with vertebral rigidity and the surface topography as constraints. Visual results as well as 2D and 3D Dice values measuring the fit between the surface topography and MRI data were obtained for real data of 3 pre-operative patients with scoliosis in order to validate the proposed method. Both qualitative and quantitative results showed a significant improvement in fitting the MRI data with the X-ray and surface topography data when compared to both rigid and simple articulated model registration. The determinant of the Jacobian of the deformation field was obtained with and without the rigidity constraint, showing higher variations in the deformation in the anterior part of the slices, and showing lower variations closer to the vertebrae in the case of the proposed method. Future work aims at incorporating tissue-specific elasticity constraints to the registration process, and at using automatically extracted higher order primitives for the articulated model registration. The precision of the obtained registration results allows us to build a complete 3D model of a patient’s trunk including soft tissue, vertebral, and trunk surface information which can be incorporated in a surgical simulator under development in order to potentially better predict the outcome of scoliosis treatments.

Competing interests

The authors of this work hereby declare that they have not received reimbursements, fees, funding, or salary from an organization that may in any way gain or lose financially from the publication of this manuscript, either now or in the future. They do not have any other financial competing interests to declare.

Authors’ contributions

RH participated in the development of the methodology, gathering patient data, and carrying out experiments, results, and writing of the manuscript. FC participated in the development of the methodology and in the writing of the manuscript. HL provided patient data and participated in writing of the manuscript. JD participated in the development of the methodology. All authors read and approved the final manuscript.

References

  1. Roach J: Adolescent idiopathic scoliosis.

    Orthop Clin North Am 1999, 30:353-365. PubMed Abstract | Publisher Full Text OpenURL

  2. Dionne O, Assi KC, Grenier S, Labelle H, Guibault F, Cheriet F: Simulation of the postoperative trunk appearance in scoliosis surgery. In to appear in: International Symposium on Biomedical Imaging (ISBI). Newark, NJ: IEEE; 2012. OpenURL

  3. III WW, Viola P, Atsumi H, Nakajima S, Kikinis R: Multi-Modal Volume Registration by Maximization of Mutual Information.

    1996, 1:35-51.

  4. Rohr K, Fornefett M, Stiehl HS: Approximating Thin-Plate Splines for Elastic Registration: Integration of Landmark Errors and Orientation Attributes. In Proc. 16th Intl. Conf.on Information Processing in Medical Imaging. London, UK: Springer-Verlag; 1999:252-265. OpenURL

  5. Rueckert D, Sonoda LI, Hayes C, Hill MDLGLeach, Hawkes DJ: Nonrigid registration using free-form deformations: Application to breast MR images.

    IEEE Trans. Med. Imag 1999, 18:712-721. Publisher Full Text OpenURL

  6. Blackall JM, Penney GP, King AP, Hawkes DJ: Alignment of sparse freehand 3–D ultrasound with preoperative images of the liver using models of respiratory motion and deformation.

    IEEE transactions on medical imaging 2005, 24(11):1405-1416. PubMed Abstract | Publisher Full Text OpenURL

  7. Skerl D, Likar B, Pernus F: Evaluation of similarity measures for non-rigid registration. In Biomedical Image Registration. Third International Workshop, WBIR 2006. Proceedings. Utrecht, Netherlands: Springer-Verlag; 2006:160-8. OpenURL

  8. Loeckx D, Maes F, Vandermeulen D, Suetens P: Nonrigid image registration using free-form deformations with a local rigidity constraint. In Medical Image Computing and Computer-Assisted Intervention - MICCAI 2002. 5th International Conference. Proceedings Part II. Berlin/Heidelberg: Springer-Verlag; 2004:461-8. OpenURL

  9. Bajcsy R, Kovačič S: Multiresolution elastic matching.

    Computer Vision, Graphics, and Image Processing 1989, 46:1-21.

    [http://www.sciencedirect.com/science/article/pii/S0734189X89800143 webcite]

    Publisher Full Text OpenURL

  10. Christensen GE, Rabbitt RD, Miller MI: 3D brain mapping using a deformable neuroanatomy.

    Phys Med Biol 1994, 39(3):609. PubMed Abstract | Publisher Full Text OpenURL

  11. Bro-Nielsen M, Gramkow C: Fast Fluid Registration of Medical Images. In SPIE Visualization in Biomedical Computing Volume 1131. Edited by Hone K, Kikinis Rlncs. Berlin/Heidelberg: Springer; 1996:267-276. OpenURL

  12. D’Agostino E, Maes F, Vandermeulen D, Suetens P: A viscous fluid model for multimodal non-rigid image registration using mutual information.

    Med Image Anal 2003, 7(4):565-575. PubMed Abstract | Publisher Full Text OpenURL

  13. Tang S, Jiang T: Fast non-rigid medical image registration by fluid model. In Asian Conference on Computer Vision. Berlin/Heidelberg: Springer, Lecture Notes in Computer Science (LNCS); 2004:27-30. OpenURL

  14. Mahapatra D, Sun Y: Nonrigid Registration of Dynamic Renal MR Images Using a Saliency Based MRF Model. In Proceedings of the 11th international conference on Medical Image Computing and Computer-Assisted Intervention - Part I , MICCAI ’08. Berlin, Heidelberg: Springer-Verlag; 2008:771-779.

    [http://dx.doi.org/10.1007/978-3-540-85988-8_92 webcite]

    PubMed Abstract OpenURL

  15. Glocker B, Komodakis N, Paragios N, Navab N: Non-rigid Registration using Discrete MRFs: Application to Thoracic CT Images. In Workshop Evaluation of Methods for Pulmonary Image Registration in conjunction with Medical Image Computing and Computer-Assisted Intervention. Berlin/Heidelberg: Springer, Lecture Notes in Computer Science (LNCS); 2010-2010. OpenURL

  16. So RWK, Tang TWH, Chung ACS: Non-rigid image registration of brain magnetic resonance images using graph-cuts.

    Pattern Recogn 2011, 44(10-11):2450-2467. Publisher Full Text OpenURL

  17. Alperin D, Pelizzari C: Retrospective registration of X-ray angiograms with MR images by using vessels as intrinsic landmarks.

    J Magn Reson Imaging 1994, 4(2):139-44. PubMed Abstract | Publisher Full Text OpenURL

  18. Chan H, Chung A, Yu S, Wells W: 2D-3D vascular registration between digital subtraction angiographic (DSA) and magnetic resonance angiographic (MRA) images vol. 1. In IEEE international symposium on biomedical imaging: nano to macro. New York, NY: Institute of Electrical and Electronics Engineers (IEEE); 2004:708-11. OpenURL

  19. Harmouche R, Cheriet F, Labelle H, Dansereau J: Three-dimensional reconstruction and fusion for multi-modality spinal imagesl.

    Computerized Medical Imaging and Graphics 2004, 28(1-2):21-31. PubMed Abstract | Publisher Full Text OpenURL

  20. Vermandel M, Betrouni N, Gauvrit J, Pasquier D, Vasseur C, Rousseau J: Intrinsic 2D/3D registration based on a hybrid approach: use in the radiosurgical imaging process.

    Cell Mol Biol 2006, 52(6):44-53. OpenURL

  21. Chung1 C, Wells WMI, Norbash A, L GWE: Multi-modal image registration by minimising Kullback-Leibler distance. In Fifth international conference on medical image computing and computer-assisted intervention (MICCAI 2002) vol. 2489. Lecture Notes in Computer Science, Springer; 2002:525-32. OpenURL

  22. Florin C, Williams J, Khamene A, Paragios N: Registration of 3D angiographic and X-ray images using sequential Monte Carlo sampling. In First international workshop on computer vision for biomedical image applications (CVBIA 2005), vol 3765. Lecture Notes in Computer Science, Springer; 2005:427-36. OpenURL

  23. Jomier J, Bullitt E, Van Horn M, Pathak C, Aylward S: 3D/2D model-to-image registration applied to tips surgery. In Ninth international conference on medical image computing and computerassisted intervention (MICCAI 2006), part 2, vol 4191. Lecture Notes in Computer Science, Springer; 2006:662-9. OpenURL

  24. Kita Y, Wilson D, Noble J: Real-time registration of 3D cerebral vessels to Xray angiograms. In First international conference on medical image computing and computer-assisted intervention (MICCAI-98), vol. 1496. Lecture Notes in Computer Science, Springer; 1998:1125-33. OpenURL

  25. McLaughlin R, Hipwell J, Hawkes D, Noble J, Byrne J, Cox T: A comparison of a similarity-based and a feature-based 2-D–3-D registration method for neurointerventional use.

    IEEE Trans Med Imaging 2005, 24(8):1058-66. PubMed Abstract | Publisher Full Text OpenURL

  26. Hipwell J, Penney G, McLaughlin R, Rhode K, Summers P, Cox T: Intensity-based 2-D–3-D registration of cerebral angiograms.

    IEEE Trans Med Imaging 2003, 22(11):1417-26. PubMed Abstract | Publisher Full Text OpenURL

  27. Qi W, Gu L, Xu J: Non-rigid 2D-3D registration based on support vector regression estimated similarity metric. In Fourth international workshop on medical imaging and augmented reality, vol. 5128. Lecture Notes in Computer Science, Springer; 2005:339-48. OpenURL

  28. Harmouche R, Cheriet F, Labelle H, Dansereau J: 3D registration of MR and X-ray spine images using an articulated model.

    Comp. Med. Imag. and Graph 2012, 36(5):410-418. OpenURL

  29. Markelj P, Tomazevic D, Likar B, Pernus F: A review of 3D/2D 638 registration methods for image-guided interventions.

    Med Image Anal 2010, 16(3):642-61. PubMed Abstract | Publisher Full Text OpenURL

  30. Little JA, Hill DLG, Hawkes DJ: Deformations incorporating rigid structures. In Proc. Workshop on Mathematical Methods in Biomedical Image Analysis. San Francisco, CA: IEEE Comput. Soc. Press; 1996:104-13. OpenURL

  31. Huesman RH, Klein GJ, Kimdon JA, Kuo C, Majumdar S: Deformable registration of multi-modal data including rigid structures. In 2002 IEEE Nuclear Science Symposium Conference Record, Volume 3. New York, NY: Institute of Electrical and Electronics Engineers (IEEE); 2003:1879-82. OpenURL

  32. Wang K, He Y, Qin H: Incorporating rigid structures in non-rigid registration using triangular b-splines. In Variational, geometric and level set methods in computer vision (VLSM). Berlin/Heidelberg: Springer, Lecture Notes in Computer Science (LNCS); 2005:235-246. OpenURL

  33. van de Kraats E, van Walsum T, Verlaan J, Oner F, Viergever M, Niessen W: Noninvasive Magnetic Resonance to Three-Dimensional Rotational X-Ray Registration of Vertebral Bodies for Image-Guided Spine Surgery.

    Spine 2004, 29(3):293-297. PubMed Abstract | Publisher Full Text OpenURL

  34. Tomazevic D, Likar B, Pernus F: 3-D/2-D registration by integrating 2-D information in 3-D.

    IEEE Trans. Med. Imaging 2006, 25:17-27. PubMed Abstract | Publisher Full Text OpenURL

  35. van de Kraats E, Penney G, van Walsum T, Niessen W: Multispectral MR to X-ray registration of vertebral bodies by generating CT-like data. In Conference on Medical Image Computing and Computer-Assisted Intervention (MICCAI 2005), Part 2. Berlin/Heidelberg: Springer, Lecture Notes in Computer Science LNCS; 2005-2005. OpenURL

  36. Markelj P, Tomazevic D, Pernus F, Likar B: Robust gradientbased 3-D/2-D registration of CT and MR to X-ray images.

    IEEE Transactions on Medical Imaging 2008., 27(12) OpenURL

  37. Tang Z, Pauli J: Fully automatic extraction of human spine curve from MR images using methods of efficient intervertebral disk extraction and vertebra registration.

    Int J CARS 2011, 6:21-33. Publisher Full Text OpenURL

  38. Cheriet F, Dansereau J, Petit Y, Labelle H de Guise: Towards the Self-Calibration of a Multi-view Radiographic Imaging System for the 3D Reconstruction of the Human Spine and Rib Cage.

    International Journal of Pattern Recognition and Artificial Intelligence 1999, 13(5):761-779. Publisher Full Text OpenURL

  39. Bookstein F: Principal warps: thin-plate splines and the decomposition of deformations.

    IEEE Trans. Pattern Analysis and Machine Intelligence 1989, 11(6):567-585. Publisher Full Text OpenURL

  40. Dice LR: Measures of the amount of ecologic association between species.

    Ecology 1945, 26(3):1037-1044. OpenURL

  41. Christensen G, Rabbitt R, Miller M: Deformable templates using large deformation kinematics.

    Image Processing, IEEE Transactions on 1996, 5(10):1435-1447. Publisher Full Text OpenURL

  42. Chevrefils C, Cheriet F, Aubin CE, Grimard G: Texture Analysis for Automatic Segmentation of Intervertebral Disks of Scoliotic Spines From MR Images.

    IEEE Transactions on Information Technology in Biomedicine 2009, 13(4):608-620. PubMed Abstract | Publisher Full Text OpenURL

  43. Liao QH: Fusion of deformable anatomical structures of human torso. M.sc.a. Ecole Polytechnique: Montreal (Canada); 2005. [Relevant - have]

Pre-publication history

The pre-publication history for this paper can be accessed here:

http://www.biomedcentral.com/1471-2342/13/1/prepub