Application of fractal theory in analysis of human electroencephalographic signals

Application of fractal theory in analysis of human electroencephalographic signals

Computers in Biology and Medicine 38 (2008) 372 – 378 www.intl.elsevierhealth.com/journals/cobm Application of fractal theory in analysis of human el...

236KB Sizes 1 Downloads 93 Views

Computers in Biology and Medicine 38 (2008) 372 – 378 www.intl.elsevierhealth.com/journals/cobm

Application of fractal theory in analysis of human electroencephalographic signals P. Paramanathan, R. Uthayakumar ∗ Department of Mathematics, Gandhigram Rural University, Gandhigram 624 302, Tamil Nadu, India Received 5 December 2006; accepted 11 December 2007

Abstract In medical discipline, complexity measure is focused on the analysis of nonlinear patterns in processing waveform signals. The complexity measure of such waveform signals is well performed by fractal dimension technique, which is an index for measuring the complexity of an object. Its applications are found in diverse fields like medical, image and signal processing. Several algorithms have been suggested to compute the fractal dimension of waveforms. We have evaluated the performance of the two famous algorithms namely Higuchi and Katz. They contain some problems of determining the initial and final length of scaling factors and their performance with electroencephalogram (EEG) signals did not give better results. In this paper, fractal dimension is proposed as an effective tool for analyzing and measuring the complexity of nonlinear human EEG signals. We have developed an algorithm based on size measure relationship (SMR) method. The SMR algorithm can be used to detect the brain disorders and it locates the affected brain portions by analyzing the behavior of signals. The efficiency of the algorithm to locate the critical brain sites (recurrent seizure portion) is compared to other fractal dimension algorithms. The K-means clustering algorithm is used for grouping of electrode positions. 䉷 2007 Elsevier Ltd. All rights reserved. Keywords: Fractal dimension (FD); K-means clustering; Critical brain sites; Size measure relationship (SMR) method; Electroencephalogram (EEG)

1. Introduction During the last two decades, we have become fascinated with the realization that some nonlinear dynamical systems have very complex behaviors. Electrical signals of electroencephalogram (EEG) from the brains of seizure patients also show complexity in the behavior. EEG monitoring systems have become important clinical tools for evaluation and treatment of seizure. Improving seizure prediction is the hot topic in epilepsy research. Detection of critical cortical site (recurrent seizure portion) is used to give treatment for the patient. The major physiological brain states that have distinct EEG and behavior correlates are concerned with the degree of vigilance. They range from extreme alertness to very deep sleep or unconsciousness (the latter may be considered as an insult-induced defensive state). Some totally different states could present apparently similar ∗ Corresponding author. Tel.: +91 451 2452371; fax: +91 451 2453071.

E-mail addresses: [email protected] (P. Paramanathan), [email protected] (R. Uthayakumar). 0010-4825/$ - see front matter 䉷 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.compbiomed.2007.12.004

EEG patterns. The sites affected in the brain area which are related to other areas may also be affected [1]. The variation of activity in the brain differs from area to area. There is no similarity in the random collection of cortical areas. Calculating the amount of chaos in the EEG signal using the methods of nonlinear dynamics consumes more time. The EEG signal is the gateway to measure the complexity of the nonlinear behavior of the brain. Human brain is considered as a dynamical system. Estimating the complexity in the signal gives the state of the dynamical system [2]. We found out the dimension for various signals obtained from different positions of the electrodes, which are needed to specify the state of a dynamical system [3]. One extreme in the spectrum of deterministic behaviors for nonlinear systems is chaos, in which the system shows no evidence of settling down to any stable state and shows extreme (exponential) sensitivity to small perturbations. The techniques from nonlinear dynamics such as Lyapunov exponents and correlation dimension can describe the complexity of an EEG data. But both work with small samples and they also need large amount of computations. A hallmark of chaotic

P. Paramanathan, R. Uthayakumar / Computers in Biology and Medicine 38 (2008) 372 – 378

2. Methods 2.1. Data acquisition 2.1.1. Synthetic signal Synthetic data were produced using the deterministic Weierstrass cosine function, given as follows: WH (t) =

M 

−iH cos(2i t),

0 < H < 1,

(1)

i=0

where  > 1, H—the Hausdorff dimension and we fixed  = 5 and M = 26. The fractal dimension of this signal is given by D equals 2-H . A set of 100 sequences, each with different FD, were generated using Weierstrass function. The FD of an experimental signal was computed using sliding window and normal approach. A sliding window of 1.25 s was used to promote stationarity in each segment analyzed, considering that our EEGs were sampled at 250 Hz. Fig. 1(a) and (b) shows the

2 FD=1.5

1.5 1 Amplitude

0.5 0 -0.5 -1 -1. 5 -2 0

100 200 300 400 500 600 700 800 900 1000 Number of Points

1.5 FD=1.2 1 0.5 Amplitude

systems is that the dimension seems fractional or fractal, that is 1.26 rather than 1 or 2. A global value that is relatively simple to compute is the fractal dimension (FD); it refers to a noninteger or fractional dimension of a geometrical object [2]. The FD can give an indication of the dimensionality and complexity of the system. Since actual living biological systems are not stable and the system complexity varies with time, one can distinguish between different states of the system by the FD; it can also determine whether a particular system is more complex than other systems [2,4]. The algorithms are developed by Higuchi [4] and Katz [6], for estimating the FD of the waveform signals, but their computational requirements are expensive. These algorithms contain more number of steps with lot of comparisons, so the time complexity is very high [4–6]. The Higuchi and Katz algorithms work only when the waveform signal contains less variance, so the algorithms do not give better results for EEG signals. The main contribution of this paper is presenting an efficient algorithm for calculating FD of EEG signal with less time and space complexity. Along with other nonlinear methods, the fractal dimension has been used in the analysis of biomedical signals exhibiting non-stationary and transient characteristics such as the EEG. The short-term analysis of EEG signals and their dimension analysis will provide the measure of rate of the complexity. We compared the results of our algorithm with the results of Higuchi and Katz algorithms. From our results the complexity decrease prior to the seizure is confirmed with fractal dimension values. It discriminates the ictal (seizure period) and pre-ictal (before seizure) portions in the brain effectively. The clustering technique (K-means algorithm) is used to locate the critical cortical site in the brain. The plan of the paper is as follows. In Section 2, we describe the method of data acquisition, the Higuchi, Katz and our proposed methods of computing fractal dimension and clustering of electrode sites. Section 3 illustrates the comparison of results obtained from three methods. Section 4 summarizes our conclusions.

373

0 -0.5 -1 -1.5 0

100 200 300 400 500 600 700 800 900 1000 Number of Points

Fig. 1. (a) and (b) Synthetic signals with known FD.

two sequences generated from Weierstrass cosine function with known fractal dimension values. 2.1.2. Spontaneous ongoing signal In general, one or more tape-recorded channels of suitably filtered and amplified EEG signal, extending over the period of interest, are digitized at a sampling frequency m to form the database. The T second long and L (=T × m) sample long digitized signal, F (n),

n = 1, . . . , (M − 1) × D + N ,

(2)

is split into t second long overlapping segments or epochs which are arranged as the columns or pattern vectors of an N × M matrix F, ⎡

F (1) ⎢ F (2) ⎢ F =⎢ · ⎣ · F (N )

⎤ F (D + 1) · · · F ((M − 1) × D + 1) F (D + 2) · · · F ((M − 1) × D + 2) ⎥ ⎥ · ··· · ⎥ , (3) ⎦ · ··· · F (D + N ) · · · F ((M − 1) × D + N )

374

P. Paramanathan, R. Uthayakumar / Computers in Biology and Medicine 38 (2008) 372 – 378

where N (=t·m) is the dimension of each pattern vector, M = ((L + D − N )/D) is the number of patterns, and D is the delay between patterns. In our calculation, we choose N = 2500 points, the length of the signal is estimated by consecutive differences of samples. Out of the EEG data observed from the seizure patient with 32 electrode channels, we have considered the most seizure exhibited channels of 14 only [7].

respect to the first point. The FD compares the actual number of units that compose a curve with the minimum number of units required to reproduce a pattern of the same spatial extent. FDs computed in this fashion depend upon the measurement units used. If the units are different, then so are the FDs. Katz’s approach solves this problem by creating a general unit or yardstick: the average step or average distance between successive points. Normalizing distances in the above equation by this average results in

2.2. Higuchi and Katz methods Consider x(1), x(2), . . . , x(N ) be the time sequence to be k as x k ={x(m), x(m+ analyzed. Construct k new time series xm m k), . . . , x(m + (N − m)/kk)} for m = 1, 2, . . . ., k, where m indicates the initial time value, k indicates the discrete time interval between points (delay) and a means integer part of k constructed, the a. For each of the curves or time series xm average length Lm (k) is computed as Lm (k) (n−m)/k |x(m + ik) − x(m + (i − 1)k)|(n − 1) = i=1 , (4) (N − m)/k where N is the total length of the data sequence and (N − 1) /(N − m)/kk is a normalization factor. An average length is computed for all time series having the same delay (or scale) k, as the mean of the k lengths Lm (k) for m = 1, 2, . . . , k. The procedure is repeated for each k ranging from 1 to kmax , yielding a sum of average lengths L(k) for each k as indicated in L(k) =

k 

Lm (k).

(6)

where L is the total length of the curve or sum of distances between successive points and d is the diameter estimated as the distance between the first point of the sequence and the point of the sequence that provides the farthest distance. Mathematically, d can be expressed as d = max(distance(1, i)).

(8)

Defining n as the number of steps in the curve, then n = L/a, and (5) can be written as D=

log (n)

10 . d + log10 (n) log10 L

(9)

Expression (9) summarizes Katz’s approach to calculate the FD of a waveform [6]. In Katz method, the algorithm got computational difficulty to estimate the chaotic behavior of the waveforms. It only depends on the value of d, but the irregular waveforms contain complexity upto the end of the sequence. So the d value cannot estimate the proper distance between the points. The Katz method suffers from the above problem only. The problems in the algorithms of Higuchi and Katz can be resolved by the Size measure relationship (SMR) algorithm which is presented in the following section. 2.3. Size measure relationship method

The total average length for scale k, L(k) is proportional to k −D where D is the FD by Higuchi’s method. In the curve of ln(L(k)) versus ln(1/k), the slope of the least squares linear best fit is the estimate of the fractal dimension [6]. In Higuchi’s algorithm there is a need to choose the value for Kmin and Kmax , that is, initial and final step of scaling factor values. But the criterion for the selection of these values is not presented, so it gives less effective results. Higuchi fixed the minimum value of K = 2 and maximum value of K = 8 for his algorithm. The FD of a curve can be defined as log10 (L) , log10 (d)

log10 (L/a) . log10 (d/a)

(5)

m=1

D=

D=

(7)

Considering the distance between each point of the sequence and the first, point i is the one that maximizes the distance with

Many other techniques of digital signal processing, the computation of fractal dimension can be undertaken in ‘real place’ (processing the data directly) [3]. Here we propose size measure relationship method based on recursive length of the signal using different measuring scales. The proposed method is very efficient and accurate with respect to time complexity, space complexity, locating seizure portion and calculating the chaos in the dynamical system. The results of this method were compared and realized that it is very sensitive comparing with all the methods so far available in the literature. This method is used to analyze the pathophysiological changes in the brain activity and postural control. In the implementation of Higuchi’s algorithm, it is necessary to determine the parameter kmax (maximum value of interval k), but in our method there is no need to select the kmax value. The Katz’s algorithm works with minimum window length and the signal contains minimum variance. The signal noise and the window length affect the dynamic range of Katz’s algorithm results. This causes new algorithm which clear the above problems. The FD of a signal F can be defined in our method as D=

log10 (L) , log10 (R)

(10)

P. Paramanathan, R. Uthayakumar / Computers in Biology and Medicine 38 (2008) 372 – 378

where L is the total length of the curve or sum of distances between successive points; L can be expressed as L = sum(abs(diff(Finor )))

(11)

and R is the range estimated as the distance between maximum value point and the minimum value point. Mathematically R can be expressed as R = distance (max(Finor ) − min(Finor )).

(12)

Sample differences depend on signal amplitude; a preprocessing step comprising signal amplitude normalization is necessary for signal sample Fi , i = 1, 2, . . . , N and N is the total number of points in the signal sample. Finor =

Fi

, 1 N | |F N j =1 j 1 N where Fi = Fi − Fj , j =1 N

(13)

i = 1, 2, . . . , N.

After normalization of amplitudes we get Finor ; it is used to calculate the fractal dimension of the signal. We can easily observe the changes in the EEG ongoing signal by observing fractal dimension in the peak portions. The perturbations in the brain make extreme sensitivity in the range of dimension [8]. FD allows us to measure the degree of complexity by evaluating how fast our measurements increase or decrease as our scale becomes larger or smaller [3]. Here we have used new algorithm to analyze the EEG signals. It is very simple with less time complexity when compared to any other FD algorithms. 2.3.1. Procedure for computing fractal dimension using SMR method Consider a given signal Fi , i = 1, 2, . . . , N, the procedure for calculating its FD would comprise the following steps: F

1. Normalize amplitudes using Finor = 1 N i  , where Fi = j =1 |Fj | N Fi − N1 N F , i = 1, 2, . . . , N. j =1 j 2. Compute the total length of the signal using L = sum(abs(diff(Finor ))). 3. Compute the range in Finor usingR = distance(max(Finor ) − min(Finor )). log (L) 4. Apply formula D = log10 (R) . 10

The procedure for SMR algorithm is implemented in MATLAB. 2.4. Clustering—selection of critical brain sites Cluster analysis has been applied in many fields. Here we have used cluster analysis to find groups of electrode sites with similar and minimal distance (or nearest neighbour) of its fractal dimensions. There are large numbers of clustering algorithms used by the researchers [9]. Among them, we have used the effective clustering algorithm known as ‘K-means algorithm’

375

for grouping the electrode sites [10]. Besides the data, input to the algorithm consists of K, the number of clusters to be constructed. The K-means algorithm consists of the following steps. 1. Begin with K clusters, each consisting of one of the first k samples. 2. For each sample, find the centroid nearest to it. Put the sample in the cluster identified with this nearest centroid. 3. If no samples changing clusters, stop. 4. Recompute the centroids of altered clusters and go to step 2. Once it is obtained the clusters can be used to classify patterns into distinct classes. Here the signals from each electrode are sampled at 250 Hz. The signals observed form electrode strips are placed over the left orbitofrontal (LOF), right orbitofrontal (ROF), left subtemporal (LST) and right subtemporal cortex (RST). Depth electrodes are placed in the left temporal depth (LTD) and right temporal depth (RTD) to record hippocampus EEG activity [11]. We have selected positions of 14 electrodes among the above locations. The EEG data in the format of 14 columns (electrode positions) and 10 rows (FD values of segments signals from electrode positions) are used for clustering. The K-means algorithm makes the clusters depending on the dimension values. Similar behavior electrodes are in the same group. 3. Results The new approach of estimating the FD of a signal yields good results with various window lengths. The reliability of the algorithm was tested with synthetic signal ranged from 1.001 to 1.099 using Weierstrass functions with known FD [9]. Fig. 2(a) shows the FD values obtained by the Higuchi and Katz algorithms and SMR algorithm against the known FDs of the synthetic data and FD values of EEG signals estimated by our method. Note that perfect reproduction of the known FDs should yield slope of a straight line equal to 1. Higuchi’s algorithm and our method provide the most accurate estimate of the FD. Katz’s method is less linear. In Katz’s method the window length affects the dynamic range of the FD yielding a dynamic range between 1 and 1.2 for window length greater than 750 and lower than 250 points. In the proposed method the window length does not affect the dynamic range of the estimated FD yielding a dynamic range between 1 and 2. It works with synthetic signal data for 8000 points and gives useful results. The FD results obtained with experimental EEG data reveal that our method is the most accurate of the three. Higuchi’s method yielded accurate result for synthetic signal and takes more time to calculate FD of the EEG signal. Katz’s method yielded the most consistent results regarding discrimination between states of the brain, but the dynamic ranges of the results are affected by window length and noise. The FD estimated with 14 channels for each method when the synthetic signal is contaminated with white noise, yielding a signal to noise ratio (SNR) of 10 db [9]. Fig. 2(b)–(d) presents the EEG signal of one epileptic patient analyzed by Katz, Higuchi and SMR

376

P. Paramanathan, R. Uthayakumar / Computers in Biology and Medicine 38 (2008) 372 – 378

FD of EEG Signal by Katz's Method FD=1.6479

20 0 -20 0

FD Estimated by three methods

4 3.5 3

Katz method

Higuchi's method

1.5 1 0.5 1

1.1

1.2

1.3

1.4

1.5

1.6

1.7

1.8

1.9

2

Theoretical FD of Synthetic Signal

0

500

1000

1500

0

500

1000

1500

0

500

1000

1500

0

500

1000

1500

20 0 -20

20 0 -20

500

1000

1500

2000

500

1000

1500

2000

500

1000

1500

2000

2500

2500

2500

FD = 1.5421

0

500

1000

1500

2000

2500

FD = 1.6125

0

500

1000

1500

20 0 -20

2000

500

1000

1500

2000

50 0 -50 20 0 -20 20 0 -20 20 0 -20

2500

1000

1500

2000

2500

2000 2500 FD=1.5645

2000

2500

2000

2500

FD = 1.5360

0

500

1000

1500

0

2500

500

1000

1500

2000

2500

FD = 1.5682

0

500

1000

1500

2000

2500

FD = 1.5707

0

500

1000

1500

2000

2500

FD = 1.5100

0

2500

2000 FD = 1.4907

500

1000

1500

20 0 -20

FD = 1.6580

0

2500

FD of EEG signals by SMR method

FD=1.2157

0

500

20 0 -20

FD =1.1156

0

2000 FD=1.6332

FD=1.5765 0

FD = 1.0016

0

2500

FD=1.6315

FD of EEG Signal by Higuchi's Method 20 0 -20

2000 FD=1.5399

20 0 -20

SMR method

20 0 -20

1500

20 0 -20

2

20 0 -20

1000

20 0 -20

2.5

50 0 -50

500

50 0 -50

2000

2500

FD = 1.4779

0

500

1000

1500

2000

2500

Fig. 2. (a) Performance of three methods for synthetic signal and (b)–(d) performance of Katz, Higuchi and SMR algorithms with EEG signal.

algorithms. The advantages of our method are greater speed and eliminating the need to choose the value for Kmax . The runtime and steps for the proposed method are less than that of other methods; window length does not affect and it works with various window lengths. The new method gives the distinguishability in results between different time points. A comparison of computational burden between the methods for different window lengths is presented in Table 1. The SMR algorithm needs less time complexity compared to Higuchi and Katz algorithms. The new algorithm separates the ictal (seizure period) and pre-ictal (before seizure period) portions of the EEG signal

Table 1 Comparison of runtime Window length

250 500 1000 2000 4000 8000

Runtime SMR method

Higuchi’s

Katz’s

0.0320 0.0470 0.0780 0.1410 0.3280 0.8750

0.1720 0.3440 0.9530 3.0160 12.9310 52.9210

0.0420 0.0680 0.950 0.1950 0.4250 1.0120

P. Paramanathan, R. Uthayakumar / Computers in Biology and Medicine 38 (2008) 372 – 378

377

Table 2 FD values of EEG data segments Electrodes

RST4 RTD2 RTD6 LST4 LOF2 RST2 LST1 LST2 LOF3 ROF2 LTD1 LTD2 LTD3 RST1

EEG data segments (N = 2500 points, sampled at 250 Hz) 1

2

3

4

5

6

7

8

9

10

1.4252 1.4692 1.4639 1.4614 1.4460 1.4303 1.4565 1.4097 1.3208 1.4199 1.4311 1.4346 1.3746 1.2459

1.4565 1.3787 1.4413 1.4583 1.4732 1.3899 1.5113 1.4899 1.3990 1.4850 1.4171 1.5420 1.5881 1.5080

1.3764 1.4339 1.3698 1.4222 1.4267 1.3955 1.5143 1.4754 1.4025 1.3100 1.3405 1.4716 1.4503 1.4754

1.4946 1.3290 1.4094 1.4495 1.3994 1.4025 1.5983 1.4544 1.4331 1.4501 1.3647 1.4143 1.4146 1.4775

1.400 1.3678 1.4042 1.2975 1.3172 1.3449 1.4280 1.5463 1.4049 1.4568 1.4822 1.5460 1.5432 1.5075

1.3837 1.3331 1.3344 1.3670 1.3658 1.4498 1.4080 1.3141 1.2633 1.2791 1.3425 1.3970 1.3738 1.1592

1.4555 1.4200 1.4574 1.3763 1.4552 1.4837 1.3548 1.3768 1.3481 1.3398 1.4288 1.5021 1.4306 1.4497

1.3180 1.3762 1.3726 1.4482 1.4567 1.4731 1.4985 1.4117 1.3157 1.3036 1.4573 1.4626 1.3628 1.5557

1.4538 1.3603 1.4281 1.4630 1.4628 1.5458 1.4099 1.3014 1.3145 1.3101 1.2875 1.3754 1.4348 1.1638

1.2680 1.3200 1.4167 1.3723 1.4667 1.4835 1.5041 1.3570 1.4048 1.3723 1.3786 1.3779 1.4985 1.5368

Table 3 Selection of electrode sites FD algorithms

Electrode groups

Proposed method Higuchi’s method Katz’s method

LST1, LST2, LTD1, LTD2, LTD3, LST4 RST1, RST4, LTD3, LST4, LTD2, RST2 LTD2, LTD3, RST4, RTD2, RTD6, RST1

Table 4 Calculation of variation FD algorithms

Coefficient of variation

Proposed method Higuchi’s method Katz’s method

2.1430 2.4820 2.3510

with respect to FD values, which are presented in Table 2. The ictal portions show greater dimension values and pre-ictal portions show lesser dimension values. The portion which has fractal dimension between 1 and 1.5 is the pre-ictal portions and which has the fractal dimension between 1.5 and 2 is known as ictal portion. The fractal dimension values discriminate the seizure and non-seizure portions. The EEG data observed from 14 electrode positions are clustered into four groups. The four groups obtained from our method and the fourth cluster contains the minimum variance electrode positions, which is presented in Table 3. The Fig. 2(b)–(d) shows the behavior of electrodes group performed by Katz, Higuchi and SMR algorithms. The electrodes sites in the fourth cluster are LST1, LST2, LTD1, LTD2 and LTD3. Computing the coefficient of variation (a relative measure of the differences of the values of various groups of electrodes from a measure of central tendency of electrode groups, while comparing the coefficient of variation of two or more series of data, the series having lesser coefficient of variation is less variable, more stable, more uniform and more consistent) among the electrode groups of three methods, the minimum variation appears in the proposed method which is presented in Table 4. The six electrodes express the related

brain activity (recurrent seizure portion) and the surface of the above electrode portions are the critical brain sites estimated by our algorithm. For records from 150 to 2500 points length, runtime for our method is very small when compared with Higuchi’s and Katz’s method. If the record length is increased to 8000 points, then our algorithm becomes more faster than Higuchi’s and Katz’s. All the three methods can be run in real time. If the window length increases up to 8000 points then our algorithm performance improves and becomes faster than Higuchi’s and Katz’s method. 4. Conclusion Investigations of brain dynamics by analysis of EEG signals and computing the fractal dimension of the brain activity are used to locate the critical brain sites of recurrent seizure portions. The proposed method described in this paper is easier than the methods from nonlinear dynamics. The treatment of epilepsy starts from the localization of brain sites. Further research will lead to the prediction of the seizure before their onset using the theory of nonlinear dynamical system and fractal analysis. Acknowledgement Authors wish to thank University Grants Commission— Special Assistance Programme (UGC-SAP), Department of Mathematics, Gandhigram Rural University, Gandhigram for their support in carrying out this research work. Conflict of interest statement None declared. References [1] M.V. Sebastian, M.A. Navascues, J.R. Valdizan, Surface Laplacian and fractal brain mapping, J. Comput. Appl. Math. 189 (1) (2004) 132–141. [2] C. Tricot, Curves and Fractal Dimension, Springer, New York, 1995.

378

P. Paramanathan, R. Uthayakumar / Computers in Biology and Medicine 38 (2008) 372 – 378

[3] J. Falconer, Fractal Geometry—Mathematical Foundations and Applications, Wiley, New York, 2003. [4] T. Higuchi, Approach to an irregular time series on the basis of the fractal theory, Physica D 31 (1988) 277–283. [5] P. Grassberger, I. Procaccia, Characterization of strange attractors, Phys. Rev. Lett. 50 (5) (1983) 346–349. [6] M. Katz, Fractals and the analysis of waveforms, Comput. Bio. Med. 18 (3) (1988) 145–156. [7] D. Iasemidis, D.-S. Shiau, J.C. Sackellares, P.M. Pardalos, A. Prasad, Dynamical resetting of the human brain at epileptic seizures: application of nonlinear dynamics and global optimization techniques, IEEE Trans. Biomed. Eng. 51 (3) (2004) 493–506. [8] A. Kalauzi, S. Spasic, Consecutive differences as a method of signal fractal analysis, Fractals 13 (4) (2005) 283–292. [9] G. Vachtsevanos, R. Esteller, A comparison of waveform fractal dimension algorithms, IEEE Trans. Circuits Syst. 48 (2) (2001) 177–183. [10] A.K. Jain, R.C. Dubes, Algorithm for Clustering Data, Prentice-Hall, Englewood Cliffs, NJ, 1988.

[11] S. Spasic, A. Kalauzi, M. Culic, G. Grbic, L. Martac, Fractal analysis of rat brain activity injury, Med. Biol. Eng. Comput. 43 (4) (2005) 345–348. Ponnaiyan Paramanathan was born in Thanjavur, India, in 1981. He received the M.Sc. degree in 2004 and B.Ed. degree in 2005 from the Gandhigram Rural University, Gandhigram, Tamil Nadu, India. Since 2006, he has been with the Gandhigram Rural University, India, where he is research scholar in Mathematics. His research interests focus on mathematical modeling, fractal theory and biomedical engineering. Ramasamy Uthayakumar was born in Dindigul, India, in 1967. He received the Ph.D. degree from the Gandhigram Rural University, Gandhigram, Tamil Nadu, India, in 2000. Since 1998, he has been with the Gandhigram Rural University, India, where he is senior lecturer in Mathematics. His research interests focus on mathematical modeling, fractal theory and biomedical engineering.