|Home | About | Journals | Submit | Contact Us | Français|
In this study, a combination of active shape model (ASM) and active appearance model (AAM) was used to segment the left and right ventricles of normal and Tetralogy of Fallot (TOF) hearts on 4-D (3-D+time) MR images. For each ventricle, a 4-D model was first used to achieve robust preliminary segmentation on all cardiac phases simultaneously and a 3-D model was then applied to each phase to improve local accuracy while maintaining the overall robustness of the 4-D segmentation. On 25 normal and 25 TOF hearts, in comparison to the expert traced independent standard, our comprehensive performance assessment showed subvoxel segmentation accuracy, high overlap ratios, good ventricular volume correlations, and small percent volume differences. Following 4-D segmentation, novel quantitative shape and motion features were extracted using shape information, volume-time and dV/dt curves, analyzed and used for disease status classification. Automated discrimination between normal/TOF subjects achieved 90%–100% sensitivity and specificity. The features obtained from TOF hearts show higher variability compared to normal subjects, suggesting their potential use as disease progression indicators. The abnormal shape and motion variations of the TOF hearts were accurately captured by both the segmentation and feature characterization.
Tetralogy of Fallot (TOF) is a common congenital heart disease characterized by four typical features: right ventricular hypertrophy, right ventricular outflow tract (RVOT) obstruction, ventricular septal defect (VSD), and overriding aorta. Initial management of infants with TOF includes palliative surgery to augment the RVOT and close the VSD, but patients often suffer from pulmonary valve insufficiency after surgery. This can over time lead to right ventricular dilation and dysfunction eventually requiring pulmonary valve replacement (PVR) in the adult years –. Although TOF patients often have similar life expectancy as healthy people with an estimated prevalence of 80000 in the U.S. , they need to be monitored throughout their lives to identify whether further surgery is required. Cardiac magnetic resonance (MR) imaging has become essential for pre- and post-operative management of TOF patients to determine the necessity and optimal timing of PVR –. It has been shown that RV dilation, RV systolic dysfunction and the left ventricle (LV) systolic dysfunction are independent predictors of progressive LV dysfunction . It has also been suggested that unfavorable LV and RV interaction may be an important indicator for poor post-operation outcome and for timing of PVR , .
The traditional functional indices measured in clinical analysis such as ventricular volumes and ejection fraction may not be sufficient to allow timely patient management and prevent long term sequelae ,  because they do not fully utilize the rich information provided by the MR data to describe the complex ventricular remodeling process. In addition, due to the inevitable and large observer variability associated with the manual data analysis, it is difficult to utilize these traditional indices in long-term patient care to produce reliable prediction of disease progression. It has been found that the post-operative TOF RV with near-normal ejection fraction exhibits unique abnormal ventricular remodeling characteristics such as larger and rounder cross-sectional shape and more basal bulging , but the quantitative relationships between these shape features and clinical outcomes are still unknown.
To improve the quality of clinical analysis for TOF and other cardiac diseases, it is essential to have quantitative indices that are capable of describing and measuring shape and motion features. Furthermore, highly-automated and accurate methods must be developed to perform 4-D (3-D+time) segmentation and reduce observer variability.
The cardiac image segmentation problem has been intensively studied by many groups in the recent years. Although there are still some feature-based methods such as using threshold and watershed , most groups acknowledged the necessity of using deformable models or parametric models to incorporate anatomical knowledge into the segmentation.
The deformable model method often starts from an atlas (template mesh). The anatomic knowledge is included in the atlas and its feasible deformation is learned from training or clinical models. During segmentation, the atlas is deformed onto the target image using registration or transform guided by image features. Lötjönen et al.  segmented all four heart chambers in 3-D MR data using landmark probability distribution, probability atlas, and nonrigid registration. Ecabert et al.  proposed a shape-constrained deformable model for 3-D CT to combine the strengths of parametric and deformable models and to segment all four heart chambers using piecewise affine transform optimization. Bistoquet et al.  introduced an incompressible myocardium model to improve the clinical feasibility of LV and RV segmentations in 4-D MR data. Sermesant et al.  combined the intensity and deformable models to achieve 3-D biventricular segmentations in multiple modalities. Zheng et al.  constructed models of all four heart chambers with holes for inflow/outflow tracts and achieved segmentation initialization through marginal space learning and steerable features.
Active shape model (ASM) and active appearance model (AAM) are popular parametric models. The statistical property of the target shape or appearance (combination of shape and texture) is derived from a training set and denoted by compact modal indices. ASM and AAM have been used for various applications such as: 2-D+time LV segmentation of ultrasound data , 3-D LV or LV+RV segmentation of ultrasound and MR data , bi-temporal LV segmentation of MR data , 3-D LV segmentation of MR and CT data , 4-D LV shape fitting of MR and CT data , and modeling four heart chambers of 4-D CT data . Andreopoulos et al.  achieved LV segmentation of 4-D MR images using 3-D AAM followed by 2-D+time ASM and found time-dependent segmentation errors—larger errors at end-systole than at enddiastole. Zhu et al.  performed LV segmentation of 4-D MR data in a phase-by-phase manner where the statistical and dynamic models were combined in 3-D segmentation. Several methods have been proposed to improve the quality and capability of parametric models using independent component analysis for modeling and segmentation , , enlarging ASM training sets , using sparse modeling with orthomax criterion , improving landmark correspondence , , and modifying ASM optimization , .
Most of the approaches to cardiac image processing, especially for MR, share several common limitations.
This study utilizes ASM and AAM to accomplish two main goals. 1) Develop an LV and RV segmentation method that is highly automated, applicable in clinical environment and capable of producing accurate 4-D segmentations of clinical-quality MR data. 2) Derive novel 4-D ventricular function indices and test their capabilities for distinguishing between normal and TOF hearts, and assessing post-operative status. A complete cardiac MR analysis pipeline was designed from defining an independent standard to disease characterization.
The unique contributions of this study are as follows.
The point distribution model (PDM)  describes the shape variation of a training set of objects. The coordinates of landmarks—a set of corresponding points representing the boundary—are concatenated to form a shape vector s, whose statistical properties are identified via principal component analysis (PCA) yielding
where is the mean shape and matrix Φ contains PCA modes (a.k.a. eigenvectors or components) ordered by their significance. After discarding less significant modes, the weight vector b (referred to as the modal indices in this paper) is used as compact alternative to approximate any shape s in the training set or synthesize a new feasible shape that is not in the training set but still within the allowed range of variations. The application of PDM in segmentation is ASM . The local intensities along scan lines associated with landmarks are modeled from the training set. During ASM segmentation, the landmarks are iteratively updated based on image features and trained scan line models, then new modal indices are computed and restricted within the allowed range.
Similarly, the intensity pattern—the texture vector t—of the object can also be modeled by a statistical texture model once the correspondence of all interior voxels is established. Combining shape modal indices (bs) and texture modal indices (bt) with appropriate weighting Ws results in an appearance vector a that can be statistically modeled by AAM
The AAM segmentation  is an iterative process often referred to as the model matching. The appearance modal indices ba synthesize model texture and model shape. The model shape is mapped to the target image as the target shape using an affine transform. Another texture, the target texture, is sampled based on the target shape. The goal of the model matching is to find the optimal model parameters d (ba and the affine transform parameters) that minimize the differences between the model texture and the target texture. The classical AAM implementation minimizes the root-mean-square of the texture difference. The model parameters are iteratively updated using knowledge acquired from model training. The training typically employs reduced-rank multivariate linear regression on a set of known model parameters to find a Jacobian matrix R that is used to predict the model parameter displacement Δd based on the texture difference r(d) as
ASM and AAM have their strengths and weaknesses. AAM is good at object tracking but it may impose too stringent of constraints and the matching can be easily trapped at local minima. ASM is good at finding local features but it requires a good initialization and can be disturbed by neighboring nontarget objects. The strengths of AAM and ASM were combined in a multistage hybrid approach proposed by Mitchell et al. . Its usage in MR ventricular segmentation yielded better results than AAM or ASM alone.
The short-axis (SA) MR images are typically traced in clinical analysis because they show clearly-identifiable cross-sectional views of the doughnut-shaped LV and the crescent-shaped RV depicted on 2-D slices. Due to the typically large slice thickness of 6–8 mm as opposed to the in-plane resolution of 1.5–2.5 mm, identifying the traceable borders in the last basal image is difficult and is the major source of volumetric measurement errors , . However, the consistency of ventricular coverage is often given a higher priority in clinical practice. In the complementary set of standard MR images acquired from four chamber long-axis (LA) planes, the ventricular bases and apexes are easy to identify, but tracing LV endocardial borders is more troublesome due to difficulty of accurately identifying papillary muscles which are commonly included in the LV blood pool. Several proposed methods combined SA and LA tracings together to produce better LV shapes ,  using simple geometric assumptions that are not feasible for the RV due to its complex shape. Cardiac MR images are also affected by the respiratory motion that can be simply modeled as misalignment of the ventricles between slices. The motion correction that uses only SA images  is feasible but not appropriate for shape analysis because it can potentially change the actual orientation and shape of the heart.
The complementary information present in the SA and LA MR data is fully utilized by our image fusion technique. The motion correction method, although similar to Lötjönen’s approach , was developed independently. In addition, the motion-corrected data from both orientations are fused into a single 4-D image with isotropic voxels and improved quality as shown below.
The number of phases (time frames) acquired per cardiac cycle in the MR scan is determined by the imaging protocol and the subject’s heart rate so it varies between 18 and 25. We first normalized our data to 16 phases per cycle using nearest neighbor (along time axis) interpolation that prevents the cardiac motion from introducing artifacts in a phase-by-phase fashion. The potential temporal motion change due to interpolation is very small in practice and is only limited to a very small portion of the cardiac cycle. Without motion correction, the end-diastolic isotropic 3-D images reconstructed by shifted linear interpolation  from SA and LA MR data are shown in Fig. 1(a) and (b).
After establishing correspondence between SA and LA MR data using imaging plane description in the DICOM file header, the motion in the SA MR data is corrected by registering each SA MR data [bottom row of Fig. 1(a)] to its counterpart in the 3-D image reconstructed from LA MR data [bottom row of Fig. 1(b)]. The registration is restricted to translation only and the criterion is minimal intensity difference. The result of motion correction is shown in Fig. 1(c). The motion correction of the LA MR data is then performed in the same way utilizing information from the motion-corrected SA MR data. The two images reconstructed from motion-corrected SA and LA MR data are fused into a single image as shown in Fig. 1(d), where the average intensity of the two reconstructed images was used in the overlapped region. The fused image has better coverage of the heart as well as improved image quality—note the difference in LA slices of Fig. 1(c) and 1(d). For SA MR data with d mm in-plane resolution, the isotropic voxel size of the fused image is d mm.
A custom-designed application is used for manual tracing of the fused image as shown in Fig. 2. This tracing environment provides sufficient 4-D context such that the ventricular apex and base locations, end-diastole and end-systole can be easily identified. The manual tracing is performed on user-selected short-axis slices as 2-D contours defined by spline control points. The control points are not used as landmarks and their number is determined by the user such that the resulting contour correctly describes the ventricular border. The contours are translated into 2-D distance maps. The 3-D shape on each phase is created by distance-based shape interpolation  and is used as the independent standard. The colored overlay shows the created shapes and provides the user with feedback regarding the tracing accuracy.
In order to automatically create landmarks on all sample shapes, a template-based general framework proposed by Frangi et al.  is adapted but implemented differently as illustrated in Fig. 3. The whole process is performed in a phase-by-phase manner treating a 4-D sample as 16 3-D shapes. It consists of the following steps.
A customized iterative closest point algorithm is used to find the optimal affine with anisotropic scaling such that later on the transformed landmarks are brought as close to the sample shape as possible and the deformation associated with landmark propagation is reduced.
The landmarks on the template shape are created by three steps using VTK-designed filters: surface meshing to create dense triangular mesh on the template surface , surface smoothing to remove possible sharp edges , and triangle decimation to reduce the number of vertices and triangles . The LV and RV are landmarked individually in this study using their own templates as shown in Fig. 4. For LV, its endocardial surface is landmarked first. To form landmarks for the epicardial surface, an exact copy of the landmarked endocardial surface is created and its vertices are moved along the surface normals to the epicardial surface. Both LV surfaces are therefore defined by the same number of vertices and triangles and exhibit the same topology.
Because the transformed land-marks are reasonably close to the sample surface, an efficient method is designed. The landmarks, pi = (xi,yi,zi), are iteratively updated by an implicit function that minimizes the cost that consists of a similarity cost Cs and a regularization constraint Cr with appropriate weightings ws and wr.
The similarity cost Cs is defined as
where D is the distance map of the target sample shape L computed by distance transform D. It can be derived that iteratively minimizing Cs is equivalent to iteratively moving pi along the gradient direction of D(pi) as
Minimizing the similarity cost forces the landmarks onto the target surface. Because it is performed on each landmark without considering neighboring structure of the mesh, geometric and topological errors could be introduced. The regularization constraint Cr is used to preserve the spatial smoothness of the deformation. It is defined as the Dirichlet integral and can be calculated for the surface mesh by
where |aj | is the length of the jth edge, αj and βj are two angles associated with the edge. Their definitions are shown in Fig. 5 where the jth edge is the edge between vertices (landmarks) pi and qj. It can be shown that when Cr is minimized, vertex pi must satisfy 
Combining the optimization of Cs and Cr together, the landmark pi is iteratively updated by
In this study, the regularization constraint is not used to produce a discrete minimal surface as originally proposed , but instead to preserve the smoothness of the deformation. The iterative deformation defined by (9) therefore cannot introduce abrupt changes to the relative sizes and shapes of the triangles associated with any mesh vertex. The proper weightings, ws and wr are empirical values related to number of iterations. In this study, 200 iterations of landmark propagation are performed with ws = 0.03, wr = 0.01.
In this study, we use virtual surfaces to restrict the size of the texture vector such that it only contains intensity information that is useful for the segmentation. Two virtual surfaces are defined for each ventricle: one by enlarging the outermost surface (scaling it with respect to its centroid by 1.2), the other by shrinking the innermost surface (using a scale of 0.8). The vertices of the virtual surfaces are not used in the shape model as landmarks, but whenever the actual ventricular surface is changed during model matching, its associated virtual surfaces are recomputed and the corresponding texture vectors are re-sampled accordingly. In the 3-D space where the landmarking template is defined, uniformly spaced voxels between the two virtual surfaces are used to define texture mapping. The voxel correspondences are established using a barycentric coordinate system, where a voxel’s barycentric coordinates are defined by the vertices (landmarks) of the tetrahedron it resides in. Because all the virtual and actual surfaces of a ventricle are designed to have the same topology, building the tetrahedrons for texture mapping becomes trivial—connecting the vertices of the two corresponding triangles from the neighboring surfaces defines three tetrahedrons.
The basic approach of the multistage hybrid ASM/AAM method proposed by Mitchell et al.  is adapted and implemented with additional new improvements. An AAM is first used to produce a shape sa. Then an ASM starts from sa and produces another shape ss. These two shapes are combined with an appropriate weight w as
The new combined shape s is then used to compute the new appearance modal indices using (2) to start a new AAM model matching until convergence. This hybrid approach provides extra momentum to bring model matching out of local minima and thus increases the chance of finding a global minimum. However, the appearance modal indices are associated with both shape and texture variations and the improvements in segmentation accuracy embedded in the combined shape s may not be fully incorporated into the resulting new appearance modal indices. In addition, an improved shape computed from updated appearance modal indices may be rejected as a better segmentation because the simultaneously updated texture modal indices may lead to increased matching error. To decouple this tie between shape and texture in the appearance model, a special segmentation fine tuning mechanism—brute-force reversed ASM—is designed as follows.
The traditional ASM first updates landmarks using image features and then applies the model constraints. Our method first creates a new shape by directly updating the shape modal indices and then tests whether the new shape fits the target image better based on the texture. Therefore, it is a reversed ASM method. Although it preforms a brute-force search in the model space and is more computationally expensive than traditional ASM, when used to fine tune the AAM segmentation in practice, it only uses 5%–10% of the total running time.
The full process of our hybrid ASM/AAM method is as follows.
Balanced steady-state free precession cardiac MR imaging was performed on 25 patients with repaired TOF and resultant pulmonic regurgitation (12 males, 13 females; ages: 18–60, median = 30). The study population also included 25 normal subjects (13 males, 12 females; ages: 22–38, median = 29). The normal subjects were scanned on a Siemens Avanto scanner (in-plane resolution 2.05 ±0.15 mm). TOF patients were scanned on a GE Signa scanner (in-plane resolution 1.48 ± 0.19 mm). The slice thickness was 6–8 mm and the slice gap was 0–2 mm. The acquired SA slices covered the whole heart and LA slices covered at least 70% of the heart. To prevent body fat from influencing motion correction, the MR data were cropped using two manually defined rectangular regions that enclosed the complete LV and RV on short-axis and long-axis views, respectively. The cropped MR data were motion corrected and fused into 4-D datasets as described in Section III-A.
The fused 4-D datasets were manually traced by a cardiologist on short-axis slices using our tracing application. Three ventricular borders were traced: LV endocardium, LV epicardium, and RV epicardium. The papillary muscles were included in the LV blood pool enclosed by the LV endocardial border. The manually traced slices were not limited to slices present in the original MR imaging planes and not all the interpolated slices were traced. When the ventricles exhibited small shape changes between slices, the distance between traced slices was about 6–8 mm. In the basal region, where the ventricles often show sudden 2-D shape changes, the traced slices were often 2–4 mm apart. After the initial pass of manual tracing, the tracings were repeatedly edited. The ventricular shapes were created and loaded into the tracing application as shown in Fig. 2, and the quantitative volume measurements of all cardiac phases were recorded. Working in consensus with a second cardiologist, the expert who performed the tracing then checked the volume measures and observed the resulting shapes to make appropriate small changes on the tracings until the ventricular shapes showed correctly identified end-diastole, end-systole, base and apex locations, and smooth 3-D shape and motion. Extra effort was made to guarantee correct depiction of the LV myocardial motion and the constancy of LV mass throughout the cardiac cycle. For normal subjects, stroke volumes of the LV and RV were expected to be approximately the same.
To minimize the loss of segmentation accuracy due to a small training set and large shape variability, the normal and TOF hearts as well as the LV and RV were modeled and segmented separately. The LV or RV was segmented using a two-step approach. In the first step, the preliminary segmentation was performed using a 4-D model. It was manually initialized by fitting the 4-D mean shape onto the target image. The fitting was visually inspected and adjusted mainly on the first cardiac phase, which corresponds to the peak of R-wave in EKG-gated MR acquisition. The full hybrid ASM/AAM process listed in Section V was applied with M = 10 and J = 10. In the second step, the preliminary segmentation was used as initialization of the final segmentation which was performed in a phase-by-phase manner using a single 3-D model that was created from 3-D shape/texture samples of all cardiac phases. The reversed ASM fine tuning was not performed in the final segmentation and M = 20 was used. In both steps, P = 1 and w = 0.35 were used.
The segmentation performance was assessed using a hold-out validation strategy. All of the 25 4-D datasets in the normal or TOF population were randomly divided into five groups each containing five datasets. The datasets in any group were segmented using models created and trained from datasets in the other four groups. The sizes of the training sets were therefore 20 for 4-D models and 320 for 3-D models. This approach guarantees that a complete independence of the training and testing sets is maintained throughout the multistage process. The ventricle (VE, either LV or RV) produced by any stage was compared to the independent standard (IS) using labeled image—a voxel is labeled as background or one of the several objects—that has the same resolution as the 4-D fused image.
The following measurements were computed on individual cardiac phases of all datasets.
The following metrics were derived to assess the performance.
To visualize the time-dependence of the performance, the following plots were used.
The objective of this task is to identify quantitative features that are unique to TOF and its progression. The shape and motion features were derived from the labeled images of the ventricles. The cardiac motion features were derived from the volumetric measurements. The independent standard and the final segmentation were analyzed to test whether they can extract the same or similar disease features.
The 4-D PDMs of the LV, RV, and LV+RV combined—by concatenating shape vectors of LV and RV—were constructed using all 50 datasets as the training set. The ventricular landmarks were generated again from the labeled volumes of the independent standard or the final segmentation. Having a single training set that contains an equal number of normal and TOF hearts is suitable to identify TOF features in the context of normal ventricles. The shape modal indices were computed and analyzed.
The nVTCs of the LV, RV, and LV+RV were used to represent the cardiac motion. Each nVTC was treated as a vector including ventricular volumes of all cardiac phases. The normalization by EDV removes the overall size variation and therefore performs a similar task to the Procrustes Analysis in PDM construction. In addition, dV/dt curves were also computed from the nVTCs and analyzed. The nVTCs and dV/dt curves do not include any shape knowledge. The nVTCs contain information about ejection fraction while dV/dt curves include no direct ejection fraction knowledge and focus on the rate of volume changes between phases. Extracting modal indices from nVTCs of all datasets was performed by PCA.
The modal indices extracted from the shape model or volume-based curve model were fed into a classifier using linear discriminative analysis to separate TOF patients from normal subjects. The classification performance was evaluated by leave-one-out validation and measured by specificity and sensitivity. In this binary classification of positive TOF and negative normal hearts, the results are represented by numbers of true positive (TP), false positive (FP), true negative (TN), and false negative (FN). The specificity and sensitivity were computed by
The classifications were performed using one to five strongest modal indices and the best combination of specificity and sensitivity achieved was reported as the final classification performance.
The motion correction and image fusion were performed on all datasets and assessed by visual inspection. The resulting fused images had no remaining motion artifact and no image quality degradation. Due to the time-consuming nature of the manual tracing, observer variability was not assessed. Again, 50 fused 4-D images were traced by a cardiologist and further inspected and approved by a senior cardiology expert. The total time spent on tracing and editing on each 4-D dataset by the two experts was about 6–8 h and at least 50% of the time was spent on the initial contour drawing.
The numbers of landmarks used for 3-D ventricular surfaces were: 400 for LV endocardium, 400 for LV epicardium, and 800 for RV epicardium. The total running time of all automated land-marking steps performed on all 50 LVs and 50 RVs was approximately one hour on a PC with Intel Core 2 Duo E6400 processor and 2 GB of RAM. The landmarking error measured by average landmark-to-surface distance was 0.00 ± 0.02 voxels, but it does not necessarily mean that the landmark-defined surfaces accurately represent the actual ventricular shapes. The performance metrics defined in Section VI-C were used to assess the land-marking quality by comparing the landmark-defined ventricles and the independent standard. The results are shown in Fig. 6(b) and Table I and Table II. The average surface positioning errors show subvoxel accuracy in mean and standard deviation. The signed errors and volume-based metrics show a weak trend towards underestimation. The volume curves in Fig. 6(b) are almost the same as those of the independent standard in Fig. 6(a) to a naked eye, and there are slightly more pronounced fluctuations in the distribution of the LV mass curves. The mean percent volume differences of LV and LV mass begin to show a time-dependent pattern even though they are almost constant within the cardiac cycle and the standard deviations are very small.
The inherent consequence of using a small number of modal indices to represent complex shape and motion is the inevitable loss of detail. When the training set is small, such loss could affect the expected best segmentation accuracy. Using the same hold-out strategy as used in the segmentation, the appearance modal indices of all datasets were computed and used to synthesize the ventricular shapes approximated by the model. This procedure can be thought of as an ideal segmentation when all the model parameters that need to be optimized are known. The synthesized shapes were compared with the independent standard and the results are shown in Fig. 6(c) and (d) in Table I and Table II.
From 4-D to 3-D modeling, the standard deviations of average positioning errors decrease by 0.1–0.3 voxels. The means of absolute errors decrease by 0.3 voxels in normal hearts and 0.5 voxels in TOF hearts. Such improved accuracy of the 3-D modeling is also shown in the volumetric metrics where the mean overlaps increase by 2%–4% and the standard deviations of percent volume differences decrease by 4%–6%. The important improvements introduced by these relatively small numbers can be seen in Fig. 6(c) and (d). The 4-D modeling results show more clustered distribution of nVTCs and wider spread of percent volume differences than for the landmarking results, while the distributions of nVTCs and percent volume differences for the 3-D modeling results are very similar to those for the landmarking results.
Training Jacobian Matrices for the 3-D or 4-D AAM of 20 hearts on the same PC used for the landmarking took approximately one hour for normal LV or RV and 90 min for TOF LV or RV. The manual segmentation initializations for LV and RV took less than 5 min, and consisted of the user adjusting the transform parameters to align the transformed mean shape onto the first phase of the target 4-D image such that the base and apex locations were correctly defined and the ventricular sizes were approximated. The preliminary LV or RV segmentation in 4-D took approximately 2–3 min. The final segmentation in 3-D took approximately 5 min for normal LV or RV and 10 min for TOF LV or RV. The total time used to segment a 4-D dataset was approximately 15 min for a normal heart and 30 min for a TOF heart. The TOF segmentation had longer running time mainly because the texture vector of the enlarged RV—even sampled at similar grid density as the normal RV—is substantially longer.
The performance assessments of preliminary and final segmentations are shown in Fig. 6(e) and (f) and Table I and Table II. Compared with their corresponding theoretical optima suggested by the modeling quality assessment, the means of the absolute positioning errors are off by 0.2–0.3 voxels for normal and TOF hearts, the standard deviations are off by 0.1–0.2 voxels for normal hearts and 0.2–0.3 voxels for TOF hearts.
The final segmentation produced better accuracy than the preliminary segmentation. The means of the absolute positioning errors are reduced by 0.1 voxels for normal hearts and 0.1–0.4 voxels for TOF hearts. The overlap ratios show increased mean and decreased standard deviation. In the volumetric metrics for LV and RV, better agreement with the independent standard is shown in increased slopes and decreased interceptions of linear regressions, increased correlation coefficients, and closer-to-zero means and decreased standard deviations of percent volume differences. The LV and RV nVTCs of the final segmentation also show improvements in motion tracking over the preliminary segmentation, especially in the rapid-fill diastole stage where the LV undergoes rapid expansion after reaching end-systole. However, the final segmentation show decreased LV mass accuracy.
Fig. 7 shows the independent standard and segmentation result on a normal heart and a TOF heart at cardiac phases 0, 3, 7, and 11. Both preliminary and final segmentations produced smooth and accurate ventricular shapes for the normal heart. The final segmentation recovered some shape details that were not captured by the preliminary segmentation, and produced better RV short-axis borders. However, the final segmentation showed reduced accuracy in myocardial thickness on phases 3 and 7 due to inaccurate LV epicardial tracking.
The TOF-affected heart in Fig. 7 has a unique motion pattern. It almost reaches end-systole early at phase 3, then expands quickly and reaches the end of rapid-filling diastole at phase 7. The preliminary segmentation could not fully capture this unique motion pattern as in Fig. 7(d). The final segmentation achieved correct cardiac motion tracking as in Fig. 7(f), where the LV endocardial surface is accurately identified. However, the final segmentation also introduced an underestimated LV myocardium at phase 3 of Fig. 7(f). Note that the user accidentally missed one apex SA slice that should have been traced on phase 3, but the AAM segmentations identified the correct apex location.
Fig. 8 shows the volume curves of the independent standard, preliminary segmentation, and final segmentation of two TOF hearts. The final segmentation’s capability of producing better cardiac motions is demonstrated in Fig. 8(a) for LV motion and Fig. 8(b) for RV motion.
To test the benefits of the reversed ASM fine tuning, the preliminary segmentation was performed without the reversed ASM step. The resulting nVTCs of the segmented volumes shown in Fig. 9 are significantly different from those in Fig. 6(c). Note that the segmentation tightly follows the mean motion pattern not allowing for subject-specific differences. Consequently, disease patterns may be masked by this strong influence of the underlying motion pattern.
Table III lists the classification performance achieved using modal indices of 4-D shapes, nVTCs and dV/dt curves derived from ventricular volumes of different stages. The 4-D shape model analyzes both ventricular shape and motion variations and captures the difference between normal and TOF ventricles even for the LV whose shape is often not severely affected by TOF. The shape model analysis of the final segmentation achieved almost the same classification performance as that of the independent standard.
When the ventricular motion patterns are partially described by the nVTCs, the difference between normal and TOF hearts begins to diminish. Lower classification performance was achieved by using nVTCs of the independent standard. The dV/dt curves measure the rate of cardiac contraction and expansion but exclude the ejection fraction knowledge. Using the independent standard, the dV/dt model achieved better classification performance than the nVTC model.
Fig. 10 shows the shape variations at phases 0 and 7 introduced by varying the value of the first modal index of the LV+RV 4-D shape model derived from the independent standard. The first shape mode is associated with variations in LV myocardium thickness, LV orientation, bulging of apex, ventricular height along the long axis, and severity of RV enlargement. In this example, increasing the value of the first modal index introduces more TOF-related shape features to the heart. Very similar shape variation patterns were observed in the shape model derived from the final segmentation.
Fig. 11 shows the distributions of the first two modal indices of 4-D shape models of LV, RV, and LV+RV derived from the independent standard and the final segmentation. The strongest mode captures the main difference between normal and TOF hearts, which also agrees with the visualization in Fig. 10. In the shape model space, the normal and TOF hearts have clear differences in clustering of modal index distribution. Although the normal and TOF clusters have slightly different distribution patterns using models derived from the independent standard and the final computerized segmentation, such small differences have almost no effect on the classification performances listed in Table III.
Fig. 12 shows the nVTC variations associated with varying the first three modal indices within [−2σ, 2σ] together with the mean nVTCs. In spite of the difference in classification performance, the patterns of nVTC variations captured by the models derived from the independent standard and the final segmentation are almost the same. The first mode is mainly associated with different diastolic motion and the second mode with systolic motion. Detailed observations showed that the ranges of variations embedded in the final segmentations were slightly smaller than those in the independent standard.
Although our image fusion and customized tracing application reduced observer variability, the manual tracings are still subject to inevitable uncertainties and inconsistencies. The results produced by varying the base and apex locations identified on LA views by one or even two voxels can still be deemed as clinically plausible due to lack of certainty about location of the atrio-ventricular valves, but the resulting fluctuation of the ventricular volume is large. The ventricular border tracings on the SA slices are affected by inconsistent criteria used across different hearts and cardiac phases. Such observer variability also exists when comparing the computer segmentation with the independent standard and contributes to imperfect accuracy indicated by validation metrics such as the overlap ratio.
Representing 3-D shapes by landmarks always introduces loss of detail. Our quality testing showed that such a loss yields only very small landmarking errors and the important ventricular shape features are intact. Therefore, the landmarked surfaces provide sufficiently accurate representations of the ventricular shapes. An alternative landmark generation approach is to generate landmarks on the first phase of the 4-D datasets and then propagate the landmarks throughout the cardiac cycle. We did not use this approach because small geometric errors introduced by one propagation will be accumulated and potentially amplified in other cardiac phases.
One consequence of a small training set is that the model-synthesized shape could be significantly different from the independent standard because the model cannot fully capture the variability of the training set and cannot approximate data that are not in the training set. Our intensive quality testing revealed that the 3-D modeling has only slightly worse accuracy than the landmarking, while it achieves better accuracy, especially in motion tracking, than 4-D modeling. Because 4-D modeling has a limited motion tracking capability as shown in Fig. 6(c), using 3-D models for final segmentation seems preferable. Fig. 6(d) also shows the potential error patterns of the segmentations such as underestimation of LV EDV, overestimation of LV ESV and reduced LV mass at end-systole.
Another consequence of a small training set is the under-trained model matching that produces imperfect segmentation. The training set contains insufficient information to produce Jacobian matrices that can generate parameter displacement with correct sign and enough magnitude to guide the model matching. The model matching therefore can be easily trapped at local minima and lead to imperfect segmentation. This under-training problem is more serious in the preliminary segmentation as shown in Fig. 9. The preliminary segmentation benefits from the reversed ASM fine tuning that provides some extra momentum during the model matching, but it still cannot achieve good cardiac motion tracking since the full motion variation is not captured by the 4-D modeling. Although it has been shown that the ASM training can be improved by an artificially enlarged training set , applying a similar technique to our approach is not appropriate because the clinical validity of the training set is critical.
To achieve accurate segmentation and minimize side effects of a small training set by reducing shape variations to be included in the models, the LV and RV as well as normal and TOF hearts are modeled with anisotropic scaling and segmented separately. Manual initialization was chosen due to its robustness and efficient utilization of expert knowledge. The combination of manual initialization and preliminary segmentation is necessary in our approach because it utilizes learned motion patterns and reduces human effort in manual initialization, which is also made more intuitive when anisotropic scaling is used. Without preliminary segmentation, the subsequent final segmentation using 3-D phase-by-phase approach either requires more human interaction or does not offer a comparable performance. When constructing the shape model used for disease characterization, isotropic scaling is used to preserve all potential disease features. The possible very small intersections and holes in the septal region—because the septum is included in both LV and RV—were corrected during segmentation and we found that they have no influence on our disease characterization.
Our final segmentation uses the same single 3-D model of all cardiac phases. Using multiple 3-D models—one model for one or several phases—is possible in principle, but such approach uses a smaller training set than the single-model approach and tracking severe irregular ventricular motion [Fig. 8(a)] and its associated appearance variations [Fig. 7(b)] becomes more difficult. The single-model final segmentation achieved improved accuracy and motion tracking compared to the preliminary segmentation. However, the segmentation performance is still bounded by the inherent limitations of the classical model matching approach because it cannot guarantee that 1) the accurate segmentation corresponds to the minimum intensity-based matching error, 2) the simplified multi-variable optimization using Jacobian matrices will find the global minimum, and 3) different but equally plausible initializations will converge to the same solution. The final LV endocardium and RV segmentations are improved because temporal motion constraints that may be too stringent are excluded from the 3-D model. The segmentation of LV epicardial surface with weak intensity evidence suffers from losing constraints on myocardium thickness and results in temporal fluctuation of LV mass. Due to this problem with the epicardium, the reversed ASM fine tuning was not performed in the final segmentation.
The 4-D ventricular shape and motion are complicated and so is the validation, where a single metric only partially describes the accuracy. For example, a segmentation that is a slightly shifted version of the independent standard but accurately captures the actual shape and motion features will correspond to reduced overlap. Our validation was performed on a large set of metrics to provide a comprehensive interpretation of segmentation accuracy. We believe that the ultimate segmentation performance assessment is reflected by the performance of disease characterization—whether the segmentation captures the same disease features as the independent standard. Our results show that the reported segmentation errors have almost no effect on the distribution and classification of the 4-D shape modal indices. The volume-based models built from segmentation produce slightly different classification results than those from the independent standard, but Fig. 12 shows that they contain the same types of motion variation with only slightly smaller ranges of modal indices.
The current running time of a full segmentation process, 15–30 min, is comparable to that of a reported 3-D ASM using fuzzy interface , if it is used for full cardiac cycle segmentation. The segmentation can be further optimized with a multi-resolution approach and the final segmentation can benefit from parallel computing. With the latest development in hardware, the human interaction can further be streamlined by utilizing virtual reality . It is expected that the fully optimized and integrated application will need less than 10 min to produce comprehensive quantitative morphology and function results for clinical analysis.
In this paper, a complete cardiac MR image processing pipeline was designed and implemented, from creating the independent standard to disease characterization. The pipeline was carefully designed to meet the requirement of clinical applicability where many practical limitations must be overcome.
Based on knowledge learned from quality assessment of modeling and model training, the segmentation was designed as a two-stage process with goal-oriented modeling and segmentation strategies. The preliminary segmentation utilizes a 4-D hybrid ASM/AAM algorithm with special reversed-ASM fine tuning. The final segmentation step was designed to achieve better cardiac motion tracking as well as improved shape details. The human intervention was minimized to manual initialization on only one cardiac phase. The comprehensively assessed segmentation accuracy showed that even when trained on a small population, the AAM-based segmentation method achieved very good performance.
In the final stage, a set of novel shape and motion features were designed using 4-D shape model and volume-based curves. The normal/TOF classification of these features achieved 90%–100% sensitivity and specificity using both segmetation result and independent standard, which further proves the high quality of the segmentation. We also showed that the features extracted from the independent standard and from the segmentation result represent very similar disease characteristics in shape and motion. In addition, the TOF hearts exhibit less clustered distributions in the feature space, suggesting these features can potentially be used as disease status and/or disease progression indicators.
In conclusion, the left and right ventricles of normal and TOF hearts can be accurately segmented by our AAM-based method, and the extracted shape model and volume-based features can correctly identify the unique characteristics of TOF and potentially its progression.
The contribution of M. Thomas during preparation of the independent standard is gratefully acknowledged.
This work was supported in part by NIH NHLBI and NIH NIBIB grants R01-HL071809 and R01-EB004640.
Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org.
Honghai Zhang, Department of Electrical and Computer Engineering, The University of Iowa, Iowa City, IA 52242 USA (Email: ude.awoiu@gnahz-iahgnoh)
Andreas Wahle, Department of Electrical and Computer Engineering, The University of Iowa, Iowa City, IA 52242 USA (Email: ude.awoiu@elhaw-saerdna)
Ryan K. Johnson, Division of Pediatric Cardiology, Department of Pediatrics, The University of Iowa, Iowa City, IA 52242 USA ; Email: ude.awoiu@nosnhoj-nayr.
Thomas D. Scholz, Division of Pediatric Cardiology, Department of Pediatrics, The University of Iowa, Iowa City, IA 52242 USA ; Email: ude.awoiu@zlohcs-samoht.
Milan Sonka, Department of Electrical and Computer Engineering, the Department of Ophthalmology and Visual Sciences, and the Department of Radiation Oncology, The University of Iowa, Iowa City, IA, 52242 USA (Email: ude.awoiu@aknos-nalim)