Keywords

1 Introduction

The recently developed 3D optical coherence tomography angiography (OCT-A) is an efficient and non-invasive imaging modality that is able to provide rich micrometer-level axial resolution with rich depth-resolved 3D microvascular information. It has been widely applied to analyze microvascular diseases like diabetic retinopathy (DR) and age-related macular degeneration (AMD) [5]. The routine way for analyzing retinal microvasculature in OCT-A images is mostly based on the 2D en face projections, while direct 3D analysis including rich depth-resolved geometric and topological information is hardly considered due to the challenges of low image quality and high capillary complexity. A comprehensive overview including clinical applications using 2D en face images can be found in [5]. Most of the analysis [5, 6] are performed on projected OCT-A images from different retinal layers to study the correlations between vessel biomarkers and disease progression. In comparison with other imaging modalities, OCT-A is proven to be a more useful tool for describing capillary loss between normal and non-proliferative diabetic retinopathy (NPDR) subjects based on the analysis of en face projections from different depth layers [3, 5].

Despite the challenges of 3D OCT-A analysis, it is valuable to investigate and analyze 3D microvascular changes by considering the rich depth-resolved information instead of just performing 2D analysis. The demands for high-quality 3D OCT-A vessel visualization and processing motivate us to build a surface-based microvascular analysis framework. A unified Reeb analysis approach [10] was proposed to detect tissue outliers on smooth and triangulated cortical surface representation [9] of brain magnetic resonance images (MRI). In this work, we propose to design a similar routine by taking advantage of the analysis on Reeb graph which is built on a delicate feature function defined on 3D OCT-A vessel surface. Reeb analysis provides the possibility to analyze local geometric and topological changes based on the description of full vascular topology. A variety of functions [2, 10] have been employed to set up Reeb graphs for removing surface outliers and identifying surface protrusions. In this work, we propose to employ the geodesic distance transform as a distinguishable feature function for precise Reeb graph construction on mesh surface.

Based on Reeb graph, we are able to extract potentially important biomarkers to perform 3D microvascular analysis on clinical data. To our knowledge, this is for the first time that a complete microvascular modeling framework is proposed to analyze geometric and topological biomarkers in 3D OCT-A images. The surface representation provides a high-quality 3D microvascular visualization to assist clinical practice. Meanwhile, the proposed Reeb analysis framework gives the possibility of analyzing 3D microvascular data in a new perspective.

2 Methodology

2.1 3D Surface Representation of OCT-A Microvasculature

To establish a high-quality 3D microvascular surface representation, we first design an effective enhancement and segmentation method to extract OCT-A vessel networks. Speckle noise reduction is achieved via 3D curvelet denoising on OCT-A volumes. Afterwards, we exploit the optimally oriented flux (OOF) filter [7] by considering its nice property of processing closely located curvilinear structures to enhance the dense OCT-A microvasculature as shown in Fig. 1(b). The OOF is a quadratic function and defines a symmetric matrix \(\mathcal {Q}(\mathbf x ,r)\) which can be decomposed as \(\mathcal {Q}(\mathbf x ,r) = \sum _{i=1}^{3} \lambda _i(\mathbf x ,r)\mathbf {v}_i(\mathbf x ,r)\mathbf {v}_i(\mathbf x ,r)^T,\) where \(\lambda _i(\mathbf x ,r)\) and \(\mathbf {v}_i(\mathbf x ,r)\) denote eigenvalues and eigenvectors, respectively. One has \(\lambda _1(\cdot ) \le \lambda _2(\cdot ) \ll \lambda _3(\cdot ) \approx 0\) inside vessels due to that the OCT-A vessel voxels have intensity values higher than background voxels. The final vesselness map for each voxel is defined by taking the largest magnitude: \(\mathcal {P}(\mathbf x ):={\text {max}}\{\underset{r}{{\text {max}}}\{-\frac{1}{r^2}\lambda _1(\mathbf x ,r)\},0\}.\) The binarized vessel segmentations are obtained via a proper thresholding. We also extract the volume of interest (VOI) including superficial and deep layers from each corresponding OCT volume using OCTExplorer [1]. A continuous OCT-A vessel surface representation is constructed by defining triangular meshes on 3D segmented vessel boundary masks. See Fig. 1(c). Hence, a complete mesh surface \(\mathcal {M}\) is represented by vertices and triangular faces, i.e. \(\mathcal {M}=\{\mathcal {V},\mathcal {T}\}\), where \(\mathcal {V}=\{\mathcal {V}_i | i=1,\cdots ,N_\mathcal {V}\}\) and \(\mathcal {T}=\{\mathcal {T}_j | j=1,\cdots ,N_\mathcal {T}\}\) denote the set of vertices and triangles, respectively.

Fig. 1.
figure 1

3D OCT-A microvascular enhancement surface reconstruction.

Full size image

2.2 Intrinsic Reeb Graph Construction Using Geodesic Distance Transform on Vessel Mesh Surfaces

Reeb Graph. Based on the surface representation, we propose to analyze 3D OCT-A vessel geometry and topology in a new perspective via Reeb graph [8]. A Reeb graph R(f) is intuitively a graph of level contours constructed by a Morse function f [4] defined on manifold. The Reeb graph on vessel surface can reflect the intrinsic topological changes of retinal microvasculature since its level contours vary only at critical points of the Morse function. For a triangular mesh \(\mathcal {M}=\{\mathcal {V},\mathcal {T}\}\), the function f is defined on each vertex in \(\mathcal {V}\). The level sets of f are sampled at a set of K values \(\xi _0< \xi _1< \cdots < \xi _k\) with the set of contours as \(\mathcal {C}=\{\mathcal {C}_k^n,0\le k \le K-1, 1 \le n \le N_k\}\), where \(N_k\) is the number of contours at level k, and where \(\mathcal {C}_k^n\) denotes the n-th contour at this level. In this work, we compute the Reeb graphs for 3D OCT-A microvascular surfaces by following the similar routine as in [10]. The edge \(E_{k,k+1}\) between neighbouring level contours are established by connecting a contour \(\mathcal {C}_k^{n_1}\) at level \(\xi _k\) and a contour \(\mathcal {C}_{k+1}^{n_2}\) at level \(\xi _{k+1}\) if they belong to the same connected component in the vascular region \(R_{k,k+1}=\{\mathcal {V}_i \in \mathcal {M}|\xi _k \le f(\mathcal {V}_i)\le \xi _{k+1}\}\). Finally, a complete Reeb graph on mesh surface \(\mathcal {M}\) with undirected edge connection of level contours as nodes is constructed to represent the 3D microvascular networks. Figure 2 shows the neighboring contours generated from a feature function on each vertex. The centroid of each contour is visualized to intuitively represent the graph nodes.

Fig. 2.
figure 2

Level contours generated from a feature function on each mesh vertex.

Full size image

Distinguishable Surface Transition via Geodesic Distance Transform. The various choices of the feature function f determines different description aspects and abilities for the same surface \(\mathcal {M}\). A proper surface-based feature representation can be used to emphasize different local characteristics while preserving the global topology. Retina-related diseases often cause pathological variations that can gradually affect the microvascular geometry and topology during disease progression. It is of vital importance to investigate and quantify 3D microvascular changes at critical points such as bifurcations and small capillaries. Hence, we propose to integrate a novel feature, i.e. the mesh geodesic distances as the Morse function, into our Reeb analysis framework to represent the topological transitions on 3D OCT-A microvascular surface. See Fig. 3(a).

The length of the minimal geodesic from vertices \(\mathcal {V}_i\) to \(\mathcal {V}_j\) is called their geodesic distance on surface \(\mathcal {M}\) and is represented by \(d(\mathcal {V}_i, \mathcal {V}_j)\). Thus, the geodesic distance transform of a point \(\mathcal {V}_p\) on the manifold \(\mathcal {M}\) can be defined as \(d_{\mathcal {V}_p} (\mathcal {V}_q)=d(\mathcal {V}_p,\mathcal {V}_q), \forall \mathcal {V}_q \in \mathcal {M}\). By calculating the geodesic distances between a given seed point \(\mathcal {V}_p\) and the other vertices \(\mathcal {V}_q\), we can conveniently obtain an intuitive representation of the microvascular topology on top of the complex linkages and hierarchical tree structures in OCT-A images. The seed point \(\mathcal {V}_p\) in our framework is automatically defined by finding the closest vertex to the volume center. The geodesic distances can reflect vascular expansion and are distinguishable at different locations. In Fig. 3(b), we can see the level contours are nicely distributed in vessel perpendicular direction with smooth transitions based on the proposed feature function f. Then by following the routine defined in [10] we build a Reeb graph to precisely represent the partition of OCT-A vessel surfaces.

Vessel Classification Based on Level-Set Contours. The OCT-A vessel topology can be easily affected by tailing artifact [12] which appears as connected “fake” blood flow signals and causes large vessel thickening below. This can result in unreliable microvascular measurements for further diseases study. Thus, we propose to analyze large vessels and small capillaries independently via a level-contour based classification step. A level contour is represented as a polyline composed of points intersecting edges of \(\mathcal {M}\) at the level \(\xi _k\). The proposed mesh geodesic distance function generates level contours that can fully represent the vessel cross-sectional boundaries as shown in Fig. 3(b). By considering the contour curve length \(L(\mathcal {C}_k^{n})\) obtained in the Reeb graph construction procedure, we can use it as a feature to categorize each level contour into large or small vessel groups. The original OCT-A vessel signal strengths can have large variations across healthy and pathological subjects, which cause differences in vessel geometry. Thus, we need define an automatic way to select suitable threshold values for different data. We first give a rough estimation of the large vessel contour size by taking the mean \(\mu _L\) of the maximum \(10\%\) level contours. Then we empirically define the automatic thresholding levels \(T_L\) at \(0.5\mu _L\) for different subjects. In Fig. 3(c), we show a typical example of the segmented level contours on 3D microvascular surface. The level contours of large vessels and small capillaries are well separated into different classes.

Fig. 3.
figure 3

Reeb graph construction and vessel classification on vessel surface.

Full size image

2.3 Geometric and Topological Feature Extraction via Reeb Analysis

The Reeb graph provides the possibility and convenience to investigate 3D geometric and topological variations in OCT-A images. It is also interesting to see the great potential of using the 3D surface representation as a novel perspective for retinal microvascular analysis. To this end, we employ Reeb analysis to extract the following geometric and topological biomarkers on Reeb graph.

Number of Bifurcations: The bifurcation number indicates the complexity of the hierarchical retinal microvascular tree. By taking each level contour \(\mathcal {C}_k^n\) as a node, we count the total number \(B_k^n\) of its 1-ring connected nodes. If \(B_k^n \ge 3\), we label it as 1 to indicate there exist a bifurcation point nearby. We label all its 1-ring connected nodes as 0 to avoid repeated counting in later iterations. The whole number of bifurcation points is obtained by summing the labels of all nodes, i.e. \(B_{total}=\sum _{k=0}^{K-1}\sum _{n=1}^{N_k}B_k^n\). In Fig. 4(a), we highlight the closest level contour where a nearby bifurcation exist.

Number of Ending Points: The ending points can approximately represent the amount of 3D small capillary segments in the deepest microvascular level. We iteratively count each node which has only one connection in its 1-ring neighborhood and obtain the ending nodes summation as a biomarker. In Fig. 4(b), we highlight all the ending level contours.

Number and Curve Length of Level Contours: The contour curve length can represent the microvascular size at each level. Since the level contours are equally sampled for different volume images, the curve length summation of all level contours can be considered as a measure to roughly reflect the amount of 3D blood flow signals of the whole volume at each moment.

Number of Microvascular Loops: The complex linkages in retinal vessel circulation can easily form multiple closed capillary loop areas, as shown in the histology study [11] in deep retinal layers. One capillary can split out from the main vessel and again merge into the same structure. This motivates us to quantify the small loop quantity in 3D OCT-A vessel topology and find its relation with disease progression. A vascular “loop” is essentially a “cycle” in the graph since we may have multiple level contours as nodes in a single vessel segment. The cycle detection is achieved via an iterative dynamic searching process. A dynamic path is iteratively updated by addling nodes along the deep-first searching direction. One cycle path is detected when the starting node appears for the 2nd time during the searching process. Then it restarts again from the beginning by tracing another neighboring path. Due to the complex vessel connections, we restrict the loop counting process by considering only the chordless cycle without having any other edge inside. In Fig. 4(c), all the contours that belong to a loop are visualized and we can see that many loops appear below the main vessels.

Fig. 4.
figure 4

Visualization of the detected feature via Reeb analysis.

Full size image

3 Experimental Results

A clinical OCT-A dataset of diabetic retinopathy (OCTA-DR) has been established for analyzing the obtained Reeb-based biomarkers. There are in total 100 OCT-A volume images from 100 subjects including 40 normal controls (NC), 20 severe non-proliferative diabetic retinopathy (NPDR) and 40 proliferative diabetic retinopathy (PDR). Only one eye of each patient has been selected for analysis to ensure data independence.

The proposed Reeb analysis framework is applied to extract geometric and topological biomarkers for each subject. We perform statistical analysis in a pair-wise way on all three groups based on Wilcoxon rank sum test with a significance level at 0.05. We calculate in total 11 biomarkers including the number of bifurcations, number of level contours, length of level contours for large, small and all vessels separately, as well as the number of ending points and loops.

The statistical analysis shows high significance (with \(p<0.05\)) on 10 biomarkers except the bifurcation numbers on small capillaries in distinguishing normal and PDR groups. This could be caused by the fact that most of the detected bifurcations in OCT-A images are split from large vessels as shown in Fig. 4(a). We provide the box and whisker plot in Fig. 5(a) to present the significant difference between each group for the number of bifurcations in large vessels, which achieves a high significance with \(p<0.001\) between normal and PDR groups. Similarly, the number of ending points and the number of level contours in small capillaries also show strong differences (with \(p<0.001\)) between this two groups as shown in Fig. 5(b)–(c). They also achieve significance (with \(p<0.05\)) when their differences between normal and severe NPDR groups are compared. These indicate that there is a small capillary loss in severe NPDR groups compared with normal controls. However, the differences between severe NPDR and PDR groups are not so significant for both biomarkers. Some of the other biomarkers also show insignificant differences when the pair-wise analysis is performed against severe NPDR. These could be caused by the limited data usage and thus further analysis should be performed by collecting more data in the study.

Fig. 5.
figure 5

Statistical analysis of Reeb graph based geometric and topological biomarkers for DR on large, small and all vessels separately.

Full size image

In Fig. 5(d)–(e), we can observe that the length of overall level contours on large vessels and small capillaries present very significant differences between normal and severe NPDR groups, as well as normal and PDR groups with a p value \({<}\,0.001\). For the difference between severe NPDR and PDR groups, the length of level contours achieves a p value \({<}\,0.001\) on large vessels and a p value \({<}\,0.05\) on small capillaries. These findings indicate that the length of level contours can be considered as a potentially important feature for the discrimination of different disease stages. Although the proposed 3D microvascular segmentation method has extracted most of the small capillaries, but it is still restricted by the low OCT-A image quality and capillary visibility. More accurate statistical performance can be obtained by exploring precise vessel segmentations with better image quality. In Fig. 5(f), we present the statistical analysis based on the number of topological loops to describe the difference between DR stages. We can see that it shows strong significance in distinguishing normal with the other two groups severe NPDR and PDR with \(p<0.05\) and \(p<0.001\), respectively. There is an increasing tendency in the amount of loops along with disease progression. This may indicate that the topological changes in the later stages of DR cause more closed capillary areas. It shows the great value of exploring more potentially important topological features in further Reeb analysis.

4 Conclusion

In this work, we have proposed an innovative 3D OCT-A microvascular analysis method by designing a Reeb analysis framework on high-quality OCT-A vessel surface representations. Several important geometric and topological vessel biomarkers have been extracted from the Reeb graph built on geodesic distances. Significant statistical performance of microvascular biomarkers has been found on DR discrimination. The proposed Reeb analysis framework opens great promise for analyzing 3D retinal microvasculature in a new perspective.