J. Cent. South Univ. Technol. (2011) 18: 834-839
DOI: 10.1007/s11771-011-0770-y
Nonlinear characteristics and functional analysis of masseter electromyography
HU Xin(胡昕)1,2, WU Xiao-bo(吴箫博)3, ZOU Bo(邹波)4,
FU Yu-li(傅予力)1, YU Yi(喻意)5, LU Guang-wen(卢广文)3
1. School of Electronic and Information Engineering, South China University of Technology,
Guangzhou 510641, China;
2. Wujing Guangdong General Hospital, Guangzhou 510507, China;
3. School of Biomedical Engineering, Southern Medical University, Guangzhou 510515, China;
4. Guanghua College of Stomatology, Sun Yat-sen University, Guangzhou 510515, China;
5. School of Foreign Languages, Guangdong University of Technology, Guangzhou 510006, China
? Central South University Press and Springer-Verlag Berlin Heidelberg 2011
Abstract: The C-C method was adopted to analyze the nonlinear characteristics of masseter electromyography (EMG) signals and the chaotic degree by the largest Lyapunov exponent (LLE) of different genders and sides. First, the embedding dimension and the delay time were obtained through this method, then the phase space was reconstructed to resume the chaotic attractor and determine the LLE. The result shows that the trajectory of attractor is denser than Chen’s attractor, and the LLE is positive, which means that not only the signal has the character of chaos, but also the chaotic degree of masseter EMG is relatively high. According to the value of the LLE, the chaotic degree of men’s masseter EMG is higher than that of women’s; when the dentition is normal, the chaotic degree of two sides is almost the same. Then, a conclusion can be deduced that if the LLE of both sides are in great difference, the unilateral mastication is likely to exist, which means that the nonlinear characteristics of masseter EMG can be applied to predict the unilateral mastication.
Key words: nonlinear characteristic; masseter electromyograph; C-C Method; phase space reconstruction; attractor; Lyapunov exponent
1 Introduction
It has been recently recognized that the EMG signal exhibits high-dimensional deterministic chaos [1-2]. Traditional approaches to analyze EMG signals always regard EMG signals as linear, smooth or smooth in time periods, so the precise qualitative analysis cannot be conducted by the data obtained by these traditional approaches. Hence, methods of nonlinear dynamics analysis have been introduced to EMG to get a better insight into the complex signals [3]. The EMG signal of masseter is the electrical manifestation of the neuromuscular muscles activation with muscle contraction. And masticatory muscles which are controlled mostly by the masseter nerves play a key part in generating a complicated and reflexible movement called chewing. Therefore, to understand the complex dynamical process of masticatory muscles while chewing, the application of nonlinear dynamical measures of masseter EMG is of great importance.
Correlation dimension and the LLE are two most important nonlinear invariant measures of nonlinear system. If we estimate the Lyapunov exponents for the unknown dynamical system, we first need to reconstruct the phase space which captures the original attractor [4]. From the attractor dynamics, we can estimate the correlation dimension and other quantities to see whether the scalar time series is chaotic or stochastic. Therefore, the attractor reconstruction is the first stage in the chaotic time series analyses [5]. A popular method to reconstruct the phase space is the delay-coordinate reconstruction. The values of the embedding dimension and the delay time should be decided first when adopting this method, and then finding an optimal delay time and a proper embedding dimension becomes the important basis for further nonlinear analysis. The C-C method based on the correlation integral [5] is a regular means to analyze the time series, which gives the delay time τd and the delay time window τw simultaneously by the application of correlation integral. Although the C-C method has no enough theoretical support, it runs well in practice and shows some advantages, such as rather simple operation, lower algorithm complexity, need for less data and better robustness [6].
By reconstructing the phase space, attractors that reflect the chaotic characteristic can be obtained; then a method from small data sets is adopted to calculate LLE, which is a measure of the rate at which nearby trajectories diverge [7]. If the LLE is positive, it suggests that the signals have the chaotic characteristic.
In the whole process, the precise calculation of the phase space reconstruction and the LLE are the two indispensible points to analyze the EMG signals of masseter while chewing.
2 Materials and methods
2.1 Acquisition of signals
Ten male and ten female volunteers, from age 19 to 25, with an average age of 22.5, whose teeth and periodontium are in good condition, were involved in this study. BL-420F biomechanism system of Chengdu TME Technology Co, Ltd serves as the equipment.
2.2 Preprocessing
Time domain and power spectrum (Fig.1) can be obtained by the preprocessing as elimination DC, elimination power line interference and lowpass filtering (LPF). The power spectrum shows that when the bilateral masseter muscles are in intercuspal position, the maximum power ranges from 125 to 250 Hz, which stands for the good conditions of myotility and the jaw stability.
Grassberger-Procaccia (G-P) algorithm [8] is used for the evaluation of the correlation dimension, which is the quantification of geometrical self-similarity of a chaotic system. The G-P method is used to calculate the correlation integral and to verify that the signal we obtained is not the noise. Figure 2 shows the correlation integral curves which descend in different dimensions exponentially. Thus, according to the method of exponential attenuation, the signal we obtained is not the noise.
2.3 Algorithm of C-C method
Correlation function used in the C-C method can estimate delay time τd and delay time window τw simultaneously [5]. Thus, embedding dimension m can be obtained according to the formula τw=(m-1)τ. τd ensures that each component of time series xi is interdependent, and does not depend on m, while τw depends on m, if τw remains fixed, τ varies with m. τw is the maximum time of data dependence.
Fig.1 Time domain (a) and power spectrum (b) of masseter EMG
Fig.2 Correlation integral curves descending in different dimensions exponentially
The correlation integral is the following function [9]:
(1)
where r is the search radius, N is the data number of the time series, M=N-(m-1)τ is the number of embedded points in m-dimensional space, θ(x) is a Heaviside function: θ(x)=0 for x<0 and θ(x)=1 for x≥0, and dij= ||Xi-Xj||∞ denotes the sup-norm.
The correlation integral is a cumulative distribution function, which denotes the probability of distance between any pairs of points in the phase space which is not greater than r. The distance is denoted by the sup-norm of the difference between two vectors.
For m=2, 3, 4, 5, ri=iσ/2, i=1, 2, 3, 4, N=3 000. Interpret the statistic S(m, N, r, t) as the serial correlation of a nonlinear time series:
(2)
(3)
(4)
(5)
Let t equal or be smaller than 200, and look for the first zero crossing of 
or the first local minimum of 
to find the first locally optimal time for independence of the data, which gives τd=tτs corresponding to the first local maximum t. Meanwhile, we look for the minimum of Scor(t) corresponding to the first independent overall maximum time window τw=tτs of time series.
2.4 Method from small data sets
When reconstructing the phase space, we can determine the LLE from the observed time series by looking at some “fiducial” orbit, which is a measure of the rate at which nearby trajectories diverge [10-11].
To each point Yi in the phase space, calculate distance dj(i):
(6)
And the estimated formula by SATO et al is
(7)
Take the logarithm of each side:
j=1, 2, …, M (8)
For general j, we have the equalizing value of y(i):
(9)
where q is the number of all non-zero dj(i); the regression line is made by the least square method, then the slope of the line is the LLE λ1.
3 Results and discussion
3.1 Embedding dimension and delay time
Figure 3 shows the delay time τd calculated by the C-C method, three dots in the figure represent the first attenuation times, which correspond to the delay times. The calculation results of embedding dimension and delay time are listed in Table 1.
Fig.3 C-C method to estimate delay time
3.2 Phase space reconstruction
When calculating the embedding dimension and delay time, the delay-coordinate method was adopted to resume the chaotic attractor:
i=1, 2, … (10)
where m is the embedding dimension and τ is the delay time. Figure 4 shows the attractor in the phase space.
In the phase space, the curves fall into a specified orbit after being stretched and folded, and this is the attractor. With the time evolution, the densities of the trajectories increase after being constantly stretched and folded [12]. When compared with Fig.4 and Fig.5, it is obvious that the trajectories of attractor as Fig.4 shows is far much denser than Chen’s attractor (Fig.5) after being stretched and folded. As it is known that the chaotic attractor is one of the typical characteristics of chaotic phenomena, the comparison shows that the masseter EMG is in a high chaotic degree.
3.3 Largest Lyapunov exponent
When combining the C-C method and the method from small data sets to calculate the LLE of EMG, it is fast, easily realized, and has good robustness for the data size, embedding dimension, time delay and the noise level [13-15]. First, the phase space is reconstructed to analyze the real masseter EMG signals, and in order to find an optimal embedding for reconstructing, the C-C method is adopted. By reconstructing the phase space, attractors that reflect the chaotic characteristic can be obtained; then a method from small data sets is adopted to calculate the LLE.
Table 1 Delay times and embedding dimensions
Fig.4 Attractor
The LLE quantitatively describes the exponential expansion and contraction speed of the adjacent tracks of the phase space, and estimates the chaos level of the systems. The value of LLE indicates the chaos level of signals. Table 2 lists the LLE of the 20 volunteers masseter EMG. The means of the LLE of male and female are calculated as follows: the right side and left side of male is 18.244±1.303 and 17.882±1.75 respectively; and of female, 14.866±1.395 and 14.839± 1.82, respectively. The LLE of men approaches to 18, and 15 of women, and both are far higher than 2.2—the LLE of Chen’s attractor, which means that the masseter EMG has a very high degree of chaos. Then, the analysis of variance is adopted to identify the effect of factors. It not only applies to the multi-level single-factor test, but also to the multi-factors test. In the analysis of variance, under the combination of each factor, the test is repeated, which is called Two-factor Repeated Analysis of Variance. Mauchly’s test of sphericity was used to judge whether there were relations among the repeated factors. If P is equal to or smaller than 0.05, then there was a correlation between factors. Greenhouse-Geisser corrected results should be taken [16]. In this experiment, the LLE is analyzed with two factors—gender and side. Table 3 indicates gender has a remarkable effect on the LLE (F=28.203, P=0.000); side doesn’t have effect on the LLE (F=0.281, P=0.603). Gender and side do not have interactive effect on the LLE (F=0.206, P=0.655).
Fig.5 Chen’s attractor
In addition, when the dentition is normal, the chaos degree of two sides is almost the same. However, for the unilateral chew, the chew ability of two sides will be different, so the chaos degree is also different, which results in the difference of the LLE. However, by the traditional potential method, such information implied in the EMG signals is difficult to be seized. So, the research shows that it is feasible to apply chaotic analysis to diagnose unilateral chew.
Table 2 Values of LLE calculated by method of small data sets
Table 3 Statistical analyses of LLE
4 Conclusions
1) Since the LLE is positive, it indicates that the signals have the chaotic characteristic; and since the LLE of both sides are almost the same, it indicates the chaos degree of two sides are almost the same. But the LLE of male is higher than female’s, that is to say, the chaos degree of men’s masseter EMG is higher than that of women’s.
2) Gender has a significant effect on the LLE of the masseter EMG, while no sufficient evidence supports that side has an appreciable effect on it. Side and gender simultaneously form no mutual interaction and affect little on the final results.
3) For the long period of unilateral chewing leads to the dissymmetric movement of the bilateral masseter, the complexity of the EMG signals is different, which definitely causes the large difference between the LLE of the two sides. Therefore, the comparative analysis of the LLE of the two sides has a great potential value to predict the unilateral chew.
References
[1] NIEMINEN H, TAKALA E P. Evidence of deterministic chaos in the myoelectric signal [J]. Electromyogarphy Clinical Neurophysiology, 1996, 36: 49-58.
[2] SWIE Y W, SAKAMOTO K, SHIMIZU Y. Chaotic analysis of electromyography signal at low back and lower limb muscles during forward bending posture [J]. Electromyogarphy Clinical Neurophysiology, 2005, 45: 329-342.
[3] XIE Hong-bo, ZHENG Yong-ping, GUO Jing-yi. Detection of chaos in human fatigue mechanomyogarphy signals [C]// The 31st Annual International Conference of the IEEE EMBS. Minneapolis: IEEE, 2009: 4379-4382.
[4] YANO M, HOMMA N, SAKAI M, ABE K. Phase-space reconstruction from observed time series using Lyapunov spectrum analysis [C]// Proceedings of the 41st SICE Annual Conference. Osaka: SICE, 2002: 701-706.
[5] KIM H S, EYKHOLT R, SALAS J D. Nonlinear dynamics, delay times, and embedding windows [J]. Physica D, 1999, 127: 48-60.
[6] CAI Wei-dong, QIN Yi-qing, YANG Bing-ru. Selection of delay time window and delay time in phase space reconstruction [C]// International Conference on Computational Intelligence and Security. Beijing: IEEE, 2007: 526-530.
[7] SHAO Chen-xi, FAN Jin-feng, LI Shao-bin. Nonlinear characteristics and qualitative analysis of sleep EEG [C]// The 1st International Conference on Bioinformatics and Biomedical Engineering. Hefei: IEEE, 2007: 794-797.
[8] GRASSBERGER P, PROCACCIA I. Characterization of strange attractors [J]. Physical Review Letters, 1983, 50(5): 346-349.
[9] GRASSBERGER P, PROCACCIA I. Measure the strangeness of strange attractors [J]. Physica D, 1983, 9: 189-208.
[10] L? Jing-hu, LU Jun-an, CHEN Shi-hua. Analysis and application of chaotic time series [M]. Wuhan: Wuhan University Press, 2002: 85-92. (in Chinese)
[11] ROSENSTEIN M T, COLLINS J J, DELUCA C J. A practical method for calculating largest Lyapunov exponents from small data sets [J]. Physica D, 1993, 65: 117-134.
[12] LI Qing-du, CHEN Shu. An improved algorithm for Lyapunov exponents of fractional-order system [C]// Proceedings of the 29th Chinese Control Conference. Beijing: 2010: 300-303. (in Chinese)
[13] WOLF A, SWIFT J B, SWINNEY H L. Determining Lyapunov exponents from a time series [J]. Physica D, 1985, 16: 285-317.
[14] SONG Wei, HOU Jian-jun, LI Zhao-hong, HUANG Liang. Chaotic system and QR factorization based robust digital image watermarking algorithm [J]. Journal of Central South University of Technology, 2011, 18(2): 116-124.
[15] ZHU Cong-xu, CHEN Zhi-gang. A novel spatial domain digital watermarking algorithm based on chaotic map [J]. Journal of Central South University: Science and Technology, 2005, 36(2): 272-276. (in Chinese)
[16] QIU Hong, JIN Guo-qin, JIN Ru-feng, ZHAO Wei-kang. Analysis of variance of repeated data measured by water maze with SPSS [J]. Journal of Chinese Integrative Medicine, 2007, 5(1): 20-31. (in Chinese)
(Edited by PENG Chao-qun)
Received date: 2010-12-08; Accepted date: 2011-03-20
Corresponding author: HU Xin, PhD candidate; Tel: +86-13226600225; E-mail: Sean_Human@hotmail.com
