Introduction

One of the most popular HARDI reconstruction techniques is Q-ball imaging [1]. It's output is the local probability density function on a sphere, the orientation distribution function (ODF). HARDI data can be beneficial in fiber tracking, where DTI-based reconstructions suffer from weak directional information, e.g., in fiber crossings, fannings or kissings. However, with multiple orientations per voxel it is not clear which is the most appropriate one for a specific fiber bundle reconstruction. Therefore, fiber tracking techniques using HARDI data must distinguish between one and more fiber populations per voxel.

Anisotropy measures provide information in terms of fiber tract integrity and can act as a stopping criteria for fiber tracking. The generalized fractional anisotropy (GFA) [1] is an important criterion for HARDI, as an extension to the fractional anisotropy (FA) in DTI

GF A = std(Ψ ) rms(Ψ ) = √ n ∑ n i=1 (Ψ (u i ) − ⟨Ψ ⟩) 2 (n − 1) ∑ n i=1 Ψ (u i ) 2(1)

Here, ⟨Ψ ⟩ is the mean of the ODF and Ψ (u i ) is the ODF value of the i-th diffusion direction u.

In addition, different approaches exist to delineate the intra-voxel fiber population. Frank et al. [2] introduced the fractional multi-fiber index (1823-02) for determining the model order, l, of the current voxel FMI =

∑ j:l≥4 |c j | 2 ∑ j:l≥2 |c j | 2(2)

Here, c are the spherical harmonics (SH) coefficients, used for ODF reconstruction.

Another approach was introduced by Chen et al. [3] and Descoteaux et al. [4], which incorporates the variance of the measurements (in the following named Chen's classifier)

R 0 = |c 0 | ∑ j |c j | R 2 = ∑ j:l=2 |c j | ∑ j |c j | R multi = ∑ j:l≥4 |c j | ∑ j |c j | (3)

If R 0 is large, the voxel's diffusion is considered to be isotropic, if R 2 is large, a one-fiber population is present in the specific voxel. A large R multi -value indicates two or more fibers' diffusion.

In this paper, we introduce the morphological fiber classification (MFC) method, which is a global heuristic to differentiate between voxels with one or more fiber populations. We will demonstrate that our approach is advantageous in challenging cases, for example acquisitions with low b-values, where other classifier fail to detect the multiple maxima.

Materials and Methods

We performed our initial experiments on a phantom dataset. This phantom was originally provided by the Laboratoire de Neuroimagerie Assistée par Ordinateur (LNAO, France) for the Fiber Cup, a tractography contest at the MICCAI conference in 2009. A ground truth of the fibers in the phantom was provided as well. We use this ground truth to evaluate the results of our proposed approach in terms of multiple fiber populations per voxel. The phantom data was acquired with two repetitions and 64 image encoding gradients, uniformly distributed over a sphere. The two repetitions were averaged before further processing. Dataset size was 64 × 64 voxels with an uniform voxel size of 3 mm. Of the different diffusion sensitizations provided, we use the dataset with b-value 2000 s/mm 2 .

To be able to characterize voxels as containing zero (isotropic), one or multiple fiber populations, first their respective diffusion profile has to be calculated. We use the Q-ball reconstruction based on SH as proposed by Descoteaux et al. [5]. The regularization parameter for the Laplace-Beltrami smoothing matrix is λ = 0.006 and the employed SH order is l = 4. This order is high enough to classify multiple fiber populations in a voxel [2,3], and low enough to avoid over-modeling perturbations due to noise in the input diffusion MRI signal [5].

The main idea of MFC is to morphologically eliminate fibers from the white matter mask so that only clusters remain. These clusters represent an estimation of voxels with multiple fiber populations. A mask is generated in the first step and separates isotropic voxels from voxels with at least one fiber population. Ideally, a mask image separates all white matter voxels from all gray matter voxels. Since the white matter mask features gaps due to the thresholding procedure, the second step is to close these gaps. We apply morphological closing with different kernel sizes for filtering. With reasonable threshold choices, the gaps in the white matter mask that need to be closed are small. Closing consists of a morphological dilation (2 × 2 × 2) followed by an erosion (4 × 4 × 4). Since our goal is to eliminate fibers, the erosion is performed with a larger kernel size. This kernel eliminates the additional white matter voxels from dilation and thins out white matter. The third step is morphological opening with a kernel size of 3 × 3 × 3. The resulting image contains clusters located at positions where multiple fibers meet. In the fourth step the median filtered mask image is combined with the cluster image to form the final result. Voxels marked in both images are characterized as containing multiple fiber populations, whereas voxels marked only in the mask image as containing only one fiber population. Due to the dilation in the third step some clusters might have been enlarged beyond the actual white matter mask. Therefore voxels marked only in the cluster image are ignored.

Results

Figure 1a shows the ground truth of the phantom dataset, containing the most challenging fiber courses: crossing, kissing, and fanning. The result obtained from MFC with a standard deviation mask image is presented for comparison in Fig. 1b. MFC is able to identify all areas of multiple fiber populations.

For our experiments thresholded mask images are calculated based on GFA, FMI, Chen's method, and the standard deviation (sDEV) of the diffusion ODF. The respective thresholds are chosen carefully to find a proper balance between gaps in the mask and false positives. The computed masks are shown with their corresponding MFC in Figure 2. GFA is not able to separate the white and gray matter of this dataset (Fig. 2a). Hence, its MFC provides no useful information. FMI and Chen's method separate the white and gray matter well, but fail to detect multiple fiber populations per voxel in this dataset (Fig. 2b, 2c, top row). However, their respective MFCs perform significantly better, albeit fail to detect some multiple fiber population areas. Also, the MFC computed from Chen's mask image detects one false positive cluster (Fig. 2c). The best separation of white and gray matter was obtained from the standard deviation (Fig. 2d). Further, the corresponding MFC detects all multiple fiber areas with no false positives.

Our experiments show that the selection of proper thresholds in step one is crucial for the success of the MFC algorithm. Further, median filtering the mask image before or instead of step two to close gaps and eliminate false positives yields worse classification results.

For further evaluation we used a human brain dataset (dataset size 128 × 128 × 60, voxel size 1.875 × 1.875 × 2 mm), which is courtesy of Poupon et al. [6]. Data was acquired with a uniform gradient direction scheme with 41 directions and a b-value of 700 s/mm 2 . We tested our method in the region of the centrum semiovale, where a known crossing of the corpus callosum, the corticospinal tract and the superior longitudinal fasciculus exist. The resulting MFC with properly adjusted kernel sizes (Fig. 3).

Discussion

In neuro-visualizations the extraction of fiber tracts connecting functional areas of the brain is of major intrest. Our approach can pose an improvement in detecting voxels traversed by more than one fiber and hence, influence fiber tracking methods for HARDI. Additionally, the MFC can be used to distinguish isotropic diffusion from multiple intra-voxel fiber population and thus, provide information about fiber tract integrity. Future work will include the integration of the proposed classifier into a fiber tracking algorithm.