consistent sulcal parcellation of longitudinal ... - Semantic Scholar

Report 2 Downloads 95 Views
CONSISTENT SULCAL PARCELLATION OF LONGITUDINAL CORTICAL SURFACES Gang Li, Yang Li, Yaping Wang, Dinggang Shen Department of Radiology and BRIC, University of North Carolina at Chapel Hill ABSTRACT Automatic consistent sulcal parcellation of longitudinal cortical surfaces is of great importance in studying morphological and functional changes of human brains. This paper proposes a novel energy-function based method for consist sulcal parcellation of longitudinal cortical surfaces. Specifically, both spatial and temporal smoothness are imposed in the energy function to obtain consistent longitudinal sulcal parcellation results. The proposed method has been successfully applied to sulcal parcellation of longitudinal inner cortical surfaces of 10 normal subjects, each with 4 scans acquired within 2 years. Both quantitative and qualitative evaluation results demonstrate the validity of the proposed method. Index Terms— sulcal parcellation, segmentation, sulcal segmentation, graph cuts

longitudinal

1. INTRODUCTION Automatic sulcal parcellation of cortical surfaces is of great importance in structural and functional mapping of human brains. Studying longitudinal change of cortical surfaces, which is important to normal development, aging, and disease progression of human brains, requires accurate cortical surface parcellation, since longitudinal cortical changes are normally very subtle, especially in aging and Alzheimer’s disease. However, existing sulcal parcellation or segmentation methods were mainly designed for working on a single cortical surface [1, 2, 3, 4]. Therefore, applying these methods independently to the cortical surface of each time point in a longitudinal study might generate longitudinal inconsistent results. To overcome this limitation, in this paper, we propose a novel method for robust and consistent sulcal parcellation of longitudinal cortical surfaces of human brains. In the method, consistent sulcal parcellation is formulated as an energy minimization problem, in which both spatial and temporal smoothness is imposed to the sulcal parcellation results. Then, the energy function is efficiently minimized by a graph cuts method [5], which can guarantee to achieve a strong local minimum for certain energy functions. The method has been successfully applied to longitudinal inner cortical surfaces of 10 normal subjects, each with 4 scans acquired within 2 years. Both quantitative and qualitative evaluation results demonstrate the validity of the proposed method. In contrast to the existing sulcal parcellation or segmentation methods [1, 2, 3, 4], the advantages of the proposed method are as follows. First, our method is very general, since it can work on both a single cortical surface and longitudinal cortical surfaces, whereas the existing watershed or flow tracking based methods are not straightforward to extend to the longitudinal cortical surfaces.

978-1-4244-4128-0/11/$25.00 ©2011 IEEE

337

Second, in our method, a clear energy function is formulated for the sulcal parcellation problem, and the energy function is efficiently minimized by a graph cuts method, whereas no clear energy function is defined in the existing methods. Third, a spatial smoothness term and a temporal smoothness term are included in the energy function to conveniently control the smoothness and consistency of segmented sulcal boundaries, whereas it is not easy to do this in existing methods. Finally, our method is very flexible to integrate other cortical attributes, such as fiber density.

2. METHODS Given the reconstructed cortical surfaces at successive time points of the same subject, we aim to consistently parcellate these longitudinal cortical surfaces into a set of sulcal regions (which are the buried cortical regions surrounding the sulcal space on cortical surfaces [2, 3, 4]), and the corresponding sulcal basins (which are the cortical regions bounded by adjacent gyral crest lines on cortical surfaces [1, 3, 4]). Figure 1(a) shows a two dimensional (2D) schematic illustration of sulcal regions and sulcal basins. As we can see from Figure 1, sulcal regions are relatively easy to be identified in contrast to sulcal basins. Therefore, in our proposed method, sulcal regions are first consistently extracted from longitudinal cortical surfaces. And then, based on segmentation results of sulcal regions, longitudinal cortical surfaces are further parcellated as a set of corresponding sulcal basins. Both sulcal region and sulcal basin segmentations are formulated as an energy minimization problem, and a graph cuts [5] method is used to efficiently solve the energy minimization problem.

2.1. Problem formulation Both consistent sulcal region and sulcal basin segmentations on longitudinal cortical surfaces can be considered as a discrete labeling problem. For sulcal region segmentation, the problem is to assign a label indicating sulcal or gyral region to each vertex. For sulcal basin segmentation, the problem is changed to assign a label indicating a sulcal basin to each vertex in gyral regions. Meanwhile both the spatial and temporal contextual information of adjacent vertices are taken into account when performing labeling. Therefore, the general energy function of sulcal region and sulcal basin segmentations is formulated as: (1) E Ed  Os Es  Ot Et where Ed is the data term, and Es and Et are spatial and temporal smoothness terms, respectively. Os and Ot are the weighting parameters for spatial and temporal smoothness terms, respectively. The data term represents the sum of the set of data cost per vertex:

ISBI 2011

Figure 1. (a) A 2D schematic illustration of sulcal regions and sulcal basins. The red colors indicate the two sulcal regions. The blue and green colors represent two adjacent sulcal basins. (b) An example of the maximum principal curvature on an inner cortical surface. The color bar is on the right. (c) The histogram distribution of the maximum principal curvature on the cortical surface in Figure 1(b).

¦ Dx (lx )

Ed

(2)

x

where Dx (lx ) indicates the cost of labeling a vertex x as l x , which indicates a possible label and will be clarified in next sections. The spatial smoothness term, representing the sum of the cost of labeling a pair of neighboring vertices in the cortical surface, is used to impose spatial smoothness on the segmentation results, and is defined as: (3) E V s (l , l )

¦

s

{ x , y}N s

x, y

x

y

x

where N s represents the set of neighboring vertex pairs in the cortical surface. In this paper, the one ring neighborhood in the triangular surface mesh is adopted. The temporal smoothness term, representing the sum of the cost of labeling a pair of temporal neighboring vertices in cortical surfaces of two successive time points, is used to impose temporal smoothness, and is defined as: (4) E V t (l , l ) t

¦

{ x , y}N t

negative values at sulcal bottoms and the large positive values at gyral crests on cortical surfaces [4]. Figure 1(b) shows an example of the maximum principal curvature map on a cortical surface. And Figure 1(c) shows the histogram of maximum principal curvatures of the cortical surface in Figure 1(b). The distributions of maximum principal curvatures in the sulcal and the gyral regions can be modeled as respective Gaussian functions, as used in [4]. In principal, when labeling a vertex as the region (sulcal or gyral region) which it belongs to, the cost should be small. Therefore, for sulcal region segmentation, Dx (l x ) is defined as: (5) Dx (l x ) 1.0  exp(  (c( x)  ml ) 2 2T l2 )

x, y

x

y

where N t represents the set of neighboring vertex pairs between two successive longitudinal cortical surfaces. In this paper, we focus on aging brains. Since the longitudinal change for the aging brains in 6 months is generally very subtle especially for the inner cortical surface, the temporal neighbors are defined as the closest vertices along the normal direction of the current vertex in the two neighboring cortical surfaces immediately before and after the current time point. When Ot is set as 0, the above energy function can be used for sulcal parcellation on a single cortical surface. The definitions of the data term and the spatial and temporal smoothness terms are quite different for sulcal region and sulcal basin segmentations, which will be detailed in sections 2.2 and 2.3, respectively.

2.2. Consistent sulcal region segmentation Before introducing the detailed terms in the energy function, the attributes on cortical surfaces used for sulcal parcellation in this paper are list as follows: the maximum principal curvature c(˜) (which is the principal curvature with the larger magnitude in the two principal curvatures at the vertex on cortical surfaces), the maximum principal direction p(˜) (which is the direction corresponding to the maximum principal curvature), and the normal direction n(˜) . To distinguish sulcal and gyral regions, we adopt the maximum principal curvatures, which are the large

338

x

where lx {0, 1} indicate sulcal and gyral regions, respectively. ml x and T l represent the mean and standard deviation of x

c(˜) in the

region l x . Thresholding the maximum principal curvature map using zero is adopted to obtain the initial sulcal and gyral region segmentation results and estimate the parameters m and T for the sulcal and gyral regions. The spatial smoothness term Vxs, y (l x , l y ) is defined as:

Vxs, y (lx , l y ) 1  G ( lx  l y )

(6)

where G is the Dirac delta function. The temporal smoothness term Vxt, y (l x , l y ) is defined as:

Vxt, y (lx , l y )

0.5 u (1  (n( x) ˜ n( y))) u (1  G ( lx  l y )

(7)

The central idea behind this setting is to force the cost of discontinuous labeling for a pair of temporal neighboring vertices with dissimilar normal direction to be small. The alpha-expansion graph cuts method [5] is adopted to efficiently solve the above energy minimization problem for sulcal region segmentation. After sulcal region segmentation, we perform connective component analysis to label each connective sulcal region as a unique value [4], which will be further used for sulcal basin segmentation as described below. Figure 2 (a) shows the longitudinal sulcal region segmentation results on a subject.

2.3. Consistent sulcal basin segmentation With the consistent sulcal region segmentation results, which provide rough location of sulcal basins, we further consistently parcellate longitudinal cortical surfaces into sulcal basins. Assuming that we have N isolated sulcal regions and each sulcal region has a corresponding sulcal basin, we aim to assign a sulcal basin label to each vertex on the gyral regions. One intuitive thought is to assign the sulcal basin label based on the closet geodesic distance to the N sulcal regions; however, due to the

(a)

(b)

(c)

Figure 2. An example of sulcal parcellation results on longitudinal cortical surfaces with 4 successive time points of a subject. The top row images show the longitudinal sulcal region segmentation results. The middle row images show the corresponding longitudinal sulcal parcellation results by the proposed method. The bottom images show sulcal parcellation results by the flow tracking method [4]. White curves indicate the boundaries of sulcal basins. Red arrows indicate several locations where our results are more consistent. non-symmetric structure of gyri, as shown in Figure 1(a), the generated boundaries of sulcal basins by this way might not be on the gyral crests. To deal with this problem, for sulcal basin segmentation, the data term is defined based on the weighted geodesic distances away from the segmented sulcal regions as: (8) Dx (lx ) 1.0  exp(E ˜ gr ( x)) x

where g r (x) is the weighed geodesic distance between vertex x

x

and sulcal region rx corresponding to sulcal basin lx  {0, ..., N - 1} . E is a weighting parameter and is set as 0.2 experimentally.

g rx (x) is calculated using the fast marching method on surface meshes and the marching speed for each vertex is defined as: F ( x) exp(  c( x) )

(9)

The central idea is that: given a sulcal region, if a vertex is not in the sulcal basin corresponding to the current sulcal region, the geodesic path from the vertex to the current sulcal region will pass through gyral crests, where marching speed are designed as small values, therefore, the weighted geodesic distance between the vertex and the current sulcal region will be large; otherwise, the weighted geodesic will be small. From each of the N segmented sulcal region, where geodesic distances are set as 0, the fast marching method is conducted to calculate the weighted geodesic distance away from the current sulcal region for each vertex on gyral regions. Therefore, each vertex on gyral regions has N weighed geodesic distance away from N sulcal regions. The spatial smoothness term Vxs, y (l x , l y ) is defined as:

Vxs, y (lx , l y )

w( x, y) u (1  G ( lx  l y ))

(10)

w( x, y) 0.5 u (1  p( x) ˜ p( y)) u (exp(c( x))  exp(c( y))) / 2 (11) where w( x, y) represents th e weight between two spatial

339

neighboring vertices. Note that the maximum principal direction p(˜) has been forced to point to the decreasing direction of c(˜) using the method in [4]. The basic idea is that two neighboring vertices belonging to two different sulcal basins normally meet at gyral crests where the maximum principal curvatures are the large positive values and the maximum principal directions point away from each other [4]. Therefore, the cost of discontinuous labeling is small around gyral crests. Herein, the temporal smoothness term is defined the same as the one for sulcal region segmentation. The alpha-expansion graph cuts method [5] is adopted to efficiently solve the above energy minimization problem for sulcal basin segmentation. The sulcal basin segmentation generates a complete parcellation of cortical surfaces, where each vertex has a unique sulcal basin label. The boundaries of sulcal basins can be considered as gyral crest lines. Figure 2 (b) shows the longitudinal sulcal basin segmentation results on a subject.

3. RESULTS The method is applied to 10 normal healthy subjects, in which each subject has been scanned 4 times with the interval of 6 months. The method CLASSIC [6] is used to obtain the consistent tissue segmentation results. All cortical surfaces are generated by the BrainSuite software [7]. For sulcal region segmentation, the parameters Os and Ot are set as 0.5 and 0.5, respectively. And for sulcal basin segmentation, the parameters Os and Ot are set as 0.8 and 0.8, respectively. The optimal setting of parameters will be investigated in future. All the longitudinal cortical surfaces are consistently segmented into anatomically meaningful sulcal regions and sulcal basins by the proposed method. Figure 2 shows an example of sulcal parcellation results on a subject with 4

successive time points. In this figure, we also compare the parcellation results with the results by the flow tracking method [4], which was developed for sulcal parcellation on a single cortical surface. As we can see, the sulcal parcellation results with the proposed method are more consistent and reasonable, since both spatial and temporal smoothness are imposed in the proposed method. To further evaluate the consistency of the sulcal parcellation results, we calculate the areas of left central sulcal basins on 10 subjects across 4 time points as shown in Figure 3. Figures 3(a) and 3(b) show the results obtained from the flow tracking method and the proposed method, respectively. Since the longitudinal change of inner cortical surfaces in 6 months for aging is normally very subtle, we expect that the areas should be consistent across time points. As we can see, the results generated by the proposed method are much more consistent. We also define a distance measurement to quantitatively evaluate the performance. Let the boundaries of sulcal basins in two successive surfaces be denoted as S1 and S 2 , the distance is defined as: d

0.5 u (

1 1 min S 2 (i)  S1 ( j ) ) ¦ min S1 (i)  S2 ( j )  m i¦ jS1 n iS1 jS2 S2

(12)

where n and m are the total numbers of points in S1 and S 2 , respectively. Figure 4 shows the distances of parcellated sulcal basin boundaries on whole cortical surfaces between each pair of successive time points on 10 subjects. Since the longitudinal change is very subtle, we expect small distance. As we can see, the distance from our results is smaller than that by the flow tracking method, indicating the consistency of the proposed method. Since our sulcal parcellation results are very consistent on longitudinal cortical surfaces, to validate the accuracy of the proposed method, we manually label the left central, postcentral and superior temporal sulcal basins on the first time point of 10 subjects and calculate the above distance error between automated and manual labeled boundaries of sulcal basins. Figure 5 shows the comparison results. Overall, the average distance error is around 1.0 mm, indicating the relative accuracy of the proposed method.

Figure 4. The distance between boundaries of parcellated sulcal basins at each pair of successive time points on whole surfaces.

Figure 5. The distance error of left central, postcentral and superior temporal sulcal basin boundaries on 10 subjects, compared to the manually labeling results.

4. CONCLUSION This paper presents a novel energy-function based method for consistent sulcal parcellation of longitudinal cortical surfaces. The method has been applied to the longitudinal images of 10 subjects and achieves good results. Currently, the longitudinal cortical surfaces are independently reconstructed from each image. In the future, we will investigate consistent reconstruction of longitudinal cortical surfaces and obtain more accurate vertex correspondence information, thus the longitudinal sulcal parcellation results could be further improved potentially. The future work also includes the further validation of our method using more longitudinal data, and also the analysis of longitudinal change of segmented sulci for aging and disease progression applications.

REFERENCES

(a)

(b) Figure 3. The areas of left central sulcal basins on the 10 subjects across 4 time points. (a) The results from the flow tracking method. (b) The results from the proposed method.

340

[1] G. Lohmann, and D.Y. von Cramon, “Automatic labelling of the human cortical surface using sulcal basins,” Med. Image Anal. 4 (3), 179-188, 2000. [2] M.E. Rettmann, X. Han, C. Xu, and J.L. Prince, “Automated sulcal segmentation using watersheds on the cortical surface,” NeuroImage 15 (2), 329-344, 2002. [3] F. Yang, and F. Kruggel, “Automatic segmentation of human brain sulci,” Med. Image Anal. 12 (4), 442-451, 2008. [4] G. Li, L. Guo, J. Nie, and T. Liu, “Automatic cortical sulcal parcellation based on surface principal direction flow field tracking,” NeuroImage 46 (4), 923-937, 2009. [5] Y. Boykov, O. Veksler, and R. Zabih, “Fast approximate energy minimization via graph cuts,” TPAMI, 23 (11), 1222-1239, 2001. [6] Z. Xue, D. Shen, and C. Davatzikos, “CLASSIC: consistent longitudinal alignment and segmentation of serial image computing,” NeuroImage 30 (2), 388-399, 2006. [7] D. Shattuck, and R. Leahy, “BrainSuite: an automated cortical surface identification tool,” Med. Image Anal. 8 (2), 129-142, 2002.