LUNet: deep learning for the segmentation of arterioles and venules in high resolution fundus images

Several attempts have been made to perform automatic A/V segmentation based on DFIs. Grisan and Ruggeri used proprietary data (2003) (Grisan and Ruggeri 2003) to extract vessels and classify A/V based on differences in the luminosity and contrast in the DFIs coupled with a vessel-tracking algorithm. The first use of convolutional neural networks for the A/V classification task was made by Welikala et al in 2017, who used 100 proprietary annotated DFIs from the UK Biobank dataset as the training dataset, and the DRIVE (Staal et al 2004) dataset for testing (Welikala et al 2017). In 2018, Hemelings et al reported on the first use of a fully convolutional network with an auto-encoder architecture, similar to a U-Net without skip connections, to perform the A/V segmentation task (Hemelings et al 2019). Two models were trained independently using the DRIVE (Staal et al 2004) and the high-resolution fundus (HRF) (Budai et al 2013) datasets. In 2021, Hu et al developed VC-Net, another U-Net variant incorporating a vessel constraint module to improve their model segmentation (Hu et al 2021). They trained their model on the DRIVE (Staal et al 2004), HRF (Budai et al 2013), and LES-AV (Orlando el al, 2018) datasets and used Tongren and Kailuan, two private datasets, for external validation (Hu et al 2021). In 2022, Galdran et al (2022) used a concatenation of two U-Net named Little W-Net to train two models for the A/V segmentation task; the first trained on the DRIVE dataset (Staal et al 2004) and the second oon the HRF dataset (Budai et al 2013). They used the LES-AV dataset (Orlando el al, 2018) as external validation of their model trained on the DRIVE dataset. Finally, Zhou et al developed BFN, which uses three U-Net models trained within an adversarial framework; one artery segmenter, one vein segmenter, and one multi-class segmenter (Zhou et al 2021). They trained three independent BFN on the DRIVE (Staal et al 2004), HRF (Budai et al 2013), and LES-AV (Orlando el al, 2018) before releasing a final BFN (Zhou et al 2022) trained on these three datasets. They evaluated their final model on the IOSTAR (Abbasi-Sureshjani et al 2015, 2016) dataset.

While previous work has focused mainly on DRIVE, HRF, and LES-AV, three public datasets that contain 40, 45, and 22 DFIs with manual reference segmentations, respectively, these datasets have different resolution, framing, field of view (FOV) and population sample. More specifically, the DRIVE dataset contains DFIs centered on the macula, with a FOV of 45° and a resolution of 584 × 565 pixels, acquired during a screening for diabetic retinopathy in the Netherlands. The HRF dataset contains DFIs centered on the macula, with a resolution of 2336 × 3504 pixels, acquired in Germany from patients with glaucoma or diabetic retinopathy and from healthy individuals. LES-AV contains DFIs centered on the optical disc, with a FOV of 30° and a resolution of 1444 × 1620. Patient age distribution is available for the LES-AV dataset only. The heterogeneity of these datasets led most researchers to use them to train an independent model for each dataset; poor generalization performance was observed for external datasets. Thus, there is a need for a robust A/V segmentation model with high performance that can generalize across external test sets with varying distribution shifts.

In this research, we focused on optic disc-centered DFIs, with a FOV of 30° and a high resolution of 1444 × 1444 pixels. In addition, a new DFI dataset, named UZLF, consisting of 240 crowd-sourced A/V segmentations performed by 15 medical students and subsequently corrected by a senior annotator (JVE) was created. LUNet, a novel robust deep learning (DL) algorithm tailored to the A/V segmentation task, is introduced. The generalization performance of LUNet was evaluated using 30 newly manually segmented A/V DFIs from UNAF and INSPIRE-AVR (Niemeijer et al 2011, Benítez et al 2021) as well as on the publicly available LES-AV dataset (Orlando el al, 2018) and HRF dataset (Budai et al 2013) which have reference A/V segmentations. LUNet was benchmarked against Little W-Net (Galdran et al 2022), VC-Net (Hu et al 2021) and BFN (Zhou et al 2021), three open-source state-of-the-art (SOTA) algorithms.

A total of five datasets were used in our experiments. The University Hospital UZ Leuven Fundus (UZLF) dataset was used for developing the models while the four other datasets were used as external test sets to evaluate the generalization performance of LUNet and benchmark algorithms. The datasets are summarized in table 1.

Table 1. Summary of the optic disc-centered DFI datasets with A/V segmentation, including the average percentage of unknown blood vessels per image. Median (Q1–Q3) are provided for Age.

Name	${{\rm{N}}}^{\underline{{\rm{o}}}}$ DFIs	Country	Age	Original FOV	Purpose	Unknown
UZLF	240	Belgium	62 (50–73)	30°	Train & test	0.2%
UNAF	15	Paraguay	—	45°	External test	0%
INSPIRE-AVR	15	United State	—	30°	External tes	0%
LES-AV	20	Belgium	71 (62–80)	30°	External test	3%
HRF (centered version)	28	Germany	—	60°	External test	13.5%

2.1.1. University Hospital UZ Leuven Fundus UZLF dataset

Human data were obtained within the context of the study "Automatic glaucoma detection, a retrospective database analysis" (study number S60649). The Ethics Committee Research UZ/KU Leuven approved this study in November 2017 and waived the need for informed consent. A total of 115 237 optic disc-centered DFIs from 13 185 unique patients, captured between 2010 and 2019 were provided by the UZL in Belgium. These DFIs were taken with a Visucam Pro NM camera with 30° FOV (Zeiss). The resolution of these DFIs was 1444 × 1444 pixels which is higher than most public DFIs datasets and enables the visualization of smaller blood vessels. The median age of the DFIs and the interquartiles (Q1–Q3) were 64 (52–75) years and 52% were female. The exclusion criteria included patients under 18 years of age and low-quality DFIs (FundusQ-Net <6) (Abramovich et al 2023) (figure 2). Active learning was used to proactively select DFIs for manual segmentation. Specifically, among the extracted DFIs, an exploration-exploitation strategy was applied to select 240 DFIs for A/V segmentation. The exploration step consisted of randomly stratified sampling of some DFIs to annotate according to the patient's sex and the imaged eye (right/left). The exploitation steps, or active learning steps, consisted of selecting a set of DFIs with low LUNet A/V segmentation performance. Low segmentation performance was defined as a lack of vessel continuity. A subset of 240 DFIs denoted UZLF, from 232 unique patients, were selected and manually segmented. The patients included in UZLF were between 18 and 90 years of age (median (interquartile): 62 (50–73) years) and 58% were female. UZLF was 57% composed of left eye DFIs. Patients who belong to the UZLF dataset were separated into the following classes: (1) normal ophthalmic findings, (2) normal tension glaucoma (NTG), (3) primary open angle glaucoma (POAG), and (4) other condition.

Zoom In Zoom Out Reset image size

2.1.2. Universidad Nacional de Asunción Fundus (UNAF) external test set

2.1.3. University of Iowa Hospitals and Clinics (INSPIRE-AVR) external test set

2.1.4. LES-AV dataset

2.1.5. HRF dataset

The INSPIRE-AVR, UZLF and UNAF datasets were manually segmented by the retinal experts of the UZ Leuven Hospital using the Lirot.ai app developed by Fhima et al (2022a) and following the protocol described in Fhima et al (2022b).

A total of sixteen annotators experienced in micro-vascular research worked between July 2021 and January 2023 to build the UZLF dataset. Annotators were divided into two groups; (1) one experienced senior annotator, an ophthalmology resident with a PhD in retinal vascular biomarkers, and (2) 15 junior annotators who were all graduate students in medicine with a research internship of >1 month in the Research Group of Ophthalmology and trained by the senior annotator. The UNAF and INSPIRE-AVR datasets were annotated by the senior annotator. For the UZLF dataset, the junior annotators performed the first segmentation of the 240 DFIs, of which 174 were later corrected by the senior annotator to avoid mistakes and to correct for different annotation styles.

A DFI was modeled by a couple (X, y), where:

• X ∈ R^w×h×3 represents a DFI with a width of w, a height of h pixels and 3 channels, i.e. red, green and blue;
• y ∈ R^w×h×3 where y=[y_a , y_v , y_ukn ] with y_a , y_v , y_ukn ∈R^w×h×1, i.e., y_a , y_v and y_ukn are binary images, where y_a is the segmented arterioles, y_v the segmented venules and y_ukn represent the blood vessel that were not distinguishable, with white pixel if part of a blood vessel, and black pixel otherwise.

Where h = w = 1444.

2.3.1. Train-validation-test set elaboration

2.3.2. Preprocessing

The attention U-Net architecture was used as the backbone for our DL algorithm (Oktay et al 2018). The A/V segmentation task is challenging due to the long-range dependencies of the blood vessels and their small width. To achieve accurate A/V segmentation, it is necessary to have a model that can identify both local and global dependencies. The local dependencies help detect the smaller blood vessels more accurately, while the global dependencies enable the reconstruction of the full blood vessel. In convolutional neural networks, the local dependencies are captured by the convolution operation, while the larger range dependencies are dependent on the receptive field of the network. The receptive field can be increased by augmenting the depth of the model and increasing the kernel size of the convolutional layers, using dilated convolution, or using some max pooling layers. The first two solutions would lead to an augmentation of the model parameters with a small increase of the receptive field, while the max pooling would not add any parameter and would lead to a much larger receptive field. Nevertheless, due to the small size of the blood vessel in a classical DFI, the max pooling could be an excessively aggressive strategy, which would result in the loss of small-blood-vessel visibility after several iterations. To tackle this problem, LUNet was designed to include several improvements, including; (1) a new double dilated convolution block with an increased receptive field, (2) a long tail that works on the full pixel resolution, (3) an increased depth, (4) an over-representation of the feature extracted by the encoder compared to the one reconstructed by the decoder at each level of depth of the LUNet auto-encoder. The LUNet architecture is shown in figure 3.

Zoom In Zoom Out Reset image size

2.4.1. Double dilated convolution block

Dilated convolution enables a larger receptive field but would lead to a less accurate feature representation of small details (due to the dilation rate). To tackle this limitation, both classical and dilated convolutions with a kernel size of 7 were used in the model. The double dilated convolution block also incorporated a spatial dropout 2D regularization (Tompson et al 2015) and a batch normalization layer with the ReLU activation function. A diagram of the double dilated convolution block is presented in figure 4.

Zoom In Zoom Out Reset image size

2.4.2. Long tail

Other approaches have described the A/V segmentation problem as a multiclass segmentation problem leading to a multiclass cross-entropy loss minimization (Hu et al 2021, Hemelings et al 2019). We formulated the problem differently, as a binary multi-label segmentation because it is common to see superimposition of arterioles and venules on several pixels of a DFI. Furthermore, some manual segmentations may contain pixels labeled as unknown. These pixels were not penalized for being classified as A or V. However, LUNet was still trained to detect these pixels and eventually classify them. A regularization term was added to the loss function to learn to detect blood vessels without A/V distinction. Accordingly, the LUNet loss function, L_LUNet, was defined as the sum of three losses computed independently on the arterioles, the venules, and the overall blood vessels

$$\begin{eqnarray}&&\begin{array}{l}{L}_{\mathrm{LUNet}}(y,\hat{y})=L({y}_{a}\times (1-{y}_{ukn}),{\hat{y}}_{a}\times (1-{y}_{ukn}))\\ +\,L({y}_{v}\times (1-{y}_{{ukn}}),{\hat{y}}_{v}\times (1-{y}_{{ukn}}))\\ +\,L({\max }({\hat{y}}_{a},{\hat{y}}_{v}),{\max }({y}_{a},{y}_{v},{y}_{{ukn}}))\end{array}\end{eqnarray} \tag{ 1 }$$

with

$$\begin{eqnarray}&&\begin{array}{l}L(y,\hat{y})={\lambda }_{1}\times {L}_{\mathrm{BCE}}(y,\hat{y})+{\lambda }_{2}\times {L}_{\mathrm{dice}}(y,\hat{y})\,+\\ {\lambda }_{3}\times {L}_{{cl}\mathrm{Dice}}(y,\hat{y})+{\rm{\nabla }}(\hat{y})\end{array}\end{eqnarray} \tag{ 2 }$$

where y is the ground truth segmentation, $\hat{y}$ is the predicted probability map, L_BCE is the binary cross-entropy loss, L_dice is the dice loss, L_clDice is the centerline dice loss (Shit et al 2021) and where the gradient of $\hat{y}$ is minimized to favor the continuity of the blood vessels by constraining the change of value in the probability map to be smoother.

LUNet hyperparameters were manually fine-tuned using the validation set. LUNet was trained for up to 1300 epochs, while retaining the model which obtains the minimum validation loss. A batch size of 8 and an Adam optimizer with 1e⁻⁴ for the learning rate were used. The L_LUNet loss function was used with λ₁ = λ₂ = 1 and λ₃ = 0.3.

For each training DFI loaded, random online data augmentation was performed with (1) horizontal and vertical flip, (2) transposition, (3) rescaling of the input and output to a lower resolution uniformly sampled between 800 × 800 and 1472 × 1472, and (4) color jittering. Test time data augmentation was used for the test set prediction (Wang et al 2018) with the following rotation angles [0, 30, 60, 90, 120, 150, 180, 210, 240, 270, 300, 330], as well as transposition of the DFI which resulted in 24 predicted segmentations for a single original DFI. The inverse transform for each prediction was performed. The final predicted probability map for the segmentation was the pixel-wise average of all the 24 segmentations.

Benchmark against SOTA algorithms: LUNet performance on the UZLF dataset was compared to that of Little W-Net (Galdran et al 2022), BFN (Zhou et al 2021) and VC-Net (Hu et al 2021), three SOTA algorithms, that were trained on the UZLF-train dataset. On the external datasets, LUNet was compared to Little W-Net (Galdran et al 2022), VC-Net (Hu et al 2021) and BFN (Zhou et al 2021) trained on the UZLF dataset as well as the publicly trained version of Little W-Net* (Galdran et al 2022) and BFN* (Zhou et al 2022).

Benchmark against junior performance: The DFIs of the UZLF-test dataset were annotated by a junior annotator and then reviewed by the senior annotator. This enables a comparison of the dice score between an individual junior annotator and the senior annotator. This dice score statistic allows to benchmark LUNet to the performance of a human junior annotator.

A dice score was computed for the arterioles segmentation (dice_a ) and for the venules segmentation (dice_v ). In order to estimate the 95% confidence interval (CI), a bootstrap method was employed by repeatedly sampling 80% of the test set with replacement, and computing the mean score of each sample. This procedure was carried out 1000 times, and the resulting distribution of means was used to determine the lower and upper bounds of the 95% confidence interval. In order to monitor the performance of LUNet with respect to the number of annotated DFIs, two learning curves one for the dice_a and one for the dice_v were computed. The reported performance was computed on the UZLF-test set as the training set size increased. For the LES-AV and HRF datasets, some of the ground truth blood vessels were annotated as unknown. The unknown pixels were not taken into account in the computation of the two dice scores. Furthermore, due to the imbalanced nature of blood vessels in DFIs, we also report the Matthews correlation coefficient (Chicco and Jurman 2020) for LUNet and the benchmark models for both the local and external test sets.

Finally, LUNet was evaluated on the local test set with respect to the ability to correctly estimate vasculature biomarkers (VBMs). VBMs were computed using the open source PVBM toolbox (https://pvbm.readthedocs.io/) (Fhima et al 2022b). The VBMs include: the area (AREA), the tortuosity index (TI), the median arc-chord tortuosity (TOR), the length (LEN), the median branching angles (BA), the number of blood vessels that intersect with the optic disc (START), the number of endpoints (END), the number of intersection points (INTER), the fractals dimension (D0, D1 and D2). Each biomarker was computed independently on the arterioles and the venules and the Pearson correlation between the ground truth and the estimated VBMs were computed. The average Pearson correlation between all the VBMs was also computed in order to provide an overall performance measure.

LUNet achieved dice scores of 81.99/84.54 for A/V segmentations on the local test set. These were higher than those of the SOTA benchmark algorithms Little W-Net (Galdran et al 2022), VC-Net (Hu et al 2021) and BFN (Zhou et al 2021) which achieved dice scores of 79.76/82.62, 78.31/81.89 and 78.14/80.79, respectively. The performance of LUNet were better for A (81.99 versus 81.50) and for V (84.54 vs 83.74) than an average junior annotator. Table 2 and figure 5 summarizes the quantitative results on the UZLF-test. Figure 6 shows examples of segmentations performed by LUNet, a junior annotator, an expert annotator, alongside those from the top-performing generalizable benchmark algorithm, VC-Net (Hu et al 2021). Specifically, the DFI with the best, intermediate and worst LUNet performance was selected. Furthermore, LUNet outperforms Little W-Net, VC-Net and BFN in estimating most VBMs (table 3). LUNet got the highest Pearson correlation score in 8 out of 22 VBMs and the second highest in 12 out of 22 VBMs. Meanwhile, the junior annotators had the best performance in 9 VBMs and were the second best in 2 VBMs. The learning curves of LUNet are shown in figure 7.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 8 presents the segmentations produced by LUNet alongside those from the top-performing generalizable benchmark algorithm, VC-Net (Hu et al 2021), as well as the ground truth for comparison for a DFI example for each of the external test sets. Table 4 summarizes the quantitative results on the external test sets.

Zoom In Zoom Out Reset image size

LUNet achieved average dice scores of 82.30/84.75 for A/V segmentation on the LES-AV dataset. This performance was superior to those of the SOTA benchmark algorithms Little W-Net (Galdran et al 2022), VC-Net (Hu et al 2021) and BFN (Zhou et al 2021) which achieved dice scores of 80.19/82.97, 78.75/82.49 and 78.68/80.34, respectively. Its performance was also higher than that of the public version of Little W-Net* (Galdran et al 2022) which achieved A/V dice scores of 65.32/70.27. The public version of BFN* (Zhou et al 2022) was not evaluated on LEV-AV because it was originally trained on it. LUNet achieved dice scores of 73.31/79.04 for A/V segmentations on the UNAF dataset. This performance was superior to those of the SOTA benchmark algorithms Little W-Net (Galdran et al 2022), VC-Net (Hu et al 2021) and BFN (Zhou et al 2021) which achieved dice scores of 66.71/68.61, 68.71/74.43 and 68.01/72.67, respectively. Its performance was also higher than those of the public version of Little W-Net* (Galdran et al 2022) and BFN* (Zhou et al 2022) which achieved A/V dice scores of 52.58/51.32 and 64.13/72.78, respectively. LUNet achieved dice scores of 73.58/77.53 for A/V segmentations on the INSPIRE-AVR dataset. These were higher than the SOTA benchmark algorithms Little W-Net (Galdran et al 2022), VC-Net (Hu et al 2021) and BFN (Zhou et al 2021) which achieved dice scores of 68.61/70.56, 71.28/75.34 and 68.61/71.92, respectively. They were also higher than those of the public version of Little W-Net* (Galdran et al 2022) and BFN* (Zhou et al 2022) which achieved A/V dice scores of 58.08/66.37 and 62.98/68.10, respectively. LUNet achieved dice scores of 78.12/80.39 for A/V segmentations on the cropped HRF dataset. This performance was superior to those of the SOTA benchmark algorithms Little W-Net (Galdran et al 2022) and, VC-Net (Hu et al 2021) and BFN (Zhou et al 2021) which achieved dice scores of 59.01/47.44, 72.90/76.99 and 75.78/79.99, respectively. Its performance was also higher than those of the public version of Little W-Net* (Galdran et al 2022) and BFN* (Zhou et al 2022) which achieved A/V dice scores of 57.93/61.42 and 70.95/75.47, respectively. Performance in terms of Matthews correlation coefficient are reported in table A1 and are consistent with our observations in terms of DICE score.

An ablation study was conducted on the validation set by removing the following architecture components, LT: long tail, CL: custom loss, DDCB: double dilated convolution block (figure 9). This ablation study demonstrates the importance of the individual components.

Zoom In Zoom Out Reset image size