Abstract
In diffusion MRI, fiber tracts, represented by densely distributed 3D curves, can be estimated from diffusion weighted images using tractography. The spatial geometric structure of white matter fiber tracts is known to be complex in human brain, but it carries intrinsic information of human brain. In this paper, inspired by studies of liquid crystals, we propose tract-based director field analysis (tDFA) with total six rotationally invariant scalar indices to quantify local geometric structures of fiber tracts. The contributions of tDFA include: (1) We propose orientational order (OO) and orientational dispersion (OD) indices to quantify the degree of alignment and dispersion of fiber tracts; (2) We define the local orthogonal frame for a set of unoriented curves, which is proved to be a generalization of the Frenet frame defined for a single oriented curve; (3) With the local orthogonal frame, we propose splay, bend, and twist indices to quantify three types of orientational distortion of local fiber tracts, and a total distortion index to describe distortions of all three types. The proposed tDFA for fiber tracts is a generalization of the voxel-based DFA (vDFA) which was recently proposed for a spherical function field (i.e., an ODF field). To our knowledge, this is the first work to quantify orientational distortion (splay, bend, twist, and total distortion) of fiber tracts. Experiments show that the proposed scalar indices are useful descriptors of local geometric structures to visualize and analyze fiber tracts.
You have full access to this open access chapter, Download conference paper PDF
Similar content being viewed by others
1 Introduction
Diffusion MRI (dMRI) provides a unique window to non-invasively reveal anatomical connections (i.e., fiber tracts) and white matter tissue properties in human brain [3]. In dMRI, a typical processing pipeline before statistical analysis is: (1) fit a diffusion model (e.g., diffusion tensor imaging (DTI) [2]) to measured diffusion signals in each voxel; (2) calculate various voxel-wise scalar indices from the model parameters; (3) calculate local fiber directions in each voxel based on the model parameters, and then perform tractography. Scalar indices from diffusion models in step 2 are to quantify different tissue properties at different levels. For example, at the voxel level, fractional anisotropy (FA), generalized FA, return-to-origin probability, mean-squared displacement, and orientation dispersion index are indices for tissue properties inside a single voxel in DTI, Q-ball imaging [11], spherical polar Fourier imaging [6], NODDI [13], etc. At the local neighborhood level, gradient norm [8], sheet probability index [10], and orientational distortion [5] are indices for tissue properties within a local neighborhood.
Fiber tracts estimated by tractography are curves densely distributed in the 3D space. The spatial geometric structure of fiber tracts is known to be very complex in human brain [3]. In a local spatial region, typical geometric structures of fiber tracts include splay (aka diverge, converge, or fanning), bend, twist, crossing, kissing, etc. Grid and sheet structures were also reported in some areas in human brain [10, 12]. See Fig. 1 for a demonstration of splay, bend, and twist of fiber tracts. Compared with various voxel-wise indices, only a few tract-based indices have been proposed. Savadjiev et al. [9] proposed the total dispersion index for each point in a tract. Curvature and torsion based on the Frenet frame along the curve are two famous features of a single tract [1, 4]. However, the Frenet frame and its features are not designed for a set of curves, thus they cannot quantify local geometric structures of fiber tracts shown in Fig. 1.
Cheng et al. [5] proposed a framework called director field analysis (DFA) to quantify orientational order, dispersion of spherical functions, and orientational distortion (splay, bend, twist, and total distortion) of a spherical function field at the local neighborhood level. Since the DFA in [5] was developed for a spherical function field, we call it as voxel based DFA (vDFA). vDFA does not work for fiber tracts. In this paper, inspired by studies of liquid crystals [7], we propose tract-based director field analysis (tDFA) which generalizes vDFA in [5] to fiber tracts. tDFA is not proposed to replace vDFA, but adopts methods and concepts in vDFA to fiber tracts. In tDFA, we define total 6 scalar indices at each point in fiber tracts, where orientational order and dispersion quantify the alignment and dispersion of fiber tracts, and splay, bend, twist and total orientational distortion quantify local orientational distortion of fiber tracts. tDFA is applied directly on fiber tracts after tractography. To our knowledge, this is the first work to quantify orientational distortion of fiber tracts (i.e., splay, bend and twist).
2 Method
2.1 Directors, Oriented Curves and Unoriented Curves
A director, borrowed from studies of liquid crystals [7], is defined as a vector \({{\varvec{v}}}\) that is equivalent to its negative \(-{{\varvec{v}}}\) [5]. A director \({{\varvec{v}}}\) can be represented as a dyadic tensor \({{\varvec{vv}}}^{T}\) without sign ambiguity. A local fiber direction in a voxel (i.e., a peak of the ODF in that voxel) is a director by definition. “Sign ambiguity” also happens for fiber tracts. A fiber tract is mathematically a differentiable curve in \(\mathbb {R}^3\) with two endpoints. A differentiable curve can be parameterized as an oriented curve \(C: [0,1] \mapsto \mathbb {R}^3\), where C(t) is differentiable and C(0) and C(1) are two endpoints. If we set \(\gamma (t)=C(1-t)\) (i.e., inverse the orientation of the curve), then \(\gamma : [0,1] \mapsto \mathbb {R}^3\) is a new parameterization with the opposite orientation, i.e., \(\gamma (0)=C(1)\) and \(\gamma (1)=C(0)\). Note that these two oriented curves actually represent the same unoriented curve in geometry, however, for a point C(t) (i.e., \(\gamma (1-t)\)) in the curve, the two tangent vectors in these two parameterization have opposite directions, i.e., \(T_{C(t)}=-T_{\gamma (1-t)}\). Thus, tangent vectors of fiber tracts are actually directors.
The sign ambiguity of fiber tracts is not a problem for some cases. For example, the curvature and torsion do not change under the above two parameterization. However, it is indeed a crucial problem if we would like to compare two tangent vectors. Thus, we will use tools of directors in [5] for calculation to avoid this problem. The difference of two directors \({{\varvec{v}}}_1\) and \({{\varvec{v}}}_2\) in the director representation, denoted as \(\text {Diff}_d\), is \(\text {Diff}_d({{\varvec{v}}}_1, {{\varvec{v}}}_2)=s_1{{\varvec{v}}}_1-s_2{{\varvec{v}}}_2\), where \(s_i=\pm 1\) such that \(s_1s_2 {{\varvec{v}}}_1^T{{\varvec{v}}}_2 \ge 0 \) [5].
2.2 Orientational Order and Dispersion of Fiber Tracts
Considering a set of fiber tracts \(\{C_i\}\) which are curves densely distributed in \(\mathbb {R}^3\), for each point \({{\varvec{x}}}\in C_i\), inspired by liquid crystals [5, 7], we define orientational order index (OO) as
where \({{\varvec{y}}}\) is a point in a curve \(C_j\) in a spatial neighborhood \(\varOmega ({{\varvec{x}}})\) of \({{\varvec{x}}}\), \(\mathbf {u}_1({{\varvec{x}}})\) and \(\mathbf {u}_1({{\varvec{y}}})\) are unit norm tangent vectors of curves at \({{\varvec{x}}}\) and \({{\varvec{y}}}\), \(w({{\varvec{y}}},{{\varvec{x}}})\) is a spatial weighting function (e.g., uniform or Gaussian weighting) and \(\sum _{{{\varvec{y}}} \in \varOmega ({{\varvec{x}}})} w({{\varvec{y}}},{{\varvec{x}}})=1\). See Fig. 2. We have \(\text {OO} \in [-0.5,1] \). Then, we define the orientational dispersion index (OD) as \(\text {OD}=1-\text {OO}\). Thus, \(\text {OD} \in [0,1.5]\). If \(\mathbf {u}_1({{\varvec{y}}})\) is parallel to \(\mathbf {u}_1({{\varvec{x}}})\), \(\forall {{\varvec{y}}} \in \varOmega ({{\varvec{x}}})\), i.e., the least dispersion case, then \(\text {OO}=1\) and \(\text {OD}=0\) at \({{\varvec{x}}}\). If \(\mathbf {u}_1({{\varvec{y}}})\) is orthogonal to \(\mathbf {u}_1({{\varvec{x}}})\), \(\forall {{\varvec{y}}} \in \varOmega ({{\varvec{x}}})\), i.e., the most dispersion case, then \(\text {OO}=-0.5\) and \(\text {OD}=1.5\) at \({{\varvec{x}}}\). Note that the sign ambiguity of tangent vectors \(\mathbf {u}_1({{\varvec{y}}})\) and \(\mathbf {u}_1({{\varvec{x}}})\) dose not change the values of OO and OD. In practice, we normally set \(\varOmega _{{\varvec{x}}}\) as a ball centered at \({{\varvec{x}}}\) with a radius \(r\,=\,4\,\mathrm{mm}\), considering typical isotropic voxel size of a DW image is 2 mm. Note that the total dispersion in [9] is defined as the mean of differences between tangent vectors in a circle and \(\mathbf {u}_1({{\varvec{x}}})\), while our definition of OD is inspired by the order parameter in liquid crystals [5, 7].
Demonstration of tract-based DFA at a point \({{\varvec{x}}}\), denoted as the red point. We use tubes, instead of traditional arrows, to denote tangent vectors because of the sign ambiguity of the unoriented curves and tangent vectors. \(\varOmega _{{\varvec{x}}}\) is a local neighborhood of \({{\varvec{x}}}\), denoted as the grey sphere. The yellow points are within \(\varOmega _{{\varvec{x}}}\), and their tangent vectors are used to calculate OO and OD, and construct the second and third directors in the local orthogonal frame. The pink plane is the orthogonal plane of red director \(\mathbf {u}_1({{\varvec{x}}})\). The red, green, and blue tubes are the 3 directors \(\{\mathbf {u}_1({{\varvec{x}}}), \mathbf {u}_2({{\varvec{x}}}), \mathbf {u}_3({{\varvec{x}}})\}\) in the local orthogonal frame. The 6 purple points \(\{{{\varvec{x}}}\pm k\mathbf {u}_i\}\) alone these 3 directors in the local orthogonal frame are used to calculate the distortion indices.
2.3 Local Orthogonal Frame for a Set of Unoriented Curves
Curvature and torsion are typical features of a single curve based on the Frenet frame [1]. However, the Frenet frame has several limitations: (1) It requires a direction (i.e., orientation) of a curve; (2) It is defined only for a single curve, not for a set of curves; (3) For a straight line, the Frenet frame is not well defined because the curvature of straight lines is 0.
We propose a local orthogonal frame at each point \({{\varvec{x}}}\) of fiber tracts. The local orthogonal frame has three directors. The first director, denoted as \(\mathbf {u}_1({{\varvec{x}}})\), is the tangent vector at \({{\varvec{x}}}\), and the other two directors are in the orthogonal plane of the first director. We project all tangent vectors of fiber tracts in a local neighborhood \(\varOmega _{{\varvec{x}}}\) of \({{\varvec{x}}}\) onto the orthogonal plane, and define the covariance matrix of the projected directors as
Then we set the second director \(\mathbf {u}_2({{\varvec{x}}})\) as the eigenvector associated with the largest eigenvalue of \(\mathbf {Q}_{{\varvec{x}}}\), i.e., the principal component in PCA. In this way, \(\mathbf {u}_2({{\varvec{x}}})\) indicates the direction of the largest change of \(\mathbf {u}_1({{\varvec{x}}})\) in the orthogonal plane. Then, the third director is the cross product of the first two directors. See Fig. 2. We prove in Proposition 1 that when \(\varOmega _{{\varvec{x}}}\) tends infinitesimally small, the local orthogonal frame will converge to the Frenet frame.
Proposition 1
Let \(\{C_i\}\) be a set of differentiable curves. Let \({{\varvec{x}}} \in C_1\) be a point in \(C_1\), and \(\mathbf {u}_1({{\varvec{x}}})\) is the tangent vector at \({{\varvec{x}}}\) in \(C_1\). Assume \({{\varvec{x}}}\) is not a crossing point. Given a local ball \(\varOmega _{{\varvec{x}}}=\{{{\varvec{y}}} \ | \ \Vert {{\varvec{y}}}-{{\varvec{x}}}\Vert \le r \}\) with a radius r, let the constructed local orthogonal frame be \(\{\mathbf {u}_1({{\varvec{x}}}), \mathbf {u}_2({{\varvec{x}}}), \mathbf {u}_3({{\varvec{x}}})\}\). For \({{\varvec{x}}} \in C_1\), let the Frenet frame be \(\{\mathbf {u}_1({{\varvec{x}}}), {{\varvec{v}}}_2({{\varvec{x}}}), {{\varvec{v}}}_3({{\varvec{x}}})\}\). The vectors in two frames all have unit norm. Then we have \(\lim _{r \rightarrow 0} (\mathbf {u}_2({{\varvec{x}}})^T{{\varvec{v}}}_2({{\varvec{x}}}))^2=1\), \(\lim _{r \rightarrow 0} (\mathbf {u}_3({{\varvec{x}}})^T{{\varvec{v}}}_3({{\varvec{x}}}))^2=1\).
2.4 Orientational Distortions (Splay, Bend, and Twist)
Based on distortion analysis of liquid crystals [7], there are 3 types of orientational distortions as showed in Fig. 1. (1) splay: bending occurs perpendicular to the director; (2) bend: bending is parallel to the director; (3) twist: neighboring directors are rotated with respect to one another, rather than aligned. Inspired by liquid crystals [5, 7], after we obtain the local orthogonal frame for each point in tracts, we define at each point three scalar indices to describe the three types of local distortions, and a total distortion index:
where \(\frac{\partial \mathbf {u}_1}{\partial \mathbf {u}_i}\), \(i=1,2,3\), is the directional derivative of \(\mathbf {u}_1({{\varvec{x}}})\) along \(\mathbf {u}_i({{\varvec{x}}})\), i.e.,
where k is small, and \(\text {Diff}_d\) means the difference of two directors in the director representation [5]. See Sect. 2.1. Since the director field is not continuous in practice, we use an interpolation method to estimate \(\mathbf {u}_1({{\varvec{x}}}+k\mathbf {u}_i)\) and \(\mathbf {u}_1({{\varvec{x}}}-k\mathbf {u}_i)\) from local neighborhoods of \({{\varvec{x}}}+k\mathbf {u}_i\) and \({{\varvec{x}}}-k\mathbf {u}_i\). Let \({{\varvec{z}}}={{\varvec{x}}}\pm k\mathbf {u}_i\), then we estimate \(\mathbf {u}_1({{\varvec{z}}})\) from its neighborhood tangent directors of tracts:
where \(I({{\varvec{z}}})\) is a neighborhood of \({{\varvec{z}}}\), \(w_s({{\varvec{y}}}, {{\varvec{z}}}) = \frac{1}{Z\Vert {{\varvec{y}}} -{{\varvec{z}}} \Vert ^2}\), \(w_b({{\varvec{y}}},{{\varvec{x}}})\) is the bundle probability weight, and Z is the normalization factor such that \(\sum _{{{\varvec{y}}} \in I({{\varvec{z}}})} w_s({{\varvec{y}}}, {{\varvec{z}}})w_b({{\varvec{y}}}, {{\varvec{z}}})=1\). we use the inverse distance weight so that if \({{\varvec{z}}}\) is a point in a tract, then \(w_s=1\) and \(\mathbf {u}_1({{\varvec{z}}})\) is just the tangent vector of the tract. It can be proved that the interpolated \(\mathbf {u}_1({{\varvec{z}}})\) is the actually the principal eigenvector of \(\sum _{{{\varvec{y}}} \in I({{\varvec{z}}})} w_{s}({{\varvec{y}}},{{\varvec{z}}}) w_{b}({{\varvec{y}}},{{\varvec{z}}})\mathbf {u}_1({{\varvec{y}}}) \mathbf {u}_1^T({{\varvec{y}}})\). We set \(I({{\varvec{z}}})\) (i.e., \(I({{\varvec{x}}}\pm k\mathbf {u}_i)\)) as a ball with a radius of 2k, and k is normally set as \(1\,mm\), half of the typical voxel size of \(2\,mm\). If one would like to calculate indices using all points \({{\varvec{y}}}\) from all fiber bundles, then \(w_b({{\varvec{y}}},{{\varvec{x}}})=1\) as a constant. If one would like to calculate indices using points \({{\varvec{y}}}\) only from the same fiber bundle as \({{\varvec{x}}}\), then we may set \(w_b({{\varvec{y}}},{{\varvec{x}}})\) as a pre-determined probability of fiber clustering, or we could simply set \(w_b({{\varvec{y}}},{{\varvec{x}}})=1\) only if the angle between \(\mathbf {u}_1({{\varvec{y}}})\) and \(\mathbf {u}_1({{\varvec{x}}})\) is smaller than a threshold \(\theta _0\) (e.g., \(45^\circ \)), and \(w_b({{\varvec{y}}},{{\varvec{x}}})=0\), otherwise. If there is no specific mention, we set \(w_b({{\varvec{y}}},{{\varvec{x}}})=1\) only if the angle is smaller than \(45^\circ \).
3 Experiments
Synthetic Data Experiments. We generated three synthetic fiber tracts (splay, bend, and twist data) which demonstrate these three types of distortions, and then calculated the proposed scalar indices for each point in the fiber tracts. These scalar indices were used to color the fiber tracts. See Fig. 3. We omit OO, because \(\text {OO}=1-\text {OD}\). The first three rows in Fig. 3 show that the proposed splay, bend, and twist indices completely separate these 3 datasets, where only one among these 3 values is non-zero in each dataset. We also combined fiber tracts in splay and bend data, then calculated these indices using interpolation with points from all fiber bundles and with points only from the same fiber bundle (angle threshold \(\theta _0=45^\circ \)). The last two rows in the figure show that in a crossing area, it is better to use fiber tracts from the same bundle to calculate bundle specific distortion indices. Note that increased OD happens in crossing areas, compared with single bundle areas, and calculation of OD does not require an interpolation.
The corpus callosum of a real subject. Left: tracts are colored by OD values in points, where low/high values are in blue/red colors. Right: tracts are colored by splay, bend, and twist values in points, where splay, bend, and twist values are set as red, green and blue color channels, respectively.
Real Data Experiments. We performed tractography on a publicly available dataset with a single subject from DIPY (dipy.org). The data has a single shell with \(b\,=\,{2000}\,\mathrm{s/mm^{2}}\) and 150 directions on the shell. The corpus callosum (CC) was extracted based on diffusion ODF and deterministic tracking. Then, all six indices were calculated from the fiber tracts of CC. Figure 4 shows the scalar indices as colors of tracts. In the right subfigure, splay, bend, and twist values are set as red, green and blue color channels, respectively. The bending areas of CC are mainly in green, which means that the bend index is higher than slay and twist indices in those bending areas. We can also see red (high splay values) in fanning areas of CC, and blue (high twist values) when tracts are twisting in the 3D space.
4 Conclusion
We propose a unified mathematical framework called tract-based director field analysis (tDFA) with six scalar indices to quantify local geometric structure of fiber tracts. OD and OO are useful to quantify the degree of alignment and dispersion of fiber tracts; The distortion indices (i.e., splay, bend, twist and total distortion) demonstrate good sensitivity in both synthetic datasets and a real dataset on CC. To our knowledge, this is the first work to quantify orientational distortion (splay, bend, and twist) of fiber tracts. The proposed indices are rotationally invariant, because they are calculated from intrinsic quantities (i.e., local orthogonal frames) of tracts.
References
Basser, P.J.: New histological and physiological stains derived from diffusion-tensor MR images. Ann. N. Y. Acad. Sci. 820, 123–138 (1997)
Basser, P.J., Mattiello, J., LeBihan, D.: MR diffusion tensor spectroscropy and imaging. Biophys. J. 66, 259–267 (1994)
Basser, P.J., Pajevic, S., Pierpaoli, C., Duda, J., Aldroubi, A.: In vivo fiber tractography using DT-MRI data. Magn. Reson. Med. 44, 625–632 (2000)
Batchelor, P.G., Calamante, F., Tournier, J.-D., Atkinson, D., Hill, D.L.G., Connelly, A.: Quantification of the shape of fiber tracts. Magn. Reson. Med. 55(4), 894–903 (2006)
Cheng, J., Basser, P.J.: Director Field Analysis (DFA): exploring local white matter geometric structure in diffusion MRI. Med. Image Anal. 43, 112–128 (2018)
Cheng, J., Ghosh, A., Jiang, T., Deriche, R.: Model-free and analytical EAP reconstruction via spherical polar fourier diffusion MRI. In: Jiang, T., Navab, N., Pluim, J.P.W., Viergever, M.A. (eds.) MICCAI 2010. LNCS, vol. 6361, pp. 590–597. Springer, Heidelberg (2010). https://doi.org/10.1007/978-3-642-15705-9_72
Collings, P.J., Hird, M.: Introduction to Liquid Crystals: Chemistry and Physics. CRC Press, Boca Raton (1997)
Pajevic, S., Aldroubi, A., Basser, P.J.: A continuous tensor field approximation of discrete DT-MRI data for extracting microstructural and architectural features of tissue. J. Magn. Reson. 154, 85–100 (2002)
Savadjiev, P., et al.: Fusion of white and gray matter geometry: a framework for investigating brain development. Med. Image Anal. 18(8), 1349–1360 (2014)
Tax, C.M.W., et al.: Sheet Probability Index (SPI): characterizing the geometrical organization of the white matter with diffusion MRI. NeuroImage 142, 260–279 (2016)
Tuch, D.S.: Q-ball imaging. Magn. Reson. Med. Off. J. Int. Soc. Magn. Reson. Med. 52, 1358–1372 (2004)
Wedeen, V.J., et al.: The geometric structure of the brain fiber pathways. Science 335(6076), 1628–1634 (2012)
Zhang, H., Schneider, T., Wheeler-Kingshott, C.A., Alexander, D.C.: NODDI: practical in vivo neurite orientation dispersion and density imaging of the human brain. Neuroimage 61(4), 1000–1016 (2012)
Acknowledgement
This work was supported by Australian Research Council Projects FL-170100117, DP-180103424, LP-150100671, Intramural Research Program of NICHD (ZIA-HD000266), the National Key R&D Program of China (2016YFF0201002), and the Beijing Municipal Science and Technology Commission (BJMSJY-160153).
Author information
Authors and Affiliations
Corresponding author
Editor information
Editors and Affiliations
Rights and permissions
Copyright information
© 2018 Springer Nature Switzerland AG
About this paper
Cite this paper
Cheng, J. et al. (2018). On Quantifying Local Geometric Structures of Fiber Tracts. In: Frangi, A., Schnabel, J., Davatzikos, C., Alberola-López, C., Fichtinger, G. (eds) Medical Image Computing and Computer Assisted Intervention – MICCAI 2018. MICCAI 2018. Lecture Notes in Computer Science(), vol 11072. Springer, Cham. https://doi.org/10.1007/978-3-030-00931-1_45
Download citation
DOI: https://doi.org/10.1007/978-3-030-00931-1_45
Published:
Publisher Name: Springer, Cham
Print ISBN: 978-3-030-00930-4
Online ISBN: 978-3-030-00931-1
eBook Packages: Computer ScienceComputer Science (R0)