Enforcing artificial neural network in the early detection of diabetic retinopathy OCTA images analysed by multifractal geometry

Diabetic retinopathy (DR) is one of the leading causes of vision loss. It causes neovascularization with blocking the regular small blood vessels. Early detection helps the ophthalmologist in patient treatment and prevents or delays vision loss. In this work, multifractal analysis has been used in some details to automate the diagnosis of diabetic without diabetic retinopathy and non-proliferative DR. Concerning using number of multifractal geometrical methods, as a necessary second step the enforcement of the sophisticated artificial neural network has been consultant in order to improve the accuracy of the obtained results. Patients and methods: Thirty normal cases’ eyes, 30 diabetic without DR patients’ eyes and 30 non-proliferative diabetic retinopathy (mild to moderate) eyes are exposed to optical coherence tomography angiography (OCTA) to get image superficial layer of macula for all cases. These images were approved in Ophthalmology Center in Mansoura University, Egypt, and medically were diagnosed by the ophthalmologists. We extract the most changeable features that associated to the morphological retinal vascular network alternations. The seven extracted features are related to the multifractal analysis results, which describe the vascular network architecture and gaps distribution. A supervised Artificial Neural Network (ANN) is used to classify the images into three categories: normal, diabetic without diabetic retinopathy and non-proliferative DR. Results: The human retinal blood vascular network architecture is found to be a fractal system. Multifractal geometry describes the irregularity and gaps distribution in the retina. We extracted seven features from the studied images. The features were the generalized dimensions D0 , D1 , D2 , α at the maximum f(α) singularity spectrum, the spectrum width, the spectrum symmetrical shift point and lacunarity. The ANN obtains a single value decision with classification accuracy 97.78%, with minimum sensitivity 96.67%. Conclusion: Early stages of DR could be noninvasively detected using high-resolution OCTA images that were analysed by multifractal geometry parameterization and implemented by the sophisticated artificial neural network with classification accuracy 96.67%. This approach could promote risk stratification for the decision of early diagnosis of diabetic retinopathy.


Introduction
Diabetic retinopathy remains a frightening prospect to patients and discomfit physicians. The pathophysiology of diabetic retinopathy is mainly microangiopathy in which the pathology of small blood vessels is affected by occlusion and leakage [1].
Though fluorescein angiography is considered gold standard for diagnosis of diabetic retinopathy [2,3]. Optical coherence tomography angiography nowadays has an upper hand in the diagnosis of early stages of non-proliferative diabetic retinopathy because of its advantages as the high resolution specially in SS-OCT [4,5]. OCT is noninvasive device and detects the early microaneurysm in the retinal microvascular network during the progression of diabetic retinopathy. Avoiding time-consuming, venipuncture, in fluorescein angiography that also reporting anaphylaxis and death cases related to fluorescein injections that have been documented, despite their rarity [6][7][8][9].
The identification of new targets and therapeutic strategies for providing the means to better manage DR, including diagnosis at early stages in the course of the disease as blindness is preventable in diabetics if it is early diagnosed. Vascular network in the superficial layer is found to be nonlinear system from mathematical point of view, as biological systems are known to be unpredictable and it is often difficult to predict the outcome or response of the body to a change in circumstances. Fractals show a system that is seemingly based on very simple principles but leads to very complex structures. Fractals are adequate to characterize the chaotic behaviour for nonlinear dynamical system as well as are effective describing biological systems [10].
As a good example, the retinal blood vascular network is found to be a fractal system. Fractal geometry describes the irregularity or fragmented shape of natural features where Euclidean geometry failed.
As the fractal dimension describes how much space is filled but it does not indicate how the space is filled, the lacunarity is a parameter that describes the distribution of gap sizes throughout the object established in an image, and it is capable of identifying different fractal structures that have approximately the same fractal dimension [10]. Previous studies have employed the lacunarity to identify the occlusion of the artery and retinal vein, they used it to diagnose amblyopic eyes [11].
In fact, the retinal blood vascular network is considered multifractal pattern as its different regions have different fractal properties, which can be characterized by a hierarchy of exponents rather than a single fractal dimension [12,13].
Artificial intelligence techniques have become the most important methods that used in medical diagnosis. Some researchers extract the retinal blood vessels with effective features using the combination of supervised and unsupervised machine learning techniques [14]. Others provide a technique for detecting various stages of DR [15]. Improving and speed-up the Convolutional Neural Network (CNN) training algorithm for medical image analysis has been done in [16]. A supervised segmentation technique used a deep neural network on a large sample [17]. The application of artificial intelligence in retina and diabetic retinopathy was found in [18,19]. Therefore, the main motive is to assess the macular vascular network using OCTA images as a multifractal system either it is normal, diabetic without DR or nonproliferative DR (NPDR) mild to moderate cases. Use the generalized dimensions (the capacity box-counting dimension, the information dimension and the correlation dimension), lacunarity and the characteristics of the singularity spectrum as numerical features for medical diagnosis using of an artificial neural network (ANN) as a classifier. The ANN summarizes the information into a single value, which is very suitable in terms of feature reduction and for statistical analysis. It is convenient for automatic assessments without requiring extensive user involvement.

Retinal images
The retrospective observational case series were appro ved by Ophthalmology Center in Mansoura University, Egypt. We have 90 cases divided into 30 eyes of healthy, 30 eyes of diabetic without DR and 30 eyes of NPDR, were included in this paper. The subjects were between 45 and 63 years old. These subjects were determined via clinical examination by retinal specialist. The subjects' images have good resolution and were carefully chosen from about 200 images.

Study design
All subjects were imaged by using Triton Topcon Swept-Source OCTA, with acquisition speed of 100,000 A-scans/s, 1050 nm wavelength, B-scan acquisition approx. 0.016 s and 3D scan acquisition time approx. 0.65 s. Imaging was performed using angiographic 3*3 mm scan pattern, consists of 4 repeated B-scans of 320 A-scans each at 320 raster positions, centred at the fovea centre.

Multifractal methods
As is well known, the irregular forms and multiple hierarchical spatial and dynamical structures found in biology and medicine, like for instance vascular network geometry of the human retina among other, are not fully understood within the framework of traditional methods. In fact, , a biological entity has a fractal geometry, hence their dimensions is always expressed by noninteger fractal dimension.
For example, in order to understand and explain the complexity of the geometry of the human retina, one can use multifractal analysis as an complementary part of a screening tool for early detection of retinal vascular diseases. The main cause for using multifractal analysis to describe human retina is the generalization of fractal concept to include intricate structures with more than one scaling exponent. There are many physical parameters which can be applied to describe the multifractal structures like, e.g. box-counting dimension, information or entropy dimension and correlation dimension. In order to appreciate the difference between them, we briefly give some short notes concerning them.

Generalized dimensions 3.1.1.1. Box-counting dimension D F . In principle,
one can calculate the box-counting dimension by superimposing a regular grid of pixels of length δ on the object and count the numbers of "occupied" pixels. This procedure is repeated using different values of δ. In other words, one can define D B as the real number, such that the number of grids N(δ) with scale δ that is needed to cover an object grows with 1 δ D B as δ → 0.
In the case of real data set, the mathematical expression for D B is given by

Information dimension D i .
As is well known, the general Shannon entropy or the information entr opy H(δ) at a scale δ is given by where N is the total number of points in a given set subject to is the number of occupied boxes of size δ and f i is the relative frequency distribution. It worth mention that for a given fractal set, one finals the following expression: where k is a given constant. Therefore, for nonuniformly distributed patterns, the information dimension will subsequently take the following form: From thermodynamics point of view, entropy is the measure of the uncertainty of a random event. This means, if an event happens very often, it provides less information leading to lower state of entropy. In general, it plays an important rule in the analysis of nonlinear dynamics systems, especially in describing the loss of information as a chaotic system evolves in time. Similar to the mathematical expression of fractal dimension, the information dimension is defined as Clearly, above expression provides that D i = D F for uniform point patterns, however, for non-uniform point patterns, one finds D i < D F .

Correlation dimension D C . Principally, this
mentioned dimension is well adapted to the characterization of a spatial clustering in non-uniform point pattern. In other words, correlation dimension measures the fractal dimension giving the relationship between two pixels inside a region under consideration. In essence, the correlation dimension is widely used to measure the fractal dimension to strange attractors for nonlinear dynamical system undergoing chaotic dynamics.
In such case as defined, the correlation function of a point pattern is calculated via the following expression of correlation interval N: where C i (δ) is given by is the number of distinction pairs in circle of radius δ, centred on the ith of N points and d ij = |x i − x j | is the Euclidean distance between the ith and jth points. In the above expression, θ (µ) is well known as Heaviside function. As is well known, θ(µ) = 0 and θ(µ) = 1 for x < 0 and x ≥ 0, respectively. From the above expression, it is easy to define D C using the following formula: Noting that, proceeding the formal definitions leading us to for uniform point patterns, while for non-uniform distribution turn to be Therefore, it should be clear that the correlation information and box-counting dimension are significantly different for human retina.

Singularity spectrum f(α)
Another way to compute the multifractal behaviour is the singularity spectrum. The multifractal spectrum f(α) shows the singularity exponent distribution. It quantifies the nature degree of nonlinearity in a compact way. From the generalized fractal dimension, the exponents τ and q can written as Therefore where α is a function of q by solving the equation If the exponents τ (q) are known, then the multifractal spectrum can be written as The generalized dimension and f(α) spectrum is shown in Figure 1.

Artificial neural network(ANN) method
Neural network takes inspiration from the learning process occurring in human brains. It is a biologically inspired programming paradigm, which enables a computer to learn from observational data. It is one type of expert system algorithms that can be used to perform nonlinear statistical modelling, the most commonly used method for developing predictive models for dichotomous outcomes in medicine.
The main element in the network architecture called neuron, which performs a specific task according to its activation function, as simple neuron can be shown in Figure 2.
This model can be expressed as The activation function F can be selected from different function varieties according to the required tasks. The network structure is organized as layer of neurons arranged in parallel layers and are connected to each other by weighted connections. The ANN is shown in Figure 3.
The basic architecture consists of three types of neuron layers: input, hidden and output layers. In feedforward networks, the signal flow is from input to output units, strictly in a feed-forward direction. Feedforward neural networks are usually trained by the original back propagation algorithm where training is usually carried out by iterative updating of weights based on the error signal. There are several other neural network architectures [27][28][29] depending on the properties and requirement of the application.
The objective of training the neural network is to minimize the cumulative error, which can be expressed as Where N is the total number of neurons, y i is the target output and F(net i ) is the actual ith neuron output. According to the error the network weights w should be modified as where α is the learning parameter and is constant, the sign may be positive or negative according to the resulted output and errors.
There are several learning methods in neural networks [30,31] and may be classified into (1) supervised learning where pair input and target vectors are presented together to the network, with calculating the errors between the desired and actual response for each node in the output layer are found, then updating the network weights. (2) Unsupervised learning, no labels are given to the learning algorithm, leaving it on its own to find structure in its input. The system is supposed to discover statistically salient features of the input population. (3) Reinforcement learning, software-defined agents to find the optimal way to attain their goals or improve the task performance. Many researchers have used neural network in classification applications [32][33][34][35][36][37], clustering [38][39][40], pattern recognition [41][42][43][44], etc.

Image processing software for OCTA images
The OCTA images are processed by using a custom developed MATLAB program (MATLAB v.9.4 for R2018a Mathworks, MA, USA) and Fiji-based Image-J software. The processing method is clear through the following two steps: Step 1: The images are cropped, with applying image processing tools for enhancing the resulted images.
Step 2: Produce information about how "tube-like" in each point in the image and emphasizes the vessels cylindrical structure. This is useful as a preprocessing step for tracing neurons or blood vessels. The "Tubeness" plugin works by Gaussian smoothing the image followed by multiplying the two lowest eigenvalues of the local Hessian matrix [45].

The diabetic retinopathy detection method
After the images have been processed, multifractal analysis and lacunarity measurement are calculated. Seven features are extracted from each image to describe the morphological structures of the retinal microvascular system. The seven features are: (1) Box-counting dimension D F (2) Information Dimension D i (3) Correlation Dimension D C (4) α at the maximum value of the singularity spectrum f(α) (5) Singularity spectrum width (6) Shift in the singularity spectrum symmetrical axis (7) Lacunarity.
The proposed methodology flow chart is shown in Figure 4.
As shown in Figure 4, the flowchart summarizes the methodology steps. The ANN was trained until it reaches the minimum trained error E min , then the system is ready for classification.

Image datasets
We analysed 90 eye images, 30 eyes for normal cases, 30 eyes for diabetic without DR and 30 eyes for NPDR cases. Table 1 shows the demographic characteristics of the studied subjects. Table 1 shows the demographic characteristics of the examined subjects. The subjects were classified and examined. Performing statistical analysis is to assess the significance of the studied subject characteristics. In our research, the data were nonsignificant if the associated P value was greater than 0.05. The average ages of the subjects and their standard deviation were shown with P > 0.41904, which is statistically nonsignificant. The female/male numbers were 42/48, which has P > 0.3729. The laterality RT/LT numbers were 47/43 with P > 0.5525. So, the studied subject characteristics were statistically nonsignificant. Therefore, it can be said that the demographic characteristics data of the cases examined have not any impact on the accuracy of the results. Figures 5-7 show four sample images for each normal, diabetic without DR and NPDR subjects respectively.
Each figure consists of four rows. The first row contains the source image, the second row contains the processed images with vessels cylindrical structure. The multifractal analysis results were illustrated in third and  fourth rows. The third row has generalized dimensions, while the fourth row contains singularity spectrum image. Figures 5-7 show the D(Q) curves in normal, diabetic without DR and NPDR cases. Both diabetic cases' curves have some slight changes but not sufficient to have a single parameter as accurate decision. Also, the figures showed the variation of entropic function f (α) with the parameter α in order to visualize the distribution of space occupied by vascular structure for all cases. Interestingly, one can infer from this mentioned Figure 8 that there is an appreciate difference between all cases under consideration. It is obvious that there is a noticeable shift to the right for all cases. Hence, the blue curves for normal subjects started earlier and ended before a given value of parameter α = 2.5. While the red curves for diabetic without DR starts nearly with normal subjects and maximum ended nearly at α = 2.5.
However, the black curves for NPDR started after the other two cases and approximately exceeded α = 2.5. Whenever the diabetic status increased, the more shift the curve to the right.
In addition, one can note a clear symmetrical behaviour for entropic function f(α) corresponding normal cases, hence the average symmetrical shift for normal case is 0.085. While the average symmetrical shift for diabetic without DR is 0.123 and 0.205 for NPDR, i.e. the normal cases spectrum is more symmetrical shape than other diabetic cases, while NPDR is the less symmetrical shift spectrum. As well as, the f(α) spectrum width is lightly wider more in NPDR cases than normal and diabetic without DR cases.
In order to classify the obtained data in more accurate manner, it is necessary to find a set of parameters, which in turn leads to decision making in the classification of pathological and normal subjects.

Classification with artificial neural network
OCTA 90 eye images for normal, diabetics without DR and NPDR subjects were processed to obtain the set of parameters that enable us to classify the subjects. The most variables should be used that contain sufficient information about the microvascular blood vessels, the homogeneity of the vessels system, the distribution of gaps through the vessels structure and the random distribution of the retinal blood vessels. Therefore, seven parameters were chosen to ensure a complete information about the microvascular system. The chosen parameters are the generalized dimensions D F , D I and D C correspond to box-counting dimension, information dimension and correlation dimension respectively, lacunarity, α at the location of maximum f(α) spectrum, the width and asymmetry of the spectrum. Table 2 shows 7 extracted features for 12 sample images, 4 for each case. Now, after obtaining these features, we need to manipulate these parameters to enable us for obtaining a final decision according to the image status. One of the most suitable classifier tool is artificial intelligence.
In our research, a classifier neural network was created as shown in Figure 9 with four layers. The input layer has seven neurons, while there are two hidden layers with five and four neurons respectively and an output layer has three neurons. The inputs were the obtained parameters from the patient image analysis D F , D I , D C , lacunarity, α at the location of maximum f(α) spectrum, the width and asymmetry of the spectrum. The classification and the diagnosis tasks were done by using ANN. A typical ANN has a multilayered structure such as feed-forward networks with backpropagation learning algorithm. We use 60 sample images in training the neural network while the other 30 images were used in network testing. The classification performance and accuracy will be discussed in the next section.

Performance measures
To find out the quality of the proposed system in terms of the accuracy of the results and its effectiveness, several performance measures can be used to assess the system effectiveness. We use in our research classification accuracy, confusion matrix, sensitivity, specificity and precision. After the network has been trained, it was tested by 30 new subject eyes, 10 for normal eye images, 10 for diabetic without DR images and 10 for NPDR images Tables 3 and 4 summarize the performance measure results.
From Tables 3 and 4, the proposed classifier app roach produces better results. Hence, the statistical parameters confirm that the minimum accuracy is within 97.78% for 90 eyes study cases. In the classification of normal versus the diabetic subjects, the sensitivity, specificity and precision measure achieved by the proposed system are 100%, 97.78% and 93.75% respectively. The system achieves low misclassification error with 2.22%. While, in the classification of DWDR  versus NPDR, the sensitivity, specificity and precision measures are equal to 96.67%, with misclassification error 3.33%. These results can easily be enhanced by increasing number of studied subjects.

Conclusion
Multifractal analysis gives us an open window to obtain a more and accurate description of the human retina network than Euclidean geometry. In principle, it can provide us an accessible and effective neat technique to characterize and investigate the human retina. This could give us a recommendation that multifractal analysis becomes an essential screening tool for the early detection of retinal diseases.
From the entropic function f(α) and correlation dimensions, we get good ideas about the presence of gaps and detection of bifurcation points. There is an appreciate difference between both cases under consideration. Seven features were chosen to ensure a complete information about microvascular system. Then, we manipulate the extracted features to enable us for obtaining a final decision according to the image status.
Artificial intelligence uses complex algorithms and software to simulate the human brain for information processing and decision making. Also, the main essence of this study is to employ artificial neural network by using the seven extracted features as network inputs for mathematically analysing the digit images. The accurate decision is obtained by testing the network with several different images between normal eyes, diabetic without DR and Non-Proliferative DR patients. The supposed system accuracy reaches 97.78% with low datasets. The obtained results could be improved by increasing the number of trained image samples. By the numerical early diagnosis, blindness will be preventable or delayed in diabetics.