Retinal OCT Texture Analysis for Differentiating Healthy Controls from Multiple Sclerosis (MS) with/without Optic Neuritis.

Multiple sclerosis (MS) is an inflammatory disease damaging the myelin sheath in the central and peripheral nervous system in the brain and spinal cord. Optic Neuritis (ON) is one of the most prevalent ocular demonstrations of MS. The current diagnosis protocol for MS is MRI, but newer modalities like Optical Coherence Tomography (OCT) are already of interest in early detection and progression analysis. OCT reveals the symptoms of MS in the Central Nervous System (CNS) through cross-sectional images from neural retinal layers. Previous works on OCT were mostly focused on the thickness of retinal layers; however, texture features seem also to have information in this regard. In this research, we introduce a new pipeline that constructs layer-stacked (LS) images containing data from each specific layer. A variety of texture features are then extracted from LS images to differentiate between healthy controls and ON/None-ON MS cases. Furthermore, the definition of texture extraction methods is tailored for this application. After performing a vast survey on available texture analysis methods, a treasury of powerful features is collected in this paper. As a primary work, this paper shows the ability of such features in the diagnosis of HC and MS (ON and None-ON) cases. Our findings show that the texture features are powerful to diagnose MS cases. Furthermore, adding information of conventional thickness values to texture features improves considerably the discrimination between most of the target groups including HC vs. MS, HC vs. MS-None-ON, and HC vs. MS-ON.

1. Introduction

Multiple sclerosis (MS) is an inflammatory disease damaging the myelin sheath in the central and peripheral nervous system in the brain and spinal cord. This disease causes the immune system to attack one or more proteins of the myelin structure and disrupts the ability of the nervous system to communicate and therefore brings about many physical signs and symptoms [1]. Those suffering from MS show neurological symptoms including disorders in the autonomic, visual, motor, and sensory nervous system [2]. Optic Neuritis (ON) is a common eye problem where inflammation or demyelination affects the optic nerve. It occurs when inflammation damages the optic nerve, a bundle of nerve fibers that transmits visual information from the eye to the brain. Signs and symptoms of ON can be the first indication of MS, or they can occur later in the course of MS. Not everyone who experiences ON goes on to develop further symptoms of MS, but a significant proportion does [3].

The current diagnosis protocol for MS is Magnetic Resonance Imaging (MRI); however, researchers are already looking for substitute methods to overcome MRI limitations like high cost, late-stage diagnosis, and inaccurate signs due to aging rather than MS [4]. The effects of MS on the Central Nervous System (CNS) make the retinal nerve fiber layer (RNFL) a proper candidate for being imaged instead of brain MRIs. The thickness of RNFL can be used to assess the existence of any damage in the CNS. Moreover, RNFL is considered as one of the main retinal layers. The role of the remaining layers is not exactly known in the case of MS and needs more investigations.

Optical Coherence Tomography (OCT) is a noninvasive imaging modality to take cross-sectional images of biological tissues. Retinal OCT provides information on symptoms of many eye diseases such as macular degeneration, glaucoma, and diabetic retinopathy and helps ophthalmologists to diagnose and treat such diseases in a timely manner [5-7]. Parisi et al.'s study in 1998 on the diagnosis of MS using retinal OCT was the first work in this field [8]. He investigated whether there is a relationship between RNFL thickness and visual pathway function in patients with MS. Since then, a great deal of research has been done on thickness changes in different retinal layers and the possibility of their use for diagnosing MS. Petzold et al. in 2017 prepared a survey covering this topic and reviewed 110 articles from 1991 to 2016 and provided a good overview of the subject [9].

In recent years, in the field of retinal OCT image processing, much attention has been paid to extracting and using texture features of layers while these types of features have not yet been used widely for MS diagnosis and there are few works addressing this issue. As an example, Varga et al. in 2015 had a study investigating the differences in texture descriptors and optical properties of retinal tissue layers in patients with MS and evaluated their usefulness in the detection of neurodegenerative changes using OCT image segmentation [10]. The term texture in image processing and machine vision refers to the amount, type, and distribution of pixel brightness throughout the image along with the texture of the image [11]. Researchers have defined it as “A texture area in an image can be constructed with an irregular and varied spatial distribution of the intensity of the brightness or color [12].” In this regard, four general categories named statistical, structural, signal processing-based, and model-based features are usually used [11]. In this study, we want to examine the texture of OCT images, and we suspect that changes in the texture of the layers must occur before the thickness changes. It seems that the deterioration of axons in the retinal nerve fiber layer and changes in the texture layers can be determined by the noninvasive OCT method, making them possible to be used as a complementary diagnostic tool in addition to the existing methods for early detection of recurrent MS-ON and MS-None-ON [13].

Here is an overview of the literature investigating texture features in OCT images. In 2007, Baroni et al. investigated the possibility of discriminating retinal OCT image layers in texture processing using Grey-Level Cooccurrence Matrix (GLCM) feature extraction [14]. In 2014, Anantrasirichai et al. presented a new method for extracting the texture of OCT retinal images in glaucoma [15]. In 2018, Sawyer et al. examined the possibility of using texture analysis to classify ovarian OCT images [16]. In 2019, Nunes et al. used texture analysis of OCT data to define new biomarkers for MS, of course, only on one specific retinal layer [17].

The rest of this paper is as follows. The proposed method for texture extraction of retinal OCT layers is described in Section 2. The performance of the method is evaluated and discussed in Section 3. Finally, Section 4 presents the conclusions of the study.

2. Material and Method

2.1. Database

The data in this study is obtained from Spectralis Heidelberg HRA+OCT device in Faiz Hospital and Sadra Ophthalmology Center, Isfahan, Iran. The size of each B-scan is 496 × 480 pixels. For some subjects, data contains 19 B-scans, and for others, it includes 25 B-scans. OCT data includes 36 health control (HC) eyes and 39 patients suffering from MS (20 eyes suffering from MS with no history of ON (MS-None-ON) and 19 eyes suffering from MS with a history of ON (MS-ON)). HC and patient populations have matched gender and age approximately. A summary data flow diagram is presented in Figure 1.

Figure 1

Detailed structure of the data.

2.2. Algorithm Flow

The workflow of the proposed method is shown in Figure 2. The first step is the preprocessing block in which the retinal delineation [18] is used to extract the layers. In the second block, layer-stacked (LS) images are created by stacking each specific layer from all B-scans of one subject. The third block is applied for masking the images as input to the next feature extraction block. Five different groups of texture features are utilized in this step. In the following, the most effective features are selected based on p value for distinguishing HC, MS-ON, and MS-None-ON population from retinal OCT layers around the fovea. Finally, in the last step, a classification between HC and abnormal population is performed. Each block of the proposed algorithm flow is elaborated below.

Figure 2

Algorithm flow of the proposed method.

The sample output of preprocessing block is shown in Figure 3. To construct layer-stacked images, we consider that data for each subject consists of a number of B-scans, and each B-scan contains 10 layers, locations of which are obtained in preprocessing step. Accordingly, we construct 10 layer-stacked images by cutting and stacking each individual layer from all B-scans of one subject (Figures 4 and 5). A sample of layer-stacked images is demonstrated in Figure 6.

Figure 3

Interretinal layers (a sample output of preprocessing block).

Figure 4

Sample of individual layers in one B-scan.

Figure 5

Construction process for layer-stacked (LS) images.

Figure 6

Layer-stacked (LS) images corresponding to each retinal layer. (a) First layers of all B-scans. (b) Second layers of all B-scans. (c) Third layers of all B-scans. (d) Fourth layers of all B-scans. (e) Fifth layers of all B-scans. (f) Sixth layers of all B-scans. (g) Seventh layers of all B-scans. (h) Eighth layers of all B-scans. (i) Ninth layers of all B-scans. (j) Tenth layers of all B-scans.

During texture calculation, boundary points in layer-stacked images have synthetic contrast which may fool the feature extraction method and lead to incorrect and outlier values. To solve this problem, an eliminating mask is developed to ignore pixels located on both sides of each individual layer.

Feature extraction is then performed on masked layer-stacked images. The features used in our work are GLCM, Local Binary Pattern (LBP), Local Directional Pattern (LDP), Local Optimal Oriented Pattern (LOOP), and fractal dimension. Finally, discriminant features are fed into Support Vector Machine (SVM) and Linear Discriminant Analysis (LDA) classifiers for differentiating between HC and MS cases.

2.3. Texture Feature Extraction

Investigating texture features is an efficient way to characterize various properties, such as structure, orientation, roughness, smoothness, or regularity of an image. Extracting features from masked layer-stacked images, we apply two categories of texture features including original and modified features.

2.3.1. Original Features

Different texture analysis methods are utilized in this research and elaborated in the next subsections. A set of features are then extracted according to Table 1.

Table 1

List of utilized features.

Texture analysis method	Features	Description	Texture analysis method	Features	Description
GLCM	Energy	Provides the sum of squared elements in the GLCM. It has values between 0 and 1	GLCM	Difference variance	Measures the dispersion (with regard to the mean) of the grey-level difference distribution of the image
	Entropy	Measure of randomness that can be used to characterize the texture of an image		Difference entropy	Measures the disorder related to the grey-level difference distribution of the image
	Contrast	A measure of intensity contrast between a pixel and its neighbor over the whole image		Maximum probability	Measures the maximum likelihood of producing the pixels of interest
	Homogeneity	Measures the closeness of the distribution of elements in the GLCM to the GLCM diagonal		IMC1	Measure of dependency between two random variables
	Correlation	A measure of how correlated a pixel is to its neighbor over the whole image		IMC2
	Sum of squares	Measures the dispersion (with regard to the mean) of the grey-level distribution	LBPLDPLOOP	Mean	Measures mean and standard deviation of histograms
	Cluster shade	Characterizes the tendency of clustering of the pixels in the region of interest
	Cluster prominence			Standard deviation
	Dissimilarity	A measure of distance between pairs of pixels in the region of interest		Dynamic range	Measures the ratio between the largest and smallest values
	Autocorrelation	Represents the degree of similarity between a given time series and a lagged version of it		Kurtosis	Measure of the “tailedness” of the probability distribution
	Sum average	Measures the mean of the grey-level sum distribution of the image		Skewness	Measure of the asymmetry of the probability distribution
	Sum entropy	Measures the disorder related to the grey-level sum distribution of the image	Fractal dimension	Mean	Measures mean and standard deviation of the box-counting method
	Sum variance	Measures the dispersion (with regard to the mean) of the grey-level sum distribution of the image		Standard deviation
	Inverse difference	A measure of local homogeneity of an image

(1) Grey-Level Cooccurrence Matrix. GLCM describes the spatial relationship between each intensity tone by considering changes between grey levels i and j at a particular displacement distance d and at a particular angle θ [15]. Here, we use a 256 quantization level and the distance is selected as one pixel with four distinct orientations (0, 45, 90, and 135 degrees). Furthermore, conditions of those pixels with 180 degrees in difference are considered to be the same.

(2) Local Binary Pattern. LBP is a method for describing the texture characteristics introduced in 1990 [19]. LBP compares the intensity of each pixel with neighboring pixels and determines the output value based on equation (1), where P is the number of neighboring points that are chosen, i.e., 8, i is the intensity of the neighborhood points, and  i is the intensity of the central point.  LBP calculates the output of LBP for P neighboring points.

(3) Local Directional Pattern. A more robust to noise modified version of LBP is LDP which computes directional components for each pixel with Kirsch kernels and provides a measure of the strength of intensity variation in those directions [20]. For each central pixel located at (x, y) with intensity i, eight rotated versions of the Kirsch edge detector should be applied on neighboring pixels with intensities i n = 0, 1, ⋯, 7. Eight corresponding responses of the Kirsch masks are m n = 0, 1, ⋯, 7. m is the kth highest Kirsch activation, and all the neighboring pixels having Kirsch response higher than m are assigned 1, and others 0. Then, the LDP value for the pixel (x, y) is given by

(4) Local Oriented Optimization Pattern. LOOP offers a nonlinear combination of LBP and LDP that overcomes their individual problems while maintaining the strengths of each. Compared to LDP, LOOP assigns an exponential weight w to each of neighboring pixels. w is a digit between 0 and 7, according to the rank of the magnitude of m among the 8 Kirsch mask outputs [21]. The value of the LOOP in (x, y) is given by

(5) Fractal Analysis. Images with self-similarity characteristics are called fractal. The box-counting analysis is an appropriate method of fractal dimension estimation for images with or without self-similarity [22]. We have a basic equation for calculating fractal dimension given by equation (4), in which N is the number of boxes that cover the pattern, and r is the magnification or inverse value of the box size.

A higher slope means that the object is more fractal, i.e., reduction in the size of the box reveals more complexity. The lower slope means that the object is closer to the straight line, i.e., less fractal, and the amount of details does not increase rapidly with increasing magnification.

2.3.2. Modified Features

Inserting zero values by masking layer-stacked images (third block in Figure 2) causes unwanted strip artifact. In order to solve this problem, we modify the output of the abovementioned texture analysis methods, to extract more accurate features. A list of used abbreviations in this paper and their explanations is shown in Table 2.

Table 2

A list of used abbreviations and their explanations.

Abbreviation	Explanation
MS	Multiple sclerosis
ON	Optic neuritis
OCT	Optical coherence tomography
CNS	Central nervous system
LSI	Layer-stacked images
MS-ON	Multiple sclerosis with optic neuritis
MS-None ON	Multiple sclerosis without optic neuritis
MRI	Magnetic resonance imaging
RNFL	Retinal nerve fiber layer
GLCM	Grey-level cooccurrence matrix
HC	Health control
LBP	Local binary pattern
LDP	Local directional pattern
LOOP	Local optimal oriented pattern
SVM	Support vector machine
LDA	Linear discriminant analysis
FD	Fractal dimension
LS image	Layer-stacked image

For GLCM, the first row and column of the output matrix (which represent unwanted zero pixels) are eliminated. The GLCM features listed in Table 1 can then be calculated fromwhere element [i, j] of the matrix is generated by counting the number of times a pixel with value i is adjacent to a pixel with value j and then dividing the entire matrix by the total number of such comparisons made. Each entry is therefore considered to be the probability that a pixel with value i will be found adjacent to a pixel of value j. μ‚ μ, σ, and σ are means and standard deviations. p and p are partial probability density functions. x and y are the coordinates (row and column) of an entry in the cooccurrence matrix, and p(i) is the probability of cooccurrence matrix coordinates summing to x + y. HX and HY are the entropies of p and p. Finally, HXY, HXY1, and HXY2 are shown in

In LBP, LDP, and LOOP methods, the features in Table 1 should be extracted from the histogram of the output. To solve the same problem of unwanted strip artifact, the first column of the histogram (which represent unwanted zero pixels) is eliminated. Finally, five statistical features including mean, standard deviation, dynamic range, kurtosis, and skewness are extracted.

The last category of texture analysis methods to be considered is fractal analysis. Here, we remove the black background above the layer-stacked images before performing the masking step. The mean and standard deviation of the fractal dimensions for each image is then reported.

2.4. Feature Selection and Classification

To handle the course of dimensionality problem caused by small number of available data compared to bunch of calculated features, more significant features are selected based on t-test and Bonferroni correction. The Bonferroni correction is an adjustment made to p values when several dependent or independent statistical tests are being performed simultaneously on a single data set. To perform a Bonferroni correction, the critical p value (α) is divided by the number of comparisons being made.

Here, considering that the majority of the features have been extracted from the GLCM matrix and this matrix has produced the features in four different angles, according to Bonferroni correction, the value of meaningful p value (p < 0.005) is divided by 4 and p < 0.001 is considered as a significant level for cut-off. After Bonferroni adjustment for multiple comparisons, features with p < 0.001 are selected as significant features. Then, two classification models, SVM and LDA, are utilized for differentiating between four possible groups including HC vs. MS, HC vs. MS-ON, HC vs. MS-None-ON, and MS-ON vs. MS-None-ON. A 10-fold cross-validation is used to evaluate accuracy, for each classification model.

3. Result

3.1. Feature Analysis

To evaluate the statistical significance of the extracted features, the t-test is used to identify which features show significant differences between healthy and MS (ON and None-ON) cases (Table 3). The p values indicate the test rejection of the null hypothesis at the 5% significance level, considering the Bonferroni correction (p value < 0.001). Frequencies of significant selected features for each retinal layer are also presented in Table 4.

Table 3

Evaluation of statistical significance of the extracted features, before feature selection. The t-test was used to identify which features show significant differences between healthy and MS (ON and None-ON) cases. The p values indicate the test rejection of the null hypothesis at 5% significance level, considering the Bonferroni correction (p value < 0.001).

Features	Layer	p (HC vs. MS-None-ON)	p (HC vs. MS-ON)	p (MS-ON vs. MS-None-ON)	p(HC vs. MS)	Features	Layer	p (HC vs. MS-None-ON)	p (HC vs. MS-ON)	p (MS-ON vs. MS-None-ON)	p(HC vs. MS)
Autocorrelation	2	<0.001	<0.001	0.966	<0.001	Fractal mean	1	<0.001	<0.001	0.832	0.051
Autocorrelation	3	<0.001	0.101	0.092	0.003	Fractal mean	4	<0.001	<0.001	0.799	0.955
Autocorrelation	4	<0.001	0.026	0.222	0.008	Fractal mean	5	<0.001	<0.001	0.862	0.896
Cluster prominence	3	<0.001	0.086	0.033	<0.001	Fractal mean	8	<0.001	<0.001	0.339	0.785
Cluster prominence	4	0.002	0.057	0.525	<0.001	Fractal mean	9	<0.001	<0.001	0.964	0.710
Cluster shade	2	0.038	<0.001	0.209	0.044	Fractal mean	10	<0.001	<0.001	0.514	0.408
Cluster shade	3	<0.001	0.025	0.054	<0.001	Fractal Std.	5	<0.001	<0.001	0.744	0.014
Cluster shade	4	<0.001	0.015	0.298	<0.001	Fractal Std.	8	0.113	0.001	0.276	<0.001
Contrast	1	0.341	<0.001	0.058	0.012	Fractal Std.	9	<0.001	<0.001	0.751	<0.001
Contrast	2	<0.001	<0.001	0.944	<0.001	Fractal Std.	10	<0.001	<0.001	1.000	<0.001
Contrast	3	<0.001	0.005	0.150	0.001	LBP mean	2	0.001	0.003	0.981	<0.001
Contrast	4	<0.001	0.001	0.398	<0.001	LBP mean	3	0.001	<0.001	0.950	0.017
Correlation	1	<0.001	<0.001	0.999	0.927	LBP mean	4	<0.001	<0.001	0.127	0.024
Correlation	2	<0.001	0.004	0.029	0.546	LBP mean	5	<0.001	<0.001	0.353	0.061
Difference entropy	2	<0.001	<0.001	0.929	0.951	LBP mean	6	<0.001	<0.001	0.941	0.335
Difference entropy	3	<0.001	<0.001	0.513	0.677	LBP mean	7	<0.001	<0.001	0.908	<0.001
Difference entropy	4	<0.001	<0.001	0.546	0.263	LBP mean	9	<0.001	<0.001	0.251	0.068
Difference entropy	5	<0.001	0.004	0.413	0.733	LBP Std.	4	<0.001	<0.001	0.409	0.319
Difference variance	2	<0.001	<0.001	0.999	0.001	LBP Std.	5	<0.001	<0.001	0.456	0.645
Difference variance	3	<0.001	0.001	0.169	<0.001	LBP Std.	6	<0.001	0.001	0.905	0.234
Difference variance	4	<0.001	<0.001	0.523	<0.001	LBP Std.	7	<0.001	<0.001	0.205	0.116
Dissimilarity	2	<0.001	<0.001	0.790	0.010	LBP Std.	9	<0.001	<0.001	0.417	0.150
Dissimilarity	3	<0.001	0.001	0.271	0.030	LBP dynamic range	3	0.525	0.002	0.088	<0.001
Dissimilarity	4	<0.001	0.001	0.380	0.013	LBP dynamic range	5	<0.001	0.004	0.567	0.164
Dissimilarity	5	<0.001	0.035	0.274	0.302	LBP dynamic range	6	<0.001	0.006	0.675	<0.001
Energy	2	<0.001	<0.001	0.999	<0.001	LBP kurtosis	2	<0.001	<0.001	0.804	0.968
Energy	3	<0.001	<0.001	0.921	<0.001	LBP kurtosis	3	<0.001	<0.001	0.941	0.048
Energy	4	<0.001	0.001	0.528	0.861	LBP kurtosis	4	0.001	<0.001	0.805	0.001
Energy	5	<0.001	0.019	0.329	0.008	LDP mean	2	0.459	0.659	0.961	<0.001
Entropy	2	<0.001	<0.001	0.896	0.903	LDP mean	4	<0.001	<0.001	0.564	0.067
Entropy	3	<0.001	<0.001	0.446	0.765	LDP mean	5	<0.001	<0.001	0.555	0.035
Entropy	4	<0.001	<0.001	0.453	0.274	LDP mean	6	<0.001	<0.001	0.121	0.799
Entropy	5	<0.001	0.029	0.222	0.871	LDP mean	7	<0.001	<0.001	0.883	0.047
Homogeneity	2	<0.001	<0.001	0.942	0.900	LDP mean	8	<0.001	<0.001	0.808	0.055
Homogeneity	3	<0.001	<0.001	0.681	0.690	LDP mean	9	<0.001	<0.001	0.219	0.037
Homogeneity	4	<0.001	<0.001	0.504	0.201	LDP mean	10	<0.001	<0.001	0.834	0.782
Homogeneity	5	<0.001	0.006	0.396	0.680	LDP skewness	2	<0.001	<0.001	0.813	0.530
IMC1	2	<0.001	<0.001	0.862	0.017	LDP skewness	3	<0.001	<0.001	0.930	0.085
IMC2	9	0.004	<0.001	0.619	0.043	LDP skewness	4	0.004	<0.001	0.633	0.002
Inverse difference moment normalized	1	0.319	<0.001	0.063	0.011	LDP Std.	2	<0.001	<0.001	0.276	0.018
Inverse difference moment normalized	2	<0.001	<0.001	0.920	<0.001	LDP Std.	4	<0.001	<0.001	0.652	0.083
Inverse difference moment normalized	3	<0.001	0.006	0.154	0.001	LDP Std.	5	<0.001	<0.001	0.882	0.025
Inverse difference moment normalized	4	<0.001	0.001	0.381	<0.001	LDP Std.	6	0.005	<0.001	0..093	0.597
Maximum probability	2	<0.001	<0.001	1.000	<0.001	LDP Std.	7	<0.001	<0.001	0.554	0.310
Maximum probability	3	<0.001	<0.001	0.876	<0.001	LDP Std.	8	<0.001	<0.001	0.989	0.870
Maximum probability	4	<0.001	0.068	0.293	0.430	LDP Std.	9	<0.001	<0.001	0.629	.882
Maximum probability	5	<0.001	0.062	0.222	0.647	LDP Std.	10	<0.001	<0.001	0.901	0.792
Sum average	2	<0.001	<0.001	0.814	0.005	LDP dynamic range	1	0.188	0.219	0.999	<0.001
Sum average	3	<0.001	0.013	0.194	0.028	LDP kurtosis	2	<0.001	<0.001	0.998	<0.001
Sum average	4	<0.001	0.006	0.383	0.048	LDP kurtosis	3	<0.001	<0.001	1.000	<0.001
Sum entropy	2	<0.001	<0.001	0.953	0.708	LDP kurtosis	4	<0.001	<0.001	0.801	0.466
Sum entropy	3	<0.001	0.002	0.402	0.938	LDP kurtosis	5	<0.001	0.003	0.781	0.001
Sum entropy	4	<0.001	0.001	0.514	0.384	LOOP mean	2	0.459	0.659	0.961	<0.001
Sum entropy	5	<0.001	0.064	0.182	0.806	LOOP mean	4	<0.001	<0.001	0.564	0.067
Sum of squares	2	<0.001	<0.001	0.528	0.001	LOOP mean	5	<0.001	<0.001	0.555	0.035
Sum of squares	3	<0.001	0.014	0.092	<0.001	LOOP mean	6	<0.001	<0.001	0.121	0.799
Sum of squares	4	<0.001	0.003	0.262	<0.001	LOOP mean	7	<0.001	<0.001	0.883	0.047
Sum variance	2	<0.001	<0.001	0.400	0.002	LOOP mean	8	<0.001	<0.001	0.808	0.055
Sum variance	3	<0.001	0.016	0.087	<0.001	LOOP mean	9	<0.001	<0.001	0.219	0.037
Sum variance	4	<0.001	0.004	0.248	<0.001	LOOP mean	10	<0.001	<0.001	0.834	0.782
LOOP Std.	2	<0.001	<0.001	0.964	0.013	LOOP skewness	2	<0.001	<0.001	0.998	<0.001
LOOP Std.	3	<0.001	<0.001	0.202	0.021	LOOP skewness	3	<0.001	<0.001	1.000	<0.001
LOOP Std.	4	<0.001	<0.001	0.800	0.218	LOOP skewness	4	<0.001	<0.001	0.818	0.641
LOOP Std.	5	<0.001	<0.001	0.830	0.019	LOOP skewness	5	<0.001	0.002	0.769	0.001
LOOP Std.	6	0.014	<0.001	0.091	0.548	LOOP Std.	9	<0.001	<0.001	0.538	0.390
LOOP Std.	7	<0.001	<0.001	0.576	0.425	LOOP Std.	10	<0.001	<0.001	0.786	0.799
LOOP Std.	8	<0.001	<0.001	0.463	0.316

Table 4

Frequency of significant selected features for each retinal layer.

Layers	Frequency
	HC vs. MS	HC vs. MS-ON	HC vs. MS-None-ON	MS-ON vs. MS-None-ON
1	1	4	2	0
2	10	22	22	0
3	10	11	21	0
4	7	15	24	0
5	0	8	18	0
6	1	5	5	0
7	1	6	6	0
8	1	5	5	0
9	1	9	8	0
10	1	5	7	0
Total	33	90	117	0

3.2. Classification Result

According to Tables 3 and 4, 15 common significant features between three groups (HC, MS-ON, and MS-None-ON) are selected as input of each classifier. However, no significant feature is found for the last group (MS-ON vs. MS-None-ON). Then, the classification step is done to discriminate between four target groups including HC vs. MS, HC vs. MS-None-ON, HC vs. MS-ON, and MS-ON vs. MS-None-ON. The accuracy results obtained from our classifiers in different groups are shown in Table 5.

Table 5

The accuracy of texture features, thickness, and combination of texture features and thicknesses.

Methods	Classifiers	HC vs. MS	HC vs. MS-ON	HC vs. MS-None-ON	MS-ON vs. MS-None-ON
Texture features	SVM	85.3	83.6	78.6	64.1
	LDA	72.0	74.6	64.3	48.8
Thicknesses	SVM	84.0	81.8	90.0	89.7
	LDA	64	69.1	73.3	82.1
Texture features & thicknesses	SVM	96.0	87.3	96.4	82
	LDA	100	98.2	96.5	56.4

In addition, to prepare a fair comparison with the previous studies, we also test the performance of the two classification models using thickness features as input. As abovementioned, utilizing texture features for our intended goal is totally novel and previous researches were only relying on thickness as discriminant features. Therefore, the thickness features are calculated as the average value of distance between two consecutive boundaries, which lead to 10 thickness values out of 11 retinal layer boundaries, and this thickness feature vector is fed also to each classifier. In summary, the following set of information are utilized as input of each two classifiers:

15 common selected texture features based on t-test and Bonferroni correction

Thickness features

Combination of I and II

As it can be found in Table 5, in cases I and II, SVM outperforms LAD. In analyzing the effect of texture and thickness features separately, it has to be mentioned that the best accuracy result for groups HC vs. MS and HC vs. MS-ON is found using texture features and SVM classifier. Meanwhile, for groups HC vs. MS-None-ON and MS-ON vs. MS-None-ON, thickness features and SVM classifier obtain the best accuracy.

Furthermore, the impressive point is that in case III and with the combination of texture and thickness features as input of the classifiers, the result improved considerably and also the performance of the LDA classifier is superior to the SVM performance in most of the conditions.

4. Conclusion

There are no specific tests for MS detection. Instead, a diagnosis of MS often relies on ruling out other conditions that might produce similar signs and symptoms, known as a differential diagnosis. Blood tests, spinal tap (lumbar puncture), evoked potential tests, and MRI are the conventional MS diagnosis methods. The first MR images of MS were produced in the early 1980s [23]. In most people with relapsing-remitting MS, the diagnosis is fairly straightforward and based on a pattern of symptoms consistent with the disease and confirmed by brain imaging scans such as MRI; however, MS diagnosis can be more difficult in patients with unusual symptoms or progressive disease.

MRI-based methods have been indeed the most successful techniques to estimate CNS damage up to the present, although it is becoming increasingly clear that due to the ability of direct visualization of retinal axons, OCT has become an extremely sensitive method for imaging neurodegeneration in MS patients. Studies show thickness reduction in retinal layers of MS patients with and without history of ON by OCT image analysis [4, 9, 24–28]. Hence, OCT is suggested as an important tool for monitoring MS and also as a complementary method for MRI-based diagnosis techniques [29-32]. However, as mentioned above, the majority of previous works are on the thickness analysis of retinal layers. Here, by combining the information of thickness and texture of retinal layers, we prepared a more comprehensive analysis of OCT imaging performance in the diagnosis of MS with or without ON.

Indeed, texture analysis is a novel strategy for studying intrinsic changes in retinal layers during neurodegenerative diseases. MS, as one of the famous neurodegenerative disorders, is investigated in this research.

After performing a vast survey on available texture analysis methods, a treasury of powerful features is collected in this paper. As a primary work, this paper shows the ability of such features in discrimination of HC and MS (ON and None-ON) cases. Even with simple classification methods, the texture features are powerful to diagnose MS cases (from HC ones) with accuracy of 85.3% and 72% with SVM and LDA classifiers, respectively.

Another valuable point is that adding information of conventional thickness values to texture features improves the discrimination between most of the target groups including HC vs. MS, HC vs. MS-None-O, and HC vs. MS-ON. It should be noted that the results of the last group (MS-ON vs. MS-None-ON) are generally weaker than other groups due to the lack of significant discriminant texture features for this group.

Furthermore, the findings show that some layers like 2, 3, and 4 carry more texture information useful in separation of HC from MS cases. Such finding can be a start point for further investigation in this area.