## Level Set Methods Incorporating Generic Constraints

We have developed a coupled surfaces approach for automatically segmenting a volumetric layer from a 3D image [43-45]. This approach uses a set of coupled differential equations, with each equation determining the evolution or propagation of a surface within a level set framework. In the case of the cortex, one surface attempts to localize the white matter/gray matter (WM/GM) inner cortical boundary and the other the gray matter/cerebrospinal fluid (GM/CSF) outer boundary. Coupling between the surfaces incorporates the notion of an approximately fixed thickness separating the surfaces everywhere in the cortex. This soft constraint helps in ensuring that the GM/CSF boundary is captured even in the deeper cortical folds in the brain. A further assumption is that across each surface there is a local difference in the gray-level values, while in between the two surfaces there is a homogeneity of gray levels. By evolving two embedded surfaces simultaneously, each driven by its own image-based information, while maintaining the coupling, we are able to achieve an automatic and robust segmentation of the cortex, and simultaneously obtain a representation of the inner and outer cortical surfaces.

FIGURE 5 Platform for initializing shapes within subcortical structures [31]. See also Plate 8.

5.1 Level Set Method

For our purposes, the volumetric layer is defined completely by its bounding surfaces and the homogeneity in between. The essential idea is to first represent each cortical surface (in this case) as a front y(t). Given a moving closed front y(t), the idea is to produce an Eulerian formulation for the motion of this surface propagating along its normal direction with speed F, where F can be a function of geometric surface characteristics and image characteristics (such as those described earlier). Such a formulation can be produced by embedding the propagating interface as the zero level set of a higher dimensional function ¥ defined by ¥(x, t) = d, where d is the signed distance from position x to y(t). The equation of evolution of ¥, inside which our surface is embedded as the zero level set, is given by ¥t + F|V¥| = 0. As ¥(x, t) evolves, it remains a function. However, the propagating hypersurface y(t) may change topology, break, merge, and form sharp corners as the function ¥ changes. This flexibility greatly assists in the convergence to complex structures. In addition, the level set approach facilitates structural measurement through the direct computation of geometric parameters, such as curvature and thickness, via the level set function. For example, the intrinsic geometric properties of the front may be computed directly and accurately from the level function ¥. At any point on the front, the outward normal vector is given by n = V¥, and the mean curvature is easily obtained from the divergence of the gradient of the unit normal vector to front, that is, K = V • (V¥/|V¥|).

5.2 Coupled Surfaces Propagation: Speed Term Design

In solving the problem of segmenting the cortical gray matter layer bounded by two surfaces, we consider two moving interfaces describing the bounding surfaces. Starting from inside the inner bounding surface (gray/white boundary), with an offset in between, the interfaces propagate along the normal direction, stopping at the desired place, while maintaining the distance between them. Embedding each surface as the zero level set in its own level function, we have two equations, ¥int + Fin|V¥in| = 0 and ¥outt + Fout|V¥out| = 0, where Fm and Fout are functions of the surface normal direction, image-derived information, and the distance between the two surfaces. The coupling is embedded in the design of Fin and Fout. Where the distance between the two surfaces is within the normal range for cortical thickness, the two surfaces propagate according to the image-based information; where the distance between the two surfaces is out of the normal range, the distance constrains the propagation. We define

Function g smoothly maps larger gray level transition probability (Eq. 1) to slower speed, i.e., as the probability gets larger, g tends to zero, while as the probability approaches zero, g tends to a constant. Function h smoothly penalizes the distance outside of the normal range. As the distance between the surfaces goes out of normal range, h goes to zero, while within the normal range, h is constant. Thus, each surface moves with constant speed along the normal direction, and slows down or stops when either the image-based information becomes strong or the distance to the other surface moves away from the normal range. Surface renderings showing the inner and outer cortical layers found from a normal 3D brain MR image using this approach are shown in Fig. 6. A 2D example showing the delineation of the endocardium and epicardium of the left ventricle from an MR heart image is shown in Fig. 7. Here, the thickness constraint is much wider to allow for the larger variation in the thickness of the myocardium.

### 5.3 Validation

For algorithm testing with ground truth information, we used brain images from an MRI simulator [22] developed at the Montreal Neurological Institute (MNI) [28]. The images are generated using a labeled source image that allows users to independently control various acquisition parameters and obtain realistic MR images. The ground truth of the phantom is provided in the form of membership functions for each voxel belonging to different tissue types, such as the skull, CSF, gray matter, and white matter.

We tested our cortical coupled surface segmentation algorithm on simulated Tl-weighted images of a normal brain from the MNI simulator. Starting from the unedited images, user interaction is required to initialize the concentric spheres within the white matter. The spheres grow out and automatically lock onto the inner and outer cortical surfaces. As long as the spheres are placed inside the white matter, the algorithm is robust to starting position. To evaluate the segmentation result we defined the following measures. TP rate for a particular tissue type T is defined as the percent of phantom defined tissue voxels that were found by the algorithm. FP rate is defined as the percentage of algorithm defined tissue voxels that were not defined by the phantom. We also define the volume ratio to be the volume of all the voxels segmented as of tissue type T by our algorithm to the total partial volume of tissue type T specified by the phantom ( partial volume voxels contribute proportionally). Table 1 shows our measurement results over four types: total brain tissue (including white and gray matter), cortical gray matter in selected slices, and frontal white matter. Since the algorithm is designed specifically for the nearly constant thickness of the cerebral cortex, it recovers only part of the gray matter in the brain stem and the cerebellum where the constant thickness constraint is not satisfied. These regions account for most of the errors in the TP rate and volume ratio for the whole brain tissue. We compare the cortical gray matter volume on 49 frontal and 56 coronal

FIGURE 6 Results of coupled level set cortical segmentation. (a) Initialization of pairs of concentric spheres in 3D MR brain images (frontal part); (b) intermediate step; (c) final result of the outer (top) and inner (bottom) cortical surfaces of the frontal lobe; (d) single vs coupled surfaces approach. (Upper) Surfaces resulting from single surface approach shown on a sagittal image (finding the inner and outer cortical surfaces separately). (Lower) Results from the coupled surfaces approach run on original 3D data overlaid on a sagittal slice of the expert tracing result. The outer cortical surface resulting from the coupled algorithm nicely fits the boundary from the expert tracing. Coupling prevents the inner surface from collapsing into CSF(*1) and the outer surface from penetrating nonbrain tissue(*2). See also Plate 9.

FIGURE 6 Results of coupled level set cortical segmentation. (a) Initialization of pairs of concentric spheres in 3D MR brain images (frontal part); (b) intermediate step; (c) final result of the outer (top) and inner (bottom) cortical surfaces of the frontal lobe; (d) single vs coupled surfaces approach. (Upper) Surfaces resulting from single surface approach shown on a sagittal image (finding the inner and outer cortical surfaces separately). (Lower) Results from the coupled surfaces approach run on original 3D data overlaid on a sagittal slice of the expert tracing result. The outer cortical surface resulting from the coupled algorithm nicely fits the boundary from the expert tracing. Coupling prevents the inner surface from collapsing into CSF(*1) and the outer surface from penetrating nonbrain tissue(*2). See also Plate 9.

FIGURE 7 From left to right, evolving coupled level set for heart segmentation showing initialization, through intermediate and final results on the epicardium (red) and endocardium (green). See also Plate 10.

FIGURE 7 From left to right, evolving coupled level set for heart segmentation showing initialization, through intermediate and final results on the epicardium (red) and endocardium (green). See also Plate 10.

slices where there is only white matter and cortical gray matter. These results show that our algorithm performs well in isolating the brain from non-brain tissues and in segmenting the cortex from simulated data.

### 5.4 Images from Human Subjects

We have applied our cortical segmentation to true image data in an effort to further evaluate our approach. A quantitative comparison of our algorithm's results and those obtained by expert manual tracing was performed using a small set (N = 14) of studies. We compared either whole brain volume or cortical gray matter volume using our algorithm with that from the expert tracing. We found an average true positive rate of 95.2% for whole brain, a false positive rate of 3.9% for whole brain, and a true positive rate of 86.6% for cortical gray matter. These rates are very good and could be improved even further by adding the missing part of the brain stem and cerebellum by hand editing within minutes.

5.5 Cortical Shape, Thickness, and Surface Area

A large battery of measurements have been recorded, including measures related to both cortical and subcortical volume and shape. Figure 8 shows the inner and outer cortical surfaces of a frontal lobe colored with their shape indices [21], enabling automatic identification of gyri (mostly ridges) and sulci

0 0