{'fraction_non_alphanumeric': 0.035432194904883745, 'fraction_numerical': 0.01429524974969407, 'mean_word_length': 4.764347547290798, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Detection of diabetic retinopathy using longitudinal self-supervised learning

Longitudinal imaging is able to capture both static anatomical structures and dynamic changes in disease progression towards earlier and better patient-specific pathology management. However, conventional approaches for detecting diabetic retinopathy (DR) rarely take advantage of longitudinal information to improve DR analysis. In this work, we investigate the benefit of exploiting self-supervised learning with a longitudinal nature for DR diagnosis purposes. We compare different longitudinal self-supervised learning (LSSL) methods to model the disease progression from longitudinal retinal color fundus photographs (CFP) to detect early DR severity changes using a pair of consecutive exams. The experiments were conducted on a longitudinal DR screening dataset with or without those trained encoders (LSSL) acting as a longitudinal pretext task. Results achieve an AUC of 0.875 for the baseline (model trained from scratch) and an AUC of 0.96 (95% CI: 0.9593-0.9655 DeLong test) with a p-value<2.2e-16 on early fusion using a simple ResNet alike architecture with frozen LSSL weights, suggesting that the LSSL latent space enables to encode the dynamic of DR progression.


Introduction
According to the International Diabetes Federation, the number of people affected by diabetes is expected to reach 700 million by 2045 [14]. Diabetic retinopathy (DR) affects over one-third of this population and is the leading cause of vision loss worldwide [8]. This happens when the retinal blood vessels are damaged by high blood sugar levels, causing swelling and leakage. In fundus retina images, lesions appear as leaking blood and fluids. Red and bright lesions are the type of lesions that can be commonly identified during DR screening. The blindness incidence can be reduced if the DR is detected at an early stage. In clinical routine, color fundus photographs (CFP) are employed to identify the morphological changes of the retina by examining the presence of retinal lesions such as microaneurysms, hemorrhages, and soft or hard exudates. The international clinical DR severity scale includes no apparent DR, mild non-proliferative diabetic retinopathy (NPDR), moderate NPDR, severe NPDR, and proliferative diabetic retinopathy (PDR), labeled as grades 0, 1, 2, 3, (illustrated in Fig.1) and 4. NPDR (grades 1, 2, 3) corresponds to the early-to-middle stage of DR and deals with a progressive microvascular disease characterized by small vessel damages and occlusions. PDR (grade 4) corresponds to the period of potential visual loss which is often due to a massive hemorrhage. Early identification and adequate treatment, particularly in the mild to moderate stage of NPDR, may slow the progression of DR, consequently preventing the establishment of diabetes-related visual impairments and blindness.
In the past years, deep learning has achieved great success in medical image analysis. Many supervised learning techniques based on convolutional neural networks have been proposed to tackle the automated DR grading task [10,11,4]. Nevertheless, these approaches rarely take advantage of longitudinal information. In this direction, Yan et al. [17] proposed to exploit a Siamese network with different pre-training and fusion schemes to detect the early stage of RD using longitudinal pairs of CFP acquired from the same patient. Further, self-supervised learning (SSL) held great promise as it can learn robust high-level representations by training on pretext tasks [1] before solving a supervised downstream task. Current self-supervised models are largely based on contrastive learning [6,3]. However, the choice of the pretext task to learn a good representation is not straightforward, and the application of contrastive learning to medical images is relatively limited. To tackle this, a self-supervised framework using lesion-based contrastive learning was employed for automated diabetic retinopathy (DR) grading [5]. More recently, a new pretext task has appeared for classification purposes in a longitudinal context. Rivail et al. [12] proposed a longitudinal self-supervised learning Siamese model trained to predict the time interval between two consecutive longitudinal retinal optical coherence tomography (OCT) acquisitions and thus capturing the disease progression. Yang et al. [18] proposed an auto-encoder named LSSL that takes two consecutive longitudinal scans as inputs. They added to the classic reconstruction term a time alignment term that forces the topology of the latent space to change in the direction of longitudinal changes. An extension of such principle was provided in [9]. To reach a smooth trajectory field, a dynamic graph in each training batch was computed to define a neighborhood in the latent space for each subject. The graph then connects nearby subjects and enforces their progression directions to be maximally aligned.
In this regard, we aim to use LSSL approaches to capture the disease progression to predict the change between no DR/mild NPDR (grade 0 or 1) and more severe DR (grade ≥2) through two consecutive follow-ups. To this end, we explore three methods incorporating current and prior examinations. Finally, a comprehensive evaluation is conducted by comparing these pipelines on the OPHDIAT dataset [7]. To the best of our knowledge, this work is the first to automatically assess the early DR severity changes between consecutive images using self-supervised learning applied in a longitudinal context.

Methods
In this work, we study the use of different longitudinal pretext tasks. We use the encoders trained with those pretext tasks as feature extractors embedded with longitudinal information. The aim is to predict the severity grade change from normal/mild NPDR to more severe DR between a pair of follow-up CFP images. Let X be the set of subject-specific image pairs for the collection of all CFP images. X contains all (x t , x s ) that are from the same subject with x t scanned before x s . These image pairs are then provided as inputs to an autoencoder (AE) structure (Fig.2c). The latent representations generated by the encoder are denoted by z t = F (x t ) and z s = F (x s ) where F is the encoder. From this encoder, we can define the ∆z = (z s − z t ) trajectory vector and then formulate ∆z (t,s) = (z s − z t )/∆t (t,s) as the normalized trajectory vector where ∆t (t,s) is the time interval between the two acquisitions. The decoder H uses the latent representation to reconstruct the input images such thatx t = H(z t ) and x s = H(z s ). E denotes the expected value. In what follows, three longitudinal self-supervised learning schemes are further described.

Longitudinal Siamese
The Siamese network takes the two image pairs (x t , x s ). These images are encoded into a compact representation (z t , z s ) by the encoder network, F . A feed forward neural network (denoted G) then predicts ∆t (t,s) , the time interval between the pair of CFP images (Fig.2a). The regression model is trained by minimizing the following L2 loss: G(z t , z s ) − ∆t (t,s) 2 2 . represents the longitudinal self-supervised learning which is composed of two independent modules, an AE and dense layers. The AE takes as input the pair of consecutive images and reconstruct the image pairs while the dense layer maps a dummy vector to the direction vector τ . The figure c) corresponds to the LNE, and takes as input the consecutive pairs and build a dynamic graph to align in a neighborhood the subjectspecific trajectory vector (∆z) and the pooled trajectory vector (∆h) that represents the local progression direction in latent space (green circle).

Longitudinal self-supervised learning
The longitudinal self-supervised learning (LSSL) exploits a standard AE. The AE is trained with a loss that forces the trajectory vector ∆z to be aligned with a direction that could rely in the latent space of the AE called τ . This direction is learned through a sub-network composed of single dense layers which map a dummy data into a vector τ ∈ Ω α , the dimension of the latent space. The highlevel representation of the network is illustrated in Fig.2b. Enforcing the AE to respect this constraint is equivalent to encouraging cos (∆z, τ ) to be close to 1, i.e. a zero-angle between τ and the direction of progression in the representation space.
Objective Function.

Longitudinal neighbourhood embedding
Longitudinal neighborhood embedding (LNE) is based on the LSSL framework.
The main difference is that a directed graph G is built in each training batch. A pair of sample (x t , x s ) serves as a node in the graph with node representation ∆z. For each node i, Euclidean distances to other nodes j = i are computed by The neighbour size (N nb ) is the closest nodes of node i form its 1-hop neighbourhood N i with edges connected to i. The adjacency matrix A for the directed graph (G) is then defined as: This matrix regularizes each node's representation by a longitudinal neighbourhood embedding ∆h pooled from the neighbours' representations. The neighborhood embedding for a node i is computed by: where O is the out-degree matrix of graph G, a diagonal matrix that describes the sum of the weights for outgoing edges at each node. They define θ ∆z,∆h the angle between ∆z and ∆h, and only incite cos(θ ∆z,∆h ) = 1, i.e., a zero-angle between the subject-specific trajectory vector and the pooled trajectory vector that represents the local progression direction in the latent space (Fig. 2c). Objective Function.

Dataset
The proposed models were trained and evaluated on OPHDIAT [ each image present in every examination, we select the image that gives the highest classification probability. We repeat this process until the selected images by the model converge to a fixed set of images per examination. From the set of selected images, we construct consecutive pairs for each patient and finally obtain 100,033 pairs of images from 26,483 patients. Only 6,690 (6.7%) pairs have severity grade changes from grade 0 or 1 to grade ≥ 2 against 93,343 (93.3%) pairs with severity changes that lie between grades 0 and 1. The resulting dataset exhibits the following proportions in gender (Male 52%, Female 48%) and diabetes type (type 2 69%, type 1 31%). This dataset was further divided into training (60%), validation (20%), and test (20%) based on subjects, i.e., images of a single subject belonged to the same split and in a way that preserves the same proportions of examples in each class as observed in the original dataset. Image pre-processing. Image registration is a fundamental pre-processing step for longitudinal analysis [15]. Therefore, using an affine transformation, we first conducted a registration step to align x t to x s . Images are then adaptively cropped to the width of the field of view (i.e., the eye area in the CPF image) and then resized to 256×256. A Gaussian filter estimates the background in each color channel to attenuate the strong intensity variations among the dataset which is then subtracted from the image. Finally, the field of view is eroded by 5% to eliminate illumination artifacts around the edges.

Experiments and Results
Implementation Details. As it was conducted in [18,9], we constructed a standard AE for all the compared methods to focus only on the presented methodology and make a fair comparison between approaches, with the hope that using advanced AE structures could lead to better encoding and generalization. In our basic architecture, we employed a stack of n pre-activated residual blocks, where n determines the depth of that scale for the encoder. In each res-block, the residual feature map was calculated using a series of three 3×3 convolutions, the first of which always halves the number of the feature maps employed at the present scale, such that the residual representations live on a lower-dimensional manifold. Our encoder comprises five levels; the first four levels are composed of two residual blocks, and the latter only one residual block. This provides a latent representation of size 64 × 4 × 4. The employed decoder is a reverse structure of the encoder. The different networks were trained for 100 epochs by the AdamW optimizer with a learning rate of 5 × 10 −4 , OneCycleLR as scheduler and a weight decay of 10 −5 , using an A6000 GPU with the PyTorch framework. The regularization weights were set to λ dir = 1.0 and λ rec = 5.0. A batch size of 64 was used for all models, and a neighbour size N nb = 5 and ∆z (t,s) were used for the LNE, as in the original paper [9].

Comparison of the approaches on the early change detection
We evaluate the LSSL encoders on detecting the severity grade change from normal/mild NPDR to more severe DR between a pair of follow-up CFP images. The classifier was constructed as the concatenation of the learned backbone (feature extractor) and a multi-layer perceptron (MLP). The MLP consists of two fully connected layers of dimensions 1024 and 64 with LeakyReLU activation followed by a last single perceptron. Receiving the flattened representation of the trajectory vector ∆z, the MLP predicts a score between 0 and 1. We compared the area under the receiver operating characteristics curve (AUC) and the accuracy (Acc) in Tab.1 for different pre-training strategies (from scratch, trained on LSSL methods, encoder from a standard AE). We also pre-trained on the OPHDIAT dateset (classification of the DR severity grade) to compare the LSSL pre-training strategies with a conventional pre-training method. The statistical significance was estimated using DeLong's t-test [13] to analyze and compare ROC curves. The results in Tab.1 and Fig.3 show the clear superiority of the LSSL encoder, with a statistical significance p-value < 2.2e-16. Due to class imbalance, the Longitudinal-siamese (L-siamese) have a high Acc while exhibiting a lower AUC than the baseline (trained from scratch).  [12] 0.8253 (0.8127-0.838) 0.9354 LSSL [18] 0.9624 (0.9593-0.9655) 0.8871 LNE [9] 0.9448 (0.9412-0.9486) 0.8646 Table 1: Comparison of the approach on the early change detection with the frozen encoder.

Norm of trajectory vector analyze
We constructed different histograms in Fig.4 representing the mean value of the norm of the trajectory vector with respect to both diabetes type and progression type. According to Fig.4, only the models with the direction alignment loss term are able to capture the change detection in the DR relative to the longitudinal representation. Therefore, we observe in the histogram that the trajectory vector (∆z) is able to dissociate the two types of diabetes (t-test p-value < 0.01) and change detection (t-test p-value < 0.01). For the diabetes type, a specific study [2] about the OPHDIAT dataset indicates that the DR progression is faster for patients with type 1 diabetes. Based on the fact that the ∆z can be seen as a Fig. 4: Mean of the trajectory vector norm for the different self-supervised method used relative speed, this observation agrees with the histogram plot of the mean of the ∆z norm represented in Fig. 4. We also observed that the norm of the ∆z vector is lower for the normal stage of the DR than for mild NPDR to more severe. This was expected because only the methods with a direction alignment term in their objective explicitly modeled longitudinal effects, resulting in more informative ∆z. This also implies that simply computing the trajectory vector itself is not enough to force the representation to capture the longitudinal change.

Discussion
We applied different LSSL techniques to encode diabetic retinopathy (DR) progression. The accuracy boost, relative to a network trained from scratch or transferred from conventional tasks, demonstrates that longitudinal pre-trained selfsupervised representation learns clinically meaningful information. Concerning the limitations of the LSSL methods, we first observe that the models with no time alignment loss perform poorly and provide no evidence of disease progression encoding. Also, we report for the LNE that the normalized trajectory vector for some pairs, that have a large time between examinations, is almost all zeros, which results in a non-informative representation. This could explain the difference between the LSSL and LNE prediction performances. Moreover, during the LSSL and the LNE training, we often faced a plateau with the direction loss alignment. Therefore, we also claim that intensive work should be done regarding the choice of the hyperparameters : constant weights for the losses, latent space size (Ω α ), neighbour size (N nb ). The results concerning quantifying the encoding of the disease progression from the models trained with a time direction alignment are encouraging but not totally clear. As it was mentioned in [16], one limitation of the LSSL approach pretains to the cosine loss (direction alignment term from equations (1,2) used to encode the longitudinal progression in a specific direction in the latent space and learned while training. The loss only focuses on the correlation with the disease progression timeline but not disentanglement of the disease progression itself. Therefore, a more in-depth analysis of the latent space is required to evaluate if the trajectory vector could be used to find a specific progression trajectory according to patient characteristics (diabetes types, age, DR severity). The pairing and the registration are critical steps in the longitudinal study. As it was previously mentioned, by using a better registration method and exploiting different fusion schemes and backbone architectures, we could get enriched latent representation and, thus, hopefully, better results. Also, the frozen encoders could be transferred to other types of longitudinal problems. In summary, LSSL techniques are quite promising: preliminary results are encouraging, and we expect further improvements.