1. Introduction
1.3. Organization of the dissertation
This dissertation consists of six chapters presenting methods for meditation-EEG re-search based on univariate and multivariate AR models. Figure 1-1 illustrates the hierarchy associating different chapters.
In the beginning, this chapter makes an introduction to this study and describes the main aim. Chapter 2 reports a computerized scheme Subband-AR EEG Viewer that provides a comprehensive view of the meditation EEG record. The scheme is mainly designed to trace
the time-varying spectral characteristics in meditation EEG and could be modified for spe-cific requirements or applications. Two algorithms are introduced for slow α-rhythm detec-tion and meditadetec-tion EEG interpretadetec-tion.
Chapter 1 Introduction
Chapter 2 Subband-AR EEG Viewer
Chapter 5 Spatiotemporal synchronization Chapter 3
Online α-rhythm detection
Chapter 4 Meditation EEG
interpreter
Chapter 6
Discussion and Conclusion
univariate AR model multivariate AR model
Figure 1-1: Chapter hierarchical structure.
A method for α-rhythm detection was implemented in a real-time manner in Chapter 3.
In addition, an alternative strategy was proposed to investigate the human brain in response to external flash stimuli during Zen meditation course. The flash-light stimulus was to be applied upon emergence of the frontal α-rhythm. This study was developed and conducted for further understanding the effect of Zen meditation on human perception.
Chapter 4 describes a software system “meditation EEG interpreter” based on the second algorithm presented in Chapter 2. This interpretation system can identify such EEG patterns
bility, especially for rejecting the artifacts like the baseline drift and EMG (electromyograph) interference. It also allows a quick overview of an enormous amount of EEG data.
Besides the univariate AR model employed in meditation EEG research, Chapter 5 fo-cuses on the multivariate autoregressive (mAR) model feasible for quantitative study of spati-otemporal behavior in multichannel EEG. The spatispati-otemporal synchronization index based on residual covariance matrix of mAR model was proposed. The index measures the degree of synchronization among neighboring channels of a local brain area. A reduction of synchroni-zation (or, significant desynchronisynchroni-zation) in the brain areas results in a relatively high index. It might provide a quick inspection of spatiotemporal characteristics for various meditation sce-narios.
Significant differences in EEG rhythmic patterns, meditation scenarios, spatiotemporal characteristics, and visual evoked potentials have been observed between experimental and control groups. The last chapter makes a summary of the results obtained from previous chap-ters.
Chapter 2－
Meditation EEG Overview Based on Subband Features Quantified by AR Model
Measurement does not necessarily mean progress. Failing the pos-sibility of measuring that which you desire, the lust for measurement may, for example, merely result in your measuring something else - and perhaps forgetting the difference - or in your ignoring some things because they cannot be measured.
~ George Udny Yule
his chapter reports a computerized scheme Subband-AR EEG Viewer that pro-vides a comprehensive view of the meditation EEG record. The scheme was mainly designed to trace the varying spectral characteristics in meditation EEG.
Following sections illustrate the main ideas and the methods adopted in the scheme. Two al-gorithms modified from the scheme are introduced for particular applications. An on-line im-plementation of the Subband-AR-EEG Viewer is also drawn in the end of this chapter.
T
2.1. Subband-AR-EEG Viewer
The scheme is focused on monitoring the time-varying characteristic frequency in
medi-tation EEG. The AR model is applied to the subband component to quantify the characteristic frequency. Consider that the EEG signal, , is generated by an autoregressive (AR(p)) process driven by unit-variance white noise [Theodoridis and Koutroumbas 1999]. A pth or-der all-pole model is formulated by
]
where is the order of the AR model. Based on this model, the spectrum of the EEGs can be obtained if the coefficients are known. Several techniques have been proposed to estimate parameters . We apply the autocorrelation method in which the AR coeffi-cients are determined by solving the autocorrelation normal equations
p
where
ε
is the modeling error andγ
x[k] is the estimated autocorrelation function defined below:As addressed previously, a higher model order (for example, p ranges from 6 to 14) is normally required to better estimate the low-frequency component in EEGs. According to (2-2), the coefficients of the AR(6) model are determined by {
γ
x[k]⏐0≤k≤6}. As demonstrated in Table 2-1, the values of the autocorrelation functionγ
x[0] -γ
x[6] for θ and δ rhythms are too close to distinguish between each other, whereas the coefficients for α and β exhibit signifi-cant deviation.Table 2-1: Values of autocorrelation function, γx[0] - γx[6], computed for δ, θ, α, and β rhythms.
γx[0] γx[1] γx[2] γx[3] γx[4] γx[5] γx[6]
δ 1.000 0.986 0.948 0.891 0.822 0.746 0.665 θ 1.000 0.985 0.946 0.889 0.823 0.752 0.676 α 1.000 0.946 0.792 0.562 0.287 −0.0004 −0.271 β 1.000 0.872 0.553 0.166 −0.155 −0.323 −0.324
In fact, even increasing the model order to p=12 cannot discriminate δ rhythm from θ rhythm. The dominant pole pairs for δ and θ modeled by AR(12) are 0.964∠±0.131 and 0.959∠±0.129, respectively, which results in a close estimate of the spectral frequencies (symbol ‘∠’ denotes the phase in radian). In addition, the computational time required by AR(12) becomes four-fold compared with that for AR(6).
A downsampling process ensures the AR modeling better characterizes the low fre-quency activities. This is revealed by the autocorrelation function values
γ
x[0] -γ
x[6] esti-mated for δ and θ rhythms downsampled by 8:γ
x[k]δ: {1.00, 0.48, −0.18, −0.48, −0.39, −0.03, 0.06}, (2-4a)γ
x[k]θ: {1.00, −0.16, 0.26, 0.06, −0.42, 0.02, −0.23}. (2-4b)Moreover, according to Gabor’s uncertainty principle [Oppenheim et al. 1998], downsampling operation improves frequency resolution which is desired for narrow-band EEG. We accord-ingly employed subband-filtering prior to the frequency analysis by AR modeling.
In summary, EEG signals are firstly decomposed into different subband components by downsampling and filtering. Then the characteristic frequency (root frequency) of each sub-band component is estimated by the AR(2) model. The entire scheme is called the Sub-band-AR EEG Viewer and it can be illustrated by the tree-structured filter banks shown in Fig.
2-1. In the Subband-AR EEG Viewer, a linear-phase lowpass FIR filter H(z) with cutoff fre-quency 30Hz is used as an anti-aliasing filter before the downsampling operation. Then the AR(2) model is applied to the decimated signal. The filtering-and-downsampling process is repeated until the equivalent cutoff frequency equals 1.875Hz.
Differing from the wavelet decomposition, the Subband-AR-EEG-Viewer only employs the lowpass linear-phase filter. A nonlinear-phase filter which is commonly employed in wavelet analysis [Vaidyanathan 1993] would distort the temporal information of EEG and consequently affect the AR model coefficients estimated by the autocorrelation method.
The Subband-AR-EEG-Viewer structure (Fig. 2-1) can be rearranged into a six-channel filter banks shown in Fig. 2-2 with decimation ratios which are powers of two. Fig. 2-3 shows the frequency responses of the lowpass filters in the filter banks. The cutoff frequencies of H1(z), …, H5(z) are, respectively, 30Hz, 15Hz, 7.5Hz, 3.75Hz, and 1.875Hz (sampling rate:
200Hz).
Classification algorithm
Figure 2-1: The structure for the Subband-AR-EEG-Viewer.
Figure 2-2: The equivalent six-channel system of the subband filter bank.
It should be noted that the cutoff frequencies approximate the upper boundaries of the four well-known EEG rhythms⎯ β (14–30Hz), α (8–13Hz), θ (4–7Hz), and δ (below 3Hz).
Therefore, changes of the characteristic frequency in meditation EEG can be traced by quan-tifying the root frequency (fr) of each subband filtered component. For example, when frs of
output1, output2 and output3 are all within the range 8–13 Hz, the dominant pattern of this windowed segment is identified, to a great degree, as the α rhythm. When frs of output1 and output2 are greater than 15Hz and fr of output4 is between 4Hz and 7Hz, the particular seg-ment most likely contains θ intermixed with β rhythm.
After the subband decomposition, the AR(2) model coefficients are computed. An AR(2) model can be expressed as
]
The model coefficients are directly computed by
⎟⎟⎠
where ]
γ
x[k is the autocorrelation function estimated by (2-3).The characteristic frequency, also called the root frequency, of outputi can be estimated from the phase of the pole, or the root of the model equation (2-5) expressed in the frequency domain. After obtaining the model coefficients, the conjugated pole pair is
2
Therefore, the root frequency fr can be formulated as
]
Because the root frequency is much smaller than the sampling frequency, the result of
sin−1x can be approximated by x. For example, the α rhythm having a higher frequency of 13Hz results in a normalized radian frequency of 0.13π (assume fs = 200Hz). The approxima-tion only causes a 2.8% deviaapproxima-tion from the true value. It should be noted that output2-output6
are the results of downsampling (Fig. 2-2), the root frequency fr,i should be further divided by 2i-1. According to equations (2-6) to (2-9), root frequency of each subband component de-pends on γx[0], γx[1], and γx[2]. Tracking the root frequency of each subband component provides an efficient way to illustrate the time evolution of the characteristic frequency in meditation EEG.
2.2. Experimental setup and protocol
The meditation EEG signals were recorded using 8-channel SynAmps amplifiers (manu-factured by NeuroScan, Inc.) connected to the Pentium MMX-166 (MHz) PC. Due to hard-ware limitation, we observed EEG characteristics in different regions with the 8-channel, unipolar recording montage involved the electrode sites at O1, O2, Oz, Cz, Pz, F3, F4 and Fz (Fig. 2-4). The common reference was the linked M1-M2 (mastoid electrodes). The EEG sig-nals were pre-filtered by a bandpass filter with passband 0.3–30 Hz, and digitized at 200 Hz sampling rate. Each recording lasted for 45 minutes, which consisted of the first 5-minute background EEG (the subject sat in a normal relaxed position with eyes closed) and the 40-minute meditation EEG. During the meditation session, the subject sat, with eyes closed, in the full-lotus or half-lotus position. Each hand formed a special mudra (called the Grand Har-mony Mudra), laid on the lap of the same side. The subject focused on the Zen Chakra and the Dharma Eye Chakra (also known as the “Third Eye Chakra”) in the beginning of meditation till transcending the physical and mental realm. The Zen Chakra is located inside the third ventricle, while the Dharma Eye Chakra is located at the hypophysis [Lo et al. 2003].
Figure 2-4: 8-channel recording montage.
2.3. Slow alpha rhythm detection (SARD) algo-rithm
The slow α-rhythm in meditation EEG was found to be related to the early stage of medi-tation. The first algorithm was designed particularly to detect the slow α-rhythm in the long-term meditation EEG record. The ideas and methods proposed can be implemented in different ways, oriented towards a specific purpose.
2.3.1. The algorithm
Frequency of the α-rhythm ranges from 8Hz to 13Hz. The slow α-rhythm is a particular pattern, normally below 10Hz, that was observed in some experimental subjects at the mind-focusing stage of meditation. The Subband-AR-EEG-Viewer can be reduced to the structure (Fig. 2-5(a)) for tracking the slow α-rhythm and the slow alpha rhythm detection (SARD) algorithm is illustrated in Fig. 2-5(b). According to analytical reasoning and practical
experience, output1 in combination with output3 highly enhances the effectiveness of slow α-rhythm detection. Note that fr,1
acts as an index to screen out the high-frequency component,
and fr,3 is employed in the classification as a major reference. The SARD algorithm depicts that the slow α-rhythm pattern is detected when both root frequencies satisfy the following criteria:fr,1 < 14Hz, and 8Hz < fr,3 <10Hz.
2.3.2. Simulation
While only examining output1 (up to 30Hz) with the criterion 8Hz < fr,1 <10Hz, it often fails to identify the noise-contaminated slow α-rhythm. Figure 2-6 demonstrates the noise-immunization capability of our scheme. When a pure 9Hz sinusoid (Fig. 2-6(a)) was partially contaminated by uniformly distributed random noise (Fig. 2-6(b)), the AR model did not recognize the noise-contaminated slow α-rhythm segment based on the criterion 8Hz<fr,1<10Hz (Fig. 2-6(c)). The results in Fig. 2-6(d) show that the SARD algorithm suc-cessfully detected the slow α-rhythm under poor environment (SNR=8dB).
To verify the performance, we first analyzed a simulated signal of 4-sec duration. The signal shown in Fig. 2-7(d) was generated by connecting three short-duration, ampli-tude-modulated sinusoids, respectively, with frequencies 9Hz, 15Hz, and 5Hz (Figs.
2-7(a)-(c)). The window length is 0.5 sec (100 samples), moving at a step of 0.25 sec. As shown in Fig. 2-7(d), the SARD algorithm effectively detected the slow α-rhythm pattern.
(a)
(b)
Figure 2-5: The Subband-AR-EEG-Viewer modified for slow α-rhythm detection: (a) structure and (b) algorithm.
1 sec
slow α
slow α slowα
(a)
(b)
(c)
(d)
SNR=8dB
Figure 2-6: Capability of noise immunization of slow α-rhythm detector. The black bold lines mark epochs identified as slow α-rhythm.
slow α
1 sec
(a)
(c) (b)
(d) 9 Hz
15 Hz
5 Hz
Figure 2-7: Slow α-rhythm detection from the time-varying rhythmic activities. The signal in (d) is simulated by adding three short-duration, amplitude-modulated sinu-soids with frequencies (a) 9Hz, (b) 15 Hz, and (c) 5 Hz.
2.3.3. Identification of slow alpha rhythm
Next, the SARD algorithm was applied to the meditation EEGs (channel O1). A 10-second segment shown in Fig. 2-8 was analyzed with the same implementing parameters as used in Fig. 2-7. Bold lines above the signal indicate the slow α-rhythm detection. As shown in Fig. 2-8(a), amplitude variation often affects the recognizability of slow α-rhythm.
It results in a crack in the first bold line. On the other hand, a transient slow α-rhythm may be of little significance. We thus designed a post-processor to further refine the result. It removed segments shorter than 0.3 sec and fused cracks shorter than 0.3 sec (Fig. 2-8(b)).
1 sec
slow α
(a)
(b)
Figure 2-8: Detection of slow α-rhythm in meditation EEG (Channel O1). Bold lines above the signal indicate the slow α-rhythm detection. (a) A crack in the first bold line is caused by variation of the α amplitude. (b) A post-processor removes the second lines with duration shorter than 0.3 sec (insignificant activity) and fuses the crack shorter than 0.3 sec.
Detection of specific EEG patterns is important in identifying various meditation states.
In addition, it may serve as a preprocessing stage in such tasks like the EEG segmentation or interpretation. Based on the Subband-AR-EEG-Viewer, we devised a strategy for meditation EEG interpretation. Details are illustrated below.
2.4. Long-term meditation EEG interpretation (MEEGI) algorithm
Changes of the characteristic frequency in meditation EEG may be a key feature for un-derstanding various states of consciousness during meditation. We therefore developed a logical strategy implemented in the computerized MEEGI algorithm to segment the EEG into sections with different frequencies. The results, illustrated as a running gray-scale chart, clearly reveal the evolution of a characteristic frequency during meditation.
2.4.1. The algorithm
In the following, we present interpretation results of the meditation EEG based on five spectral features frequently observed during meditation: (1) The χ features which is a slow waveform intermixing with high-frequency rhythms, (2) δ, (3) θ, (4) α, and (5) β features.
The χ feature mostly appears at the transition from one EEG rhythm to another. To provide a long-term visual record, the five spectral features are displayed by different gray tones. The gray tones from the darkest to the lightest indicate, respectively, the χ, δ, θ, α, and β features.
In this task, the structure of Subband-AR-EEG Viewer can be reduced to that shown in Fig.
2-9.
Figure 2-9: The Subband-AR-EEG Viewer modified for meditation EEG interpretation.
The MEEGI algorithm examines each windowed segment to check the following criteria in order:
Criterion-χ: fr,3<7Hz<fr,1 and ⏐p3⏐>0.8;
Criterion-δ: fr,1<7Hz and fr,4<3.5Hz;
Criterion-θ: fr,1<7Hz;
Criterion-α: 7Hz< fr,1<14Hz and 7Hz<fr,3;
Criterion-β: 7Hz<fr,1;
where p3 is the AR(2) pole of output3. The criteria checkup is ordered according to a sound logic realizing the subband filtering method. The root frequency fr,1 is used to differentiate between 0-7Hz and 7-30Hz EEG bands, while the fr,3 , fr,4 and ⏐p3⏐ is employed in the sub-sequent discrimination process. The length of p3 can be considered as an indication of the sig-nificance of the root frequency. Because the χ wave represents an intermixed signal composed of both low- and high-frequency components, we impose restrictions on the range of ⏐p3⏐ to ensure the significance of the low frequency component. The flowchart of the MEEGI algo-rithm is shown in Fig. 2-10.
2.4.2. Simulation
To verify the effectiveness of feature recognition, the MEEGI algorithm was first applied to a simulated signal. As displayed in Fig. 2-11(e), the signal was formed by connecting five segments of δ, θ, χ, α, and β patterns. We assumed the sampling rate is 200Hz. This signal can be simulated by the pole placement method, that is, by placing each pole in the corre-sponding frequency band (Table 2-2) and adding Gaussian noise. The transition from θ to β normally results in a compound pattern like χ. The running gray-scale chart (Fig. 2-11(e)) above the simulated sequence successfully signals the temporal patterns.
χ
α
θ
β δ
Figure 2-10: Meditation EEG interpretation (MEEGI) algorithm.
1 sec
(a)
(b)
(c)
(d)
(e) (1)
(2) (3) (5) (4)
(1) χ, (2) δ, (3) θ, (4) α, and (5) β.
Figure 2-11: Simulation results (a) δ activity (b) θ activity (c) β activity (d) α activity (e) simulated signal and classification results.
Table 2-2: Locations of poles of the simulated signal
δ θ α β
Pole
locations 0.98∠0.04 0.98∠0.16 0.98∠0.4 0.88∠0.63
The above simulation demonstrates the feasibility of the scheme and algorithm in Fig.
2-9 and 2-10 for automatically identifying different EEG rhythms and revealing time-varying schema of the simulation. In empirical data, as more complex rhythmic patterns are involved, discrepancies between experienced EEG interpreters may occur. Methodology development thus focused on reliable recognition of some key features in meditation EEG analysis. Figure 2-12 demonstrates the robustness of the Subband-AR-EEG Viewer for identifying even the little jittering of β rhythms embedded in the high-amplitude slow activity.
(2) (2) (2)
Figure 2-12: Subband-AR-EEG Viewer applied to the EEG signal for feature recogni-tion.
2.4.3. Long-term meditation EEG interpretation
When applied to the long-term meditation EEG, this algorithm is particularly robust for automatic interpretation with no need to determine the implementing parameters. Fig. 2-13 displays three running gray-scale charts for two experimental subjects (Fig. 2-13(a) and (b)) and one control subject (Fig. 2-13(c)). Since Watts [Watts 1957] has reported the variation of alpha amplitude and frequency in the frontocentral region during meditation, we thus selected channel F3 for further analysis. The error rate was approximately 8.7% in comparison with the results of naked-eye examination by an experienced EEG interpreter. Both meditators have been practicing Zen Buddhist meditation for more than eight years. The control subject sat in a normal, relaxed position with the eyes closed. During the meditation session, the two
meditators exhibited different meditation scenarios. Meditation EEG of subject 2k1019p was apparently dominated by β rhythm, sometimes transforming into a short-duration α rhythm.
According to the post-experimental interview, the subject did not always stay in the Alaya consciousness and occasionally went back to normal consciousness.
The chart in Fig. 2-13(a) reveals this scenario for 20-minute. In Fig. 2-13(b), subject 2k0830a exhibited a large portion of χ activities. In our meditation EEG experiment, EEG signals of some meditators indeed were found to be characterized by large-amplitude, slow-drifting rhythms interwoven with high-frequency tiny jiggles. Meditators with this kind of EEG characteristic normally have their meditation process wandering among normal con-sciousness, subconsciousness (subliminal consciousness), and Alaya consciousness [Lo et al.
2003]. Compared with the experimental group, EEGs collected from the control subject are normally dominated by α rhythm, as illustrated in Fig. 2-13(c). Note that this subject was drowsy in the experiment, resulting in occurrence of θ and δ rhythms.
2.4.4. On-line implementation of the Subband-AR-EEG Viewer
Due to its simplicity, the algorithms provides a robust tool for on-line processing re-quired for meditation EEG interpretation as well as the biofeedback scheme in the BCI (brain-computer interface) research. Currently, we have implemented the MEEGI algorithm under Simulink (MathWorks, Inc., Natick, MA) with Real-Time Workshop on a Pentium-M 1.4 (GHz) notebook (Fig. 2-14(a)). By generating a real-time code with Real-Time Workshop, the algorithm can be downloaded to the kernel and run in a real-time manner under Windows [Guger 2001]. The classification results were displayed on a monitor. As shown in Fig.
2-14(b), the height of each bar reflects the power percentage of the corresponding EEG rhythm within a 2-sec frame. We can thus monitor the meditator’s state in a real-time manner that enables the development of a more subtle correlation between the EEG characteristics
and the meditation scenario.
χ δ θ α β
300 sec
600 sec
900 sec
1200 sec
(a)
χ δ θ α β
300 sec
600 sec
900 sec
1200 sec
(b)
χ δ θ α β
300 sec
600 sec
(c)
Figure 2-13: Running gray-scale charts (Channel F3) for two meditators (a) subject 2k1019p, and (b) subject 2k0830a, and (c) one non-meditator (control) subject
(a)
(b)
Figure 2-14: On-line implementation of the Subband-AR-EEG Viewer. (a) Computer 1 executes the MEEGI algorithm (Fig. 2-10), and Computer 2 displays the classification results, as illustrated in (b). The height of each bar reflects the power percentage of the corresponding EEG rhythm within a 2-sec frame.
Chapter 3－
Investigation of Visual Perception under Zen-Meditation
If the brain were so simple we could understand it, we would be so simple we couldn't.
~Lyall Watson
ne topic of interest in the meditation study is the evoked potentials (EP) (or event-related potentials, ERP) of a practitioner under meditation, which includes auditory evoked potential (AEP), somatosensory evoked potential (SEP), visual evoked potential (VEP), and so on. Each parameter is meaningful to the respective perception function. Due to unusual perceptions often experienced during meditation, human brain in
ne topic of interest in the meditation study is the evoked potentials (EP) (or event-related potentials, ERP) of a practitioner under meditation, which includes auditory evoked potential (AEP), somatosensory evoked potential (SEP), visual evoked potential (VEP), and so on. Each parameter is meaningful to the respective perception function. Due to unusual perceptions often experienced during meditation, human brain in