Classification of transcranial Doppler signals using their chaotic invariant measures

Classification of transcranial Doppler signals using their chaotic invariant measures

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180 journal homepage: www.intl.elsevierhealth.com/j...

843KB Sizes 0 Downloads 26 Views

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

journal homepage: www.intl.elsevierhealth.com/journals/cmpb

Classification of transcranial Doppler signals using their chaotic invariant measures Ali Ozturk a,∗ , Ahmet Arslan b,1 a b

Selcuk University, Electrical and Electronics Engineering Department, Alaeddin Keykubat Kampusu, Konya, Turkey Selcuk University, Computer Engineering Department, Alaeddin Keykubat Kampusu, Konya, Turkey

a r t i c l e

i n f o

a b s t r a c t

Article history:

In this study, chaos analysis was performed on the transcranial Doppler (TCD) signals

Received 13 June 2006

recorded from the temporal region of the brain of 82 patients as well as of 24 healthy people.

Received in revised form

Two chaotic invariant measures, i.e. the maximum Lyapunov exponent and the correlation

31 January 2007

dimension, were calculated for the TCD signals after applying nonlinearity and stationarity

Accepted 17 February 2007

tests to them. The sonograms obtained via Burg autoregressive (AR) method demonstrated that the chaotic invariant measures represented the unpredictability and complexity levels

Keywords:

of the TCD signals. According to the multiple linear regression analysis, the chaotic invari-

Transcranial Doppler signals

ant measures were found to be highly significant for the regression equation which fitted

Nonlinear analysis

to the data. This result suggested that the chaotic invariant measures could be used for

Chaotic measures

automatically differentiating various cerebrovascular conditions via an appropriate classi-

Classification

fier. For comparison purposes, we investigated several different classification algorithms.

NEFCLASS

The k-nearest neighbour algorithm outperformed all the other classifiers with a classifica-

Decision trees

tion accuracy of 94.44% on the test data. We used the receiver operating characteristic (ROC)

K-nearest neighbour

curves in order to assess the performance of the classifiers. The results suggested that the

Multilayer perceptron

classification systems which use the chaotic invariant measures as input have potential in detecting the blood flow velocity changes due to various brain diseases. © 2007 Elsevier Ireland Ltd. All rights reserved.

1.

Introduction

Transcranial Doppler (TCD) study of the adult intracerebral circulation is generally used to evaluate intracranial stenosis, cerebral arteriovenous malformations, cerebral vasospasm and cerebral hemodynamics [1]. To make intra-brain blood flow measurements requires sophisticated imaging methods. The Doppler principle utilized in medicine is based on transmitting a signal into the body and observing the frequency changes between that transmitted signal and the received signals reflected from the targets. Under these conditions, there



is a shift in the ultrasound frequency given by fd = ft − fr = 2ft v cos

 c

(1)

In the equation above, ft and fr are the transmitted and received ultrasound frequencies, respectively, v the velocity of the target, c the velocity of the sound in the medium, and  is the angle between ultrasound beam and the direction of motion of the target. Diagnoses made with TCD are based on the detection of increased or decreased blood flow velocity, absence of blood

Corresponding author. Tel.: +90 3322232000; fax: +90 3322410520. E-mail addresses: [email protected] (A. Ozturk), [email protected] (A. Arslan). 1 Tel.: +90 3322232000; fax: +90 3322410520. 0169-2607/$ – see front matter © 2007 Elsevier Ireland Ltd. All rights reserved. doi:10.1016/j.cmpb.2007.02.004

172

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

flow or changes in resistivity. TCD examination technique is used in the diagnosis of stenosis and occlusion of arteries by measuring the blood flow velocities in brain arteries. This is an alternative non-invasive diagnosis method to the angiography and magnetic resonance methods. However, the difficulty of properly enabling the expert medical staff to interpret TCD signals prevents their wider and effective usage in the clinics. Chaos theory provides very powerful methods to understand the dynamics of the complex systems like the human brain which depend on many variables. These methods are especially valuable in the analysis of one-dimensional time series data. According to the Takens theorem [2], a system’s behaviour can be modelled by using the time series generated by only one of its variables even if it is a multidimensional system. A TCD signal may be regarded as a series of discrete values varying in time. This kind of series is called a time series. If a time series produced by a nonlinear system is analyzed by linear methods such as power analysis, linear transformations or parametric linear modelling, the critical features of the system remain undetected and most parts of this time series can be incorrectly evaluated as random noise. Chaos theory (nonlinear time series analysis) provides some tools necessary to quantitatively analyze an experimental time series such as a TCD signal. The values such as maximum Lyapunov exponent and correlation dimension are some of the quantitative measures of chaos. For example, if the maximum Lyapunov exponent is found to be positive, then we can certainly say that there is an underlying chaotic behaviour in the system that generates the nonlinear time series. This means that the system has a chaotic attractor embedded in the phase space which includes nearby orbits that diverge from each other in time. According to the Takens theorem [2], if the real dimension of the attractor is DA , then choosing the embedding dimension DE as DE > 2DA prevents the orbits from crossing each other resulting in a non-periodic flow. In chaotic attractors, since the attracting set is very irregular, the correlation dimension is not integer. On the other hand, in order to apply nonlinear time series analysis (chaos theory) methods, it must be proven that the TCD signals are both nonlinear and stationary. In the literature, the features extracted from Doppler signals via spectral analysis methods (e.g. peak systolic velocity, end diastolic velocity, resistive index, pulsality index, etc.) were widely used for automatic medical diagnosis [3–6]. For example in [4], only two features (i.e. resistivity index and pulsatility index) extracted from sonograms were applied successfully to an adaptive neuro-fuzzy system for the analysis of internal carotid arterial Doppler signals. Two different artificial neural networks were compared with each other in classifing the TCD signals using the parameters such as systole, diastole and resistive index obtained from their spectral curve in [6]. Also, there are studies on the analysis of Doppler signals using chaos theory methods [7–10]. Keunen et al. have suggested that there is an underlying nonlinear dynamics in the TCD signals of healthy subjects [8], whereas Visee et al. recognized that the nonlinear phenomena were lost and uniformity appeared in ischemic cerebrovascular territory compared to the noncompromised side in patients with occlusive cerebrovascular disease [9]. Furthermore, in

[11,12] the spectrum of Lyapunov exponents was estimated from local Jacobi matrices and was applied to a multilayer perceptron neural network (MLPNN) in order to detect the variability of the Doppler signals. In our study, chaos analysis was performed on the TCD signals recorded from the temporal region of the brain of 82 patients as well as of 24 healthy people. Among 82 patients, 20 of them had cerebral aneurysm, 10 had brain hemorrhage, 22 had cerebral oedema and the remaining 30 had brain tumor. All the TCD signals were passed from nonlinearity and stationarity tests. The two strongest chaos indicators, i.e. the maximal Lyapunov exponent and the correlation dimension, were calculated for each signal. For these two chaos measures, the intra-class similarity for each subject group could be easily captured just by inspection. To uncover the possible interclass discrimination quantitatively, we applied multiple linear regression. The result of this statistical analysis suggested that an implicit dependency between the chaotic variables and output values existed. Relying on this observation, we concluded that the chaotic invariant measures could be used as input features to train an appropriate classification algorithm for the automatic diagnosis of the mentioned brain diseases. Several different classifiers (NEFCLASS, Decision tree, K-nearest neighbour (kNN) and MLPNN) were evaluated. Receiver operating characteristic (ROC) curves were used to assess the classification performance of the machine learning methods. The resultant graphs indicated that all of the algorithms had excellent correct classification rates. With the highest accuracy, the kNN algorithm outperformed all the other methods. These findings suggest that the chaotic invariant measures extracted from the TCD signals may be applied as input features to an appropriate classifier to obtain a successful classification system for automatic diagnosis of various cerebrovascular diseases.

2.

Materials and method

2.1.

Hardware

The hardware of the system used for this study involves a 2 MHz ultrasound transducer, analog Doppler unit (multi Doppler transducer XX, DWL Gmb, Uberlingen, Germany), analog/digital interface board (Sound Blaster Pro-16), and PIII 600 MHz microprocessor PC with printer. The Doppler unit is also equipped with imaging software that makes it possible to focus the sample volume at a desired location in the temporal region. The analog Doppler unit can work in both continuous and pulse wave modes. The sample volume sizes can be adjusted according to vessel diameter for recording different blood flow rates. The signal obtained from the blood vessel is transferred to a PC via a 16-bit sound card on an analog/digital interface board.

2.2.

Nonlinearity and stationarity analysis

In this paper, the surrogate data sets method was applied in order to detect the possible nonlinearity. The method of iterative amplitude adjusted Fourier transform (IAAFT)

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

discussed in [13] was used to generate surrogate data sets. The Fourier based surrogates depend on the idea of creating constrained realizations. In the constrained-realization approach, the measurable properties of the time series are considered rather than the underlying model equations. The AAFT was originally proposed by Theiler and Prichard [14]. In [13], Schreiber and Schmitz proposed a method which iteratively corrects deviations of the measured data in spectrum and distribution from the goal set. The surrogate is filtered towards the correct Fourier amplitudes and rankordered to the correct distribution. Almost all methods of time series analysis, linear or nonlinear, must assume a kind of stationarity. Nonlinear statistics include higher order correlations, dimensions, Lyapunov exponents and binned probability distributions [15]. The method used in this study is based on the similarity between parts of the time series themselves, rather than the similarity of parameters derived from the times series by local averages. The nonlinear cross-predictions error [16], i.e. the predictability of one segment using another segment as a reference, was applied to detect and analyse the nonstationarity of the TCD time series. This method provides more detailed, time resolved information than just saying that non-stationarity has been found.

In order to estimate the correlation dimension, we have to provide the time delay information between embedding vectors. Fraser and Swinney [20] suggested a method called mutual information to obtain a reasonable delay. The time delay where the mutual information takes the first minimum value is chosen as the optimum delay.

2.4.

Npairs

N N  

(ε − |sj − sk |)

(2)

j=m k
Here, si are m-dimensional delay vectors, Npairs = (N − m + 1)(N − m − w + 1)/2 is the number of pairs of points covered by the sums,  the Heaviside step function and w is the Theiler window [18]. Since successive elements of a time series are not usually independent, Theiler suggested excluding these temporally correlated points from the pair counting operation. This is achieved by simply ignoring all pairs of points whose time indices differ by less than w. On sufficiently small length scales and when the embedding dimension m exceeds the box-dimension of the attractor [19]: C(m, ε) ∝ εD2

  sn+t − sn+1 

(4)

n0=1

If an attractor is chaotic, the neighbouring trajectories that are initially close will diverge exponentially in time and the points on these trajectories will be dynamically uncorrelated. Since these points are located on the attractor they are correlated in phase space, although they appear randomly scattered points. The idea behind the dimension quantifiers is that the weight p(ε) of a typical ε-ball covering part of the invariant set scales has its diameter given by p(ε) = εD , where the value for D depends also on the weight definition. Using the square of the probability pi of finding a point of the set inside the ball, the dimension is called the correlation dimension, which is computed most efficiently by the correlation sum [17]: 1



1 1 ln |ϑn | N N

Correlation dimension

C(m, ε) =

The maximal Lyapunov exponent

Chaos arises from the exponential growth of infinitesimal perturbations, together with global folding mechanisms to guarantee boundedness of the solutions. This exponential instability is characterized by the spectrum of Lyapunov exponents [21]. Lyapunov exponents are invariant under smooth transformations and are thus independent of the measurement function or the embedding procedure. In order for an attractor to be chaotic, the maximal Lyapunov exponent must be greater than zero which means that the neighbouring trajectories of the attractor are locally diverging in at least one direction of the phase space. According to Rosenstein et al. [22], the following S value is computed for different N values:

S(t) =

2.3.

173

(3)

Since the box-dimension is not known a priori, the convergence of the estimated D2 values in m must be checked.

If S exhibits a linear increase with an identical slope, then this slope can be taken as an estimate of the maximal Lyapunov exponent.

2.5.

The classification algorithms

A NEFCLASS system has a three-layer feedforward architecture that is derived from a generic fuzzy perceptron [23]. Neuro-fuzzy classification approaches aim at creating fuzzy classification rules from data. NEFCLASS is able to learn fuzzy rules and fuzzy sets by simple heuristics. The learning algorithm of NEFCLASS has two stages: rule (structure) learning and parameter learning. Rule learning is achieved via a variation of the approach by Wang and Mendel [24]. Each input space corresponding to a feature is partitioned by a given number of initial fuzzy sets. Then, by processing the training data only once, the antecedents for the rules are created by checking which areas of the input space contain data. An evaluation procedure then creates a rule base by assigning appropriate consequents (class labels) to the discovered antecedents and selects only a certain number of rules with good performance [25–27]. Decision trees are used to classify an example starting at the root of the tree and moving through it until a leaf node, which provides the classification of the instance. Most algorithms developed for learning decision trees are variations on a core algorithm that employs a top-down, greedy search through the space of possible decision trees. The algorithm constructs the tree from a set of training cases. The ID3 algorithm designed by Quinlan [28] searches through the attributes of the training instances and extracts the attribute that best seperates the given examples. C4.5 is a software extension of the basic ID3

174

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

algorithm to address some issues not dealt with by ID3 [29,30]. K-nearest neighbour algorithms handle missing values, are robust to outliers, and can be good predictors. They tend to only handle numeric variables which make their usage suitable for this study. K-nearest neighbour classifier, relying on a distance function, is sensitive to noise and to irrelevant features, because such features have the same influence on the classification as do good and highly predictive features. A solution to this is to pre-process the data to weight features so that irrelevant and redundant features have a lower weight. The kNN is simple, but could be very powerful in some applications [31,32]. There have been several practical applications of the kNN method [33]. The MLPNNs are non-parametric techniques for performing a wide variety of detection and estimation tasks [34–36]. They have features such as the ability to learn and generalize, smaller training set requirements, fast operation after training and ease of implementation. The backpropagation algorithm, which is relatively easy to implement, searches the error surface for points with minimal error via the gradient descent approach. This algorithm computes the error at the output nodes of the network and adjusts the weights between layers in a backward fashion to reduce the overall error of the network to a predefined level. The learning behaviour of the MLPNN may be violated due to higher weight values. The learning rate should be lowered to prevent this. On the other hand, a lower learning rate slows down the learning. Thus, a trade-off must be made by taking the above discussion into account. One way to increase the learning rate without leading to oscillation is to modify the backpropagation algorithm by including a momentum term.

2.6.

Statistical analysis

Multiple linear regression attempts to model the relationship between two or more explanatory variables and a response variable by fitting a linear equation to the observed data. The significance level for the F-test shows whether there is a linear relationship between the input variables and the output variable. The coefficient of determination value indicates how well the obtained linear equation fits to the observed data. For the regression equation, the P-values of the coefficients are estimated to find the statistical significance of the input variables. By means of Bland–Altman plots which show the difference between the methods (A − B) against (A + B)/2 the average, it is much easier to assess the magnitude of disagreement (both error and bias), spot outliers, and see whether there is any trend, for example an increase in A − B for high values. Plotting the data in this way is a powerful technique for displaying the results of a method comparison study [37]. An ROC curve represents the trade-offs between sensitivity and specificity. A good diagnostic test is the one that has small false positive and false negative rates across a reasonable range of cut off values. The evidence of a good test is an ROC curve climbing rapidly towards the upper left hand corner of the graph. This means that the 1-false negative rate is high and the false positive rate is low. The same conclusion can be reached by measuring the area under the curve. A larger area

under the ROC curve means a better diagnostic test. The area may change between 0.5 and 1.0 [38].

3.

Experimental results

For this study, the TCD signals belonging to 106 subjects were recorded. According to the examination which has taken their sonograms into consideration, among 82 patients, 20 of them had cerebral aneurysm, 10 had brain hemorrhage, 22 had cerebral oedema and the remaining 30 had brain tumor. The group consisted of 64 males and 42 females with ages ranging from 3 to 71 years and a mean age of 40.5 ± 0.5 years. The group suffering from cerebral aneurysm consisted of 12 males and 8 females with a mean age 59.5 ± 0.5 years (range 55–65), the group suffering from brain hemorrhage consisted of six females and four males with a mean age 27.0 ± 0.5 years (range 21–36), the group suffering from cerebral oedema consisted of 11 females and 11 males with a mean age 25.0 ± 0.5 years (range 3–40), the group suffering from brain tumor consisted of 18 females and 14 males with a mean age 29.5 ± 0.5 years (range 12–41) and the healthy subjects were 9 females and 15 males with a mean age 31.5 ± 0.5 years (range 23–65). The patients with cerebral aneurysm had an abnormal dilation or balloning of their brain artery due to an acquired or congenital weekness in the wall of the vessel. The brain hemorrhage patients were suffered from either intracerebral hemorrhage or subarachnoid hemorrhage (SAH) due to a rupture of a blood vessel within the head. The patients with brain tumor had changes in their cerebrovascular resistance. There is a significant negative correlation between cerebral vascular resistance and cerebral blood flow, where higher values of cerebral vascular resistance were associated with lower blood flow levels, and vice versa. A representative Doppler sonogram for each patient group is given in Fig. 1. The sonograms were obtained using the Burg autoregressive (AR) method. In the sonogram of the patient with cerebral aneurysm, an elevated mean velocity can be seen, whereas the primary characteristic of the sonogram of the cerebral oedema patient is the diminished end-diastolic velocity. However, the chaotic invariant measures provide some quantitative information about the complexity of a signal which can be visually sensed by the sonograms. Higher values of these measures show randomness and unpredictability of the signal, whereas lower values indicate periodicity and uniformity existing in the signal. For example, the correlation dimension of the TCD signal of the healthy subject whose sonogram is seen in Fig. 1(e) is 3.51, which is the highest among all the samples. Likewise, the D2 of the cerebral aneurysm patient whose sonogram is in Fig. 1(a) is 3.42. The obvious random apperance of these two sonograms also suggests the truth of the assumption of positive correlation between the complexity and a high correlation dimension value. However, despite their randomness there are some repeating peaks and valleys on the sonogram in Fig. 1(a) causing a slightly lower D2 . There is a very prominent periodicity in the sonogram of the brain tumor patient in Fig. 1(d). The D2 value found for the TCD signal of this patient is 2.15 which is the lowest value among all the sample subjects. There appears periodicity along with some randomness in the sonograms in Fig. 1(b) and (c), whereas the D2 values

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

175

Fig. 1 – Representative Doppler sonograms obtained for each of the patient groups: (a) cerebral aneurysm, (b) brain hemorrhage, (c) cerebral oedema, (d) brain tumor, and (e) healthy.

of these two sample patients were 2.92 and 2.36, respectively. The mean values of 1 for the TCD signals of the brain tumor patients and cerebral oedema patients were about 1.58 and 1.83, respectively. The decline in 1 indicates the loss of the chaotic behaviour and implies a more prominent periodicity which can also be detected from the sonograms in Fig. 1(c) and (d), respectively. On the other hand, the 1 found for the TCD signal of the healthy subject was 3.63. Higher values of 1 imply randomness and unpredictability in the long term which can also be observed from the sonogram in Fig. 1(e). According to our chaos analysis, it has been shown that all the TCD time series are not only nonlinear, but also stationary to an acceptable level. This meant that we could apply the chaos theory methods in order to estimate the invariant measures like correlation dimension (D2 ) and maximum Lyapunov

exponent (). The 1 was found to be positive whereas the D2 was found to be greater than 2 and fractional for all the TCD signals. A descriptive statistical summary of these estimations is given in the following Table 1. These two chaotic invariant measures were nearly the same for the patients having the same disease. This means that there was an obvious intra-class similarity for the chaotic invariants. However, an evident inter-class discrimination cannot be easily captured by inspection only. To uncover the implicit dependency between the chaotic variables and the output values, we applied multiple linear regression to make a further quantitative analysis on our data as a very first step. According to the significance level for the F-test which was found to be 0.012, there was a linear relationship between the independent variables, i.e. chaotic invariant measures, and

176

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

Table 1 – The statistical summary of the chaotic invariant measures for different subject groups Subject group Cerebral aneurysm Brain hemorrhage Cerebral oedema Brain tumor Healhty

Mean D2 3.58 2.67 2.46 2.12 3.71

Min D2 3.23 2.41 2.23 2.08 3.42

the output. The coefficient of determination value was found to be 0.0829 which meant that the obtained linear regression equation was not good on the data. For the regression equation, the P-values of the two coefficients were found to be 0.003 and 0.2676, respectively. This means that although the correlation dimension coefficient is statistically significant and has a greater effect on the output values than the Lyapunov exponent coefficient, we cannot completely ignore the latter in a possible linear model either.

3.1.

Nonlinearity and stationarity analysis

In order to perform nonlinearity analysis on the TCD time series, the surrogate data were generated with IAAFT method. The root mean squared errors for the original time series and the surrogate time series were produced and sorted optionally in ascending or descending order. If the roots mean squared prediction error corresponding to the original time series was the last or the first row depending on the sort option, then it could be said that the original time series was nonlinear. In order to decide the number of surrogates to be produced, we use a residual probability ˛ of a false rejection, corresponding to a level of significance (1 − ˛) × 100%. Then, we need 1/˛ number of surrogates (including the original time series) for a one-sided test. Thus, for a significance level of 95%, we had to produce at least 19 surrogates. For stationarity analysis, the total time series of 10,000 points was divided into 40 segments each of which involved 250 points. For the TCD time series of a cerebral aneurism patient, the mean cross-prediction error was 6.94, the standard deviation (S.D.) was 0.72, the maximum cross-prediction error was 9.8 and the minimum cross-prediction error was 5.0. As a comparison, the diagonal cross-prediction errors for the Ikeda map data of 20,000 points and Henon map data of 20,000 points were evaluated. Each of these time series was also divided into 40 segments having 500 points. The Ikeda and Henon map time series are known to be chaotic time series representing low-dimensional determinism. For the Ikeda map time series, the mean cross-prediction error was 6.90, the S.D. was 0.31, the maximum cross-prediction error was 8.0 and the minimum cross-prediction error was 6.2. For the Henon map time series, the mean cross-prediction error was 2.22, the S.D. was 0.27, the maximum cross-prediction error was 2.8 and the minimum cross-prediction error was 1.8. These values show that for a time series, being chaotic does not mean low cross-prediction errors between segments, as is the case for the Ikeda map data. On the other hand, for artificial chaotic time series like the Ikeda and Henon maps, the S.D. of prediction errors is lower than the S.D. of prediction errors for the real world time series such as TCD signals. The truth of this assumption is also demonstrated by the

Max D2 3.84 2.92 2.62 2.32 3.92

Mean  3.28 3.12 1.84 1.47 3.68

Min 

Max 

2.96 2.71 1.74 1.42 3.60

3.38 3.26 1.98 1.58 3.83

Figs. 2 and 3. The Bland–Altman plot shows the level of agreement between the cross-prediction errors of the Ikeda and Henon map data. According to this plot, the upper and lower 95% confidence limits of 1.96 S.D. of the differences between the cross-prediction errors are 51.5 and 42.3, respectively. It can also be seen in this plot that the mean of the differences between the cross-prediction errors is 46.9.

3.2.

Correlation dimension estimation

To estimate the correlation dimension, the local slopes of the correlation sums for different embedding dimensions were obtained and written into a file. Then we plotted the contents of the file in x–y coordinates where the ε values correspond to x coordinate and the local slopes of the correlation sums correspond to y coordinate. The y value corresponding to the plateau which the curves for different ε values exposed was the correlation dimension of the chaotic attractor. In Fig. 4, the correlation dimension estimation for the TCD time series of a patient with cerebral aneurism is given as an example. For all the nonlinear time series investigated in this study, the time delay  was found to be 1 according to the mutual information method. In other words, the mutual information takes the first minimum value at time delay  = 1 for all of the TCD time series.

3.3.

The maximal Lyapunov exponent estimation

The program which computes the maximal Lyapunov exponent for the TCD time series implements the method of Rosenstein et al. [22] which is known to be reliable in small data sets. The program produces a file with two columns, first of which shows the iteration number and the second one shows the logarithm of the stretching factor. Maximal Lyapunov exponent was obtained by plotting these values using linear scales for both of the axes. The first slope of the resulting curve gives the maximal Lyapunov exponent. In Fig. 5, it is shown how this estimation was made for the TCD time series of a patient with cerebral aneurism.

3.4.

The classification results

The data set of 106 subjects was divided into two subsets for training and testing, containing 52 and 54 subjects, respectively. We intuitively assumed and empirically observed that this kind of division of the data set would provide generalization while avoiding overfitting for the machine learning algorithms used in this study. The reason of using a major part of the data set for testing was to observe the performance of the classifiers in case of many different unseen subjects.

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

177

Fig. 2 – Comparison of the diagonal cross-prediction errors: (a) Ikeda map, (b) Henon map, and (c) the TCD signal of a patient with cerebral aneurism.

The NEFCLASS system was trained for 1000 epochs with automatic pruning. The learning rate was chosen as 0.01. The overall accuracy for the test set after pruning was found to be 88.88% with six misclassifications out of 54 test subjects, whereas the accuracy before pruning was found to be 87% with seven misclassifications. The paired samples t-test statistics were made on the error values both before and after prun-

Fig. 5 – Maximum Lyapunov exponent estimation for the TCD signal of a patient with cerebral aneurism.

Fig. 3 – Bland–Altman plot for diagonal cross-prediction errors of Ikeda and Henon maps data.

Fig. 6 – Bland–Altman plot for before and after pruning error values of NEFCLASS.

Fig. 4 – Correlation dimension estimation for the TCD signal of a patient with cerebral aneurism.

ing. A P-value of 0.1351 was found as the result of this test. This means that we cannot reject the null hypothesis that the average of the differences between the paired error values is zero. A Bland–Altman plot was also given in Fig. 6 showing the level of agreement between the error values before and after pruning.

178

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

Among 15 cerebral oedema patients, three subjects were incorrectly classified as having cerebral aneurism whereas among 10 cerebral aneurism patients, three subjects were incorrectly classified as healthy. The implemented decision tree had five leaf nodes each of which corresponds to a different output class. For the test set, the overall classification accuracy of the decision tree was found to be 92.60% with four misclassifications out of 54 subjects. The mean absolute error and the root mean squared error were found to be 0.0296 and 0.1721, respectively. The misclassified subjects were the ones with brain tumor. Only three subjects among 17 brain tumor subjects were incorrectly classified as having brain aneurism and one of them was incorrectly classified as healthy subject. This instance-based algorithm was implemented using the four nearest neighbours around a central point. After the training stage, the overall classification accuracy was found to be 94.44% with three misclassifications out of 54 subjects. The mean absolute error and the root mean squared error were found to be 0.031 and 0.1369, respectively. The misclassified

subjects were the ones with brain tumor. Only two subjects among 17 brain tumor subjects were incorrectly classified as having brain aneurism and one of them was incorrectly classified as healthy. The MLPNN which was used in this study contained three layers. The input layer had two nodes for the chaotic invariant measures. The output layer had five nodes corresponding to the classes to which the subjects belonged whereas there were three nodes in the hidden layer. The MLPNN was trained for 1000 epochs with a learning rate of 0.3 and a momentum value of 0.2. Initially small random values were assigned to the weights. After the training, the mean absolute error and the root mean squared error were 0.017 and 0.027, respectively. The overall classification accuracy of the MLPNN for the test set was found to be 92.60% with four misclassifications out of 54 subjects. After the application of the test set, the mean absolute error and the root mean squared error of the MLPNN were found to be 0.042 and 0.1445, respectively. The misclassified subjects were the ones with brain tumor. Only

Fig. 7 – The ROC curves indicating the performance of the classification methods (a) NEFCLASS, (b) Decision Tree, (c) K-nearest neighbour, and (d) MLPNN.

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

three subjects among 17 brain tumor subjects were incorrectly classified as having brain aneurism and one of them was misclassified as a healthy subject. According to the classification results of each method, the ROC curves which describe the performance of the classification systems were obtained. These plots provided a view of the whole spectrum of sensitivities and specificities, because all possible sensitivity and specificity pairs for a particular method were graphed. The ROC curves given in Fig. 7 indicate that the classification performances of the machine learning algorithms are excellent. The areas under the ROC curves are (a) 0.944, (b) 0.929, (c) 0.945, and (d) 0.929, respectively.

4.

Discussions

This paper demonstrated that the two chaotic invariant measures, i.e. maximum Lyapunov exponent (1 ) and correlation dimension (D2 ), represented the unpredictibility and complexity levels of the TCD signals, respectively. Since the Doppler signals are conventionally interpreted by analyzing their spectral content, we investigated the sonograms obtained by the Burg AR method in order to verify this finding. We noticed that the more periodic and uniform the sonogram is due to a blood flow velocity change, the lower the 1 is, whereas the more complex the sonogram is as in healthy subjects, the higher the D2 is. These findings confirm the previous studies in the literature which found a positive correlation between the complexity and unpredictibility levels of the biomedical signals and the chaotic invariant measures. Moreover, according to the result of the multiple linear regression analysis, both of the chaotic measures were found to be highly significant and had a direct effect on the dependent variable of the regression equation. Based on this observation, we concluded that these two chaotic invariant measures could be used for automatic diagnosis of blood flow changes due to brain diseases by means of an appropriate classification algorithm. The whole data set was divided into train and test subsets. For comparison purposes, NEFCLASS, Decision Tree, kNN and MLPNN algorithms were trained using the training set. The classification accuracy of 94.44% of the kNN algorithm was the best, while the accuracy of 88.88% of the NEFCLASS was the worst among them. The performance of the C4.5 decision tree algorithm and the MLPNN had the same accuracy level of 92.59%. The ROC curves also suggested that the classification performance of the classifiers was very promising, indicating the success of the feature extraction method. In summary, we used two chaotic invariant measures (i.e. maximum Lyapunov exponent and correlation dimension) which were extrated from TCD signals as the input features to implement an automatic diagnosis system. The results suggest that our diagnostic system successfully detects the blood flow velocity changes in brain arteries caused by various cerebrovascular diseases. The classification result of our proposed method can be used as supporting clinical information in addition to the sonograms and spectral waveform analysis parameters in order to improve the overall accuracy of the medical diagnosis.

179

references

[1] D.H. Evans, W.N. McDicken, R. Skidmore, J.P. Woodcock, Doppler Ultrasound: Physics. Instrumentation and Clinical Applications, Wiley, Chichester, 1989. [2] F. Takens, Detecting Strange Attractors in Turbulence. Lecture Notes in Math, vol. 898, Springer, New York, 1981, pp. 226–381. [3] U. Ergun, S. Serhatlioglu, F. Hardalac, I. Guler, Classification of carotid artery stenosis of patients with diabetes by neural network and logistic regression, Comp. Biol. Med. 24 (2004) 389–405. [4] E.D. Ubeyli, I. Guler, Adaptive neuro-fuzzy inference systems for analysis of internal carotid arterial Doppler signals, Comp. Biol. Med. 35 (2005) 687–702. [5] I. Guler, F. Hardalac, N. Barisci, Application of FFT analyzed cardiac Doppler signals to fuzzy algorithm, Comp. Biol. Med. 32 (2002) 435–444. [6] S. Serhatlioglu, F. Hardalac, I. Guler, Classification of transcranial Doppler signals using artificial neural network, J. Med. Syst. 27 (2003) 205–214. [7] J.H.R. Vliegen, C.J. Stam, S.A.R. Rombouts, R.W.M. Keunen, Rejection of the ‘filtered noise’ hypothesis to explain the variability of transcranial Doppler signals: a comparison of original TCD data with Gaussian-scaled phase randomized surrogate data sets, Neurol. Res. 18 (1996) 19–24. [8] R.W. Keunen, J.H. Vliegen, C.J. Stam, D.L. Tavy, Nonlinear transcranial Doppler analysis demonstrates age-related changes of cerebral hemodynamics, Ultrasound Med. Biol. 22 (1996) 383–390. [9] H.F. Visee, R.W. Keunen, H.C. Pijlman, J.H. Vliegen, D.L. Tavy, K.J. Stam, C.A. Giller, The physiological and clinical significance of nonlinear TCD waveform analysis in occlusive cerebrovascular disease, Neurol. Res. 17 (1995) 384–388. [10] R.W. Keunen, H.C. Pijlman, H.F. Visee, J.H. Vliegen, D.L. Tavy, K.J. Stam, Dynamical chaos determines the variability of transcranial Doppler signals, Neurol. Res. 16 (1994) 353–358. [11] E.D. Ubeyli, I. Guler, Determining variability of ophthalmic arterial Doppler signals using Lyapunov exponents, Comp. Biol. Med. 35 (2005) 405–420. [12] I. Guler, E.D. Ubeyli, Detecting variability of internal carotid arterial Doppler signals by Lyapunov exponents, Med. Eng. Phys. 26 (2004) 763–771. [13] T. Schreiber, A. Schmitz, Improved surrogate data for nonlinearity tests, Phys. Rev. Lett. 77 (1996) 635–638. [14] J. Theiler, D. Prichard, Constrained-realization Monte-Carlo method for hypothesis testing, Phys. D. 94 (1996) 221–235. [15] H. Isliker, J. Kurths, A test for stationarity: finding parts in time series apt for correlation dimension estimates, Int. J. Bifurcation Chaos 3 (1993) 1573–1579. [16] T. Schreiber, Detecting and analyzing non-stationarity in a time series with nonlinear cross-predictions, Phys. Rev. Lett. 78 (1997) 843–846. [17] P. Grassberger, I. Procaccia, Characterization of strange attractors, Phys. Rev. Lett. 50 (1983) 346–349. [18] J. Theiler, Estimating fractal dimension, J. Opt. Soc. Am. A7 (1990) 1055–1073. [19] T. Sauer, J. Yorke, How many delay coordinates do you need? Int., J. Bifurcation Chaos 3 (1993) 737–744. [20] A.M. Fraser, H.L. Swinney, Independent coordinates for strange attractors from mutual information, Phys. Rev. A 33 (22) (1986) 1134–1140. [21] J.P. Eckmann, D. Ruelle, Ergodic theory of chaos and strange attractors, Rev. Mod. Phys. 57 (1985) 617–656.

180

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 8 6 ( 2 0 0 7 ) 171–180

[22] M.T. Rosenstein, J.J. Collins, C.J. De Luca, A practical method for calculating largest Lyapunov exponents from small data sets, Phys. D 65 (1993) 117–134. [23] D. Nauck, R. Kruse, Designing Neuro-Fuzzy Systems through Backpropagation. Fuzzy Modelling: Paradigms and Practice, Kluwer, 1996. [24] L.X. Wang, J.M. Mendel, Generating fuzzy rules by learning from examples, IEEE Transactions on System, Man Cybern. 22 (6) (1992) 1414–1427. [25] D. Nauck, U. Nauck, R. Kruse, NEFCLASS for Java – new learning algorithms, in: Proceedings of the IEEE 18th International Conference of the North American Fuzzy Information Processing Society (NAFIPS99), New York, NY, 1999, pp. 472–476. [26] D. Nauck, R. Kruse, A neuro-fuzzy method to learn fuzzy classification rules from data, Fuzzy Sets Syst. 89 (1997) 277–288. [27] D. Nauck, U. Nauck, R. Kruse, Generating classification rules with the neuro-fuzzy system NEFCLASS, in: Proceedings of the Biennial Conference of the North American Fuzzy Info. Proceedings Society NAFIPS’96, Berkeley, CA, June, 1996. [28] J.R. Quinlan, Induction of decision trees, Machine Learn. 1 (1986) 81–106.

[29] J.R. Quinlan, C4. 5: Programs for Machine Learning, Calif.: Morgan Kaufmann, San Mateo, 1993. [30] J. Han, M. Kamber, Data Mining Concepts and Techniques, The Morgan Kaufmann Series in Data Management Systems, Jim Gray, Series Editor Morgan Kaufmann Publishers, 2000. [31] T. Mitchell, Machine Learning, McGraw Hill, 1997. [32] Z. Zhang, Association Rule Mining, Springer, 2002. [33] A.W. Moore, Fast, Robust Adaptive Control by Learning only Forward Models. Advances in Neural Information Processing Systems, Morgan Kaufmann, 1992. [34] S. Haykin, Neural Networks: A Comprehensive Foundation, Macmillan, New York, 1994. [35] I.A. Basheer, M. Hajmeer, Artificial Neural Networks: Fundamentals, Computing, Design and Application, J. Microbiol. Meth. 43 (2000) 3–31. [36] B.B. Chaudhuri, U. Bhattacharya, Efficient training and improved performance of multilayer perceptron in pattern classification, Neurocomputing 34 (2000) 11–27. [37] J.M. Bland, D.G. Altman, Statistical methods for assessing agreement between two methods of clinical measurement, Lancet 1 (1986) 307–310. [38] M.H. Zweig, G. Campbell, Receiver-operating characterictic (ROC) plots: a fundamental evaluation tool in clinical medicine, Clin. Chem. 39 (4) (1993) 561–577.