This wiki pages describes the research effort for improving the subcortical segmentation accuracy for data acquired across different scanner platforms or imaging sequences.

The subcortical segmentation method within FreeSurfer is an atlas-driven method. It uses a probabilistic brain atlas to encode prior intensity models for individual brain structures at every atlas location, as well as prior models for the spatial distribution of these brain structures. As a result, the accuracy of the method often degrades when processing data acquired using scanner plaforms and/or imaging sequences that differ from the data used for the atlas training, due to changes in underlying image contrast. To resolve this problem and reduce the sensitivity of the method to changes in scanner platforms or acquisition parameters, we have developed an intensity renormalization procedure that automatically adjusts the prior intensity model of the atlas, in order to better fit it to the input data.

The intensity renormalization is performed for individual structures. The basic idea is to try to get an estimation of the intensity distribution of each structure in the input image, and then the prior atlas intensity model can be properly normalized to match the given data for each structure.

The intensity renormalization consists of a multi-linear image registration and a histogram matching step. It starts with a global linear registration that grossly aligns the atlas to the input image as in the original method. This global linear registration has only 12 degrees of freedom, and thus the registration accuracy is much less affected by cross-sequence image contrast changes.

Next, a local linear registration is being sought for each individual struture. The registration result will help define a (rough) delineation of the particular structure on the input image, so that the intensity distribution of the input image for that structure can be estimated. Since the linear registration in this step is computed for each structure independently (although all use the global linear registration as the starting point), the final registration result is not a single linear transform, and thus comes the name "multi-linear" registration. For the local registration of an individual structure, a template (sub-)image for it is first constructed from the atlas by building the maximum likelihood volume for the structure. The maximum likelihood volume is simply a gray-scale image where the intensity value at each voxel location is equal to the mean intensity of the maximum-likely tissue label (the label with the highest prior probability) at that location. The maximum likelihood volume for a particular structure is built in the same way but all the other labels are considered as background (intensity value set to zero) except the structure under study. Note that only non-background voxels (intensity value > 0) are considered in the image regitration algorithm. In practice, we also keep a small neighborhood of non-zero values around the given structure to provide sufficient contextual information to help the registration.

The matching criterion for the local registration is a variant of the standard SSE criterion (which is simply sum-of-squared-error or sum-of-squared-intensity-differences). From our problem setting, it is clear that the matching criterion has to take into account of possible image contrast difference. As a result, we modified the original sum-of-squared-error criterion to allow a linear mapping between the source and template image intensities. The criterion is similar to what is used in the following paper:

S. Periaswamy and H. Farid, "Elastic Registration in the Presence of Intensity Variations", IEEE T. Med. Imag., 22(7):865-874, 2003

Of cause, other criteria such as the normalized correlation ratio can also be used, but the optimization process will become more complicated and time-consuming.

After the matching criterion or objective function is defined, many different optimization algorithms can be applied to find the optimal linear registration mapping that optimizes (minimize in our case) the objective function. One example is the Powell's method, which aims to find a global optimal solution. Unfortunately, we found this method not robust enough for the particular problem. We currently settled on a gradient-descent optimization approach, which finds a local optimal solution around the initial global linear mapping.

After the registration is computed for a particular structure, the maximum likely label of the atlas is then transformed to the input image to define a mask for that structure. At this stage, the structure mask or segmentation cannot be very accurate, and thus the intensity normalization will rely on global features of the intensity distribution. In particular, two intensity histograms are built for the structure under consideration, one from the masked region in the input image, one from the atlas (using the intensity value of the maximum-likelihood volume for that structure). Next, a linear intensity transformation (a scaling plus an offset) is sought that maximizes the cross-correlation between the two histograms (currently, we force the offset to be zero and thus only optimize for the unknown scaling factor). The intensity transformation is then applied to the intensity model of the given structure at each atlas location.

After the intensity normalization is finished for every structure, the new atlas created is then applied in the sequential steps of the subcortical segmentation procedure as in the original method. In a newer version, we also use the multi-linear registration results to initialize the non-linear atlas morphing instead of using the original global linear registration. The advantage of this new initilization approach is not conclusive, however.

We note that the key component for the atlas intensity renormalization procedure is the local structure-wise linear registration. This type of part-to-all registration is much harder than the whole image registration, and thus robust and accurate results are hard to obtain. To take into account the potential difficulties, we deliberately reduce the degree-of-freedom in the problem definition. For example, as mentioned earlier, in the histogram matching, we only search for a scaling factor instead of both a scaling and an offset. One good thing is that the histogram matching based scaling factor computation is not sensitive to registration errors.

Based on the same consideration, the local linear registration is computed for a subset of structures (small structures like vessels are not considered; see the program for details). Based on the results from these structures, we can compute an average normalization factor for each of three major tissure classes (WM, GM, and CSF), and the average factor is assigned to the remaining structures.

Another detail not mentioned earlier involves some constraints we added to detect and get rid of obviously wrong registration or histogram matching results. For example, a registration result is dropped is the determinant of the registration matrix is too small or too big (to avoid unreasonable expansion or contraction of the structure). Also, a histogram matching result is dropped if the overlap of the final histograms is low. When this happens, the intensity normalization factor for the failed struture is set according to the average normalization factors computed using other structures. For exmaple, the average GM normalization factor will be assigned to hippocampus if the normalization for hippocampus initially failed.

Future directions:

1. More reliable structure-to-image registration method can be sought. One approach is to incorporate prior models for admissible structure registration, as a resort to improve the robustness of the registration method. The 1997 NeuroImage paper by Ashburner et al is a good reference for such type of approach.

2.