141
790

 

.THS'

cQODq

 

LIBRARY
Michigan State
University

 

 

 

This is to certify that the
thesis entitled

Measuring the Phase Synchrony of Brain Signals Using Time—

MS.

Frequency Distributions

presented by

Westley Evans

has been accepted towards fulfillment
of the requirements for the

degree in Electrical Engineering

 

 

Major Professor’s Signature

mm 000;;

 

Date

MSU is an Affirmative Action/Equal Opportunity Employer

 

 

PLACE IN RETURN BOX to remove this checkout from your record.
TO AVOID FINES return on or before date due.
MAY BE RECALLED with earlier due date if requested.

 

DATE DUE DATE DUE DATE DUE

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

5/08 K1IProj/Acc8PreleIRC/DateDue‘indd

 

MEASURING THE PHASE SYNCHRONY OF BRAIN
SIGNALS USING TIME—FREQUENCY DISTRIBUTIONS

By

Westley Evans

A THESIS

Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of

MASTER OF SCIENCE
Electrical Engineering

2008

ABSTRACT

MEASURING THE PHASE SYNCHRONY OF BRAIN SIGNALS
USING TIME-FREQUENCY DISTRIBUTIONS
By

Westley Evans

Quantifying the phase synchrony of brain signals is important for the study of
large—scale interactions in the brain. Current methods for computing phase syn-
chrony are limited to amplitude-dependent correlation methods, frequency-domain
coherence methods, and wavelet transform methods. However, many of these meth-
ods fail to take the time—varying nature of real life signals into account. To address
this issue, we propose two measures of phase synchrony based on Cohen’s class
of time—frequency distributions. The first proposed method measures the phase
synchrony between a pair of signals using a time-varying estimate of the phase dif—
ference between the two signals. The phase difference estimate is computed using
a reduced interference complex time-frequency distribution. The second proposed
method measures the phase synchrony among a group of signals by quantifying the
frequency locking among the signals using a time-frequency based estimation of the
instantaneous frequency. The instantaneous frequency maps of individual signals
are combined to obtain the instantaneous frequency histogram as an estimate of the
amount of frequency locking across the signals. This analysis is then extended to
the estimation of frequency locking across multiple electrodes and multiple trials.
Results are shown for both methods using synthetic signal models and electroen—

cephalogram (EEG) data collected from control and schizophrenic subjects.

TABLE OF CONTENTS

LIST OF TABLES ................................. v
LIST OF FIGURES ................................ vi
KEY TO SYMBOLS AND ABBREVIATIONS ................. ix
CHAPTER 1
Introduction ..................................... 1
1.1 Problem Statement ............................ 1
1.2 Contributions of Thesis .......................... 2
1.3 Organization of Thesis .......................... 3
CHAPTER 2
Background ..................................... 4
2.1 Electroencephalogram .......................... 4
2.2 Functional Integration .......................... 5
2.3 Time-Frequency Distributions ...................... 6
2.3.1 Cohen’s Class of TFDS and the Wigner Distribution ...... 8
2.3.2 Reduced Interference Distributions ............... 8
2.3.3 Rihaczek Distribution ...................... 10
2.3.4 Complex Continuous Wavelet Transforms ............ 11
CHAPTER 3
Bivariate Phase Synchrony ............................. 16
3.1 Background on Phase Synchrony Measures ............... 16
3.2 Reduced Interference Rihaczek Distribution (RID-Rihaczek) ..... 19
3.3 Time-Varying Phase Difference ..................... 20
3.4 Measuring Bivariate Phase Synchrony .................. 21
3.4.1 Single Trial Phase Synchrony .................. 22
3.4.2 Phase Locking Value ....................... 22
3.4.3 Comparison Between the Wavelet and Rihaczek based Phase
Synchrony Measures ....................... 23
CHAPTER 4
Multivariate Phase Synchrony ........................... 27
4.1 Instantaneous Frequency ......................... 27
4.2 Estimating IF of Multicomponent Signals ................ 28
4.3 Instantaneous Frequency Histograms .................. 30
4.4 Multiple Trial IFHs ............................ 31
4.5 Quantifying IFHs ............................. 34
4.6 IFH of Synthetic Signals ......................... 34

iii

CHAPTER 5

Results ....................................... 36
5.1 Gamma Band Synchrony in Schizophrenic Subjects .......... 36
5.2 Bivariate Phase Synchrony in Schizophrenic Subjects ......... 37
5.3 Multivariate Phase Synchrony in Schizophrenic Subjects ....... 44

CHAPTER 6

Conclusions and Future Work ........................... 59
6.1 Contributions ............................... 59
6.2 Future Work ................................ 60

6.2.1 Multivariate Method Improvements ............... 60
6.2.2 Volume Conduction ........................ 61
6.2.3 Directional Measures ....................... 62
6.2.4 Functional Networks ....................... 62
BIBLIOGRAPHY ................................. 64

iv

LIST OF TABLES

Table 2.1 The frequency bands of an EEG signal ................ 5

Table 5.1 The correlation averages in five different frequency bands for the
four most significant centroids of each control subject. ...... 57

Table 5.2 The correlation averages in five different frequency bands for the
four most significant centroids of each schizophrenic subject. . . . 58

Table 5.3 The weighted average of the correlation averages for the four most
significant centroids of each control and schizophrenic subject in
five different frequency bands. The final row contains the dif-
ference between the average control subject’s correlation average
and the average schizophrenic subject’s correlation average. . . . 58

Figure 2.1

Figure 2.2

Figure 2.3

Figure 2.4

Figure 3.1

Figure 3.2

Figure 3.3

Figure 4.1

Figure 4.2

Figure 4.3

Figure 5.1

Figure 5.2

Figure 5.3

LIST OF FIGURES

a) The real component of a complex exponential chirp signal that
increases in normalized frequency from 0.05 to 0.35, b) the mag—
nitude of the FT for the chirp signal, and c) the Wigner TFD of

the same chirp signal. ........................ 12
a) The non-uniform time-frequency resolution of the wavelet
transform and b) the uniform time-frequency resolution of Co-
hen’s class of TFDs. ......................... 13
a) Wigner distribution and b) Choi—Williams distribution for the
sum of two Gabor logons ....................... 14
The ambiguity domain for the sum of two Gabor logons ...... 15
a) Magnitude of Rihaczek distribution and b) Magnitude of re—
duced interference Rihaczek distribution .............. 24
a) Real component of the two complex exponentials and b) Com-
parison of the theoretical and estimated time-varying phase dif-
ference of the two complex exponentials. .............. 25
Average PLV for two chirps with a constant phase difference: a)
wavelet distribution, b) Rihaczek distribution. ........... 26

Example of estimating IF on a multicomponent signal. Black
indicates the expected IF, and gray indicates the estimated IF. . . 31

Synthetic signals (a) 1‘1(t) and (b) 172(t) ............... 35

The IFH of a pair of synthetic signals with common IF from 0 to
1 seconds and from 2 to 3 seconds.

Average PLV between electrode pairs in the 7 band of the P300
time window for the four control subjects. The diagonal of the
figure, which is the PLV of an electrode with itself, is set to zero
rather than one to improve scaling. ................. 38

Average PLV between electrode pairs in the '7 band of the P300
time window for the four schizophrenic subjects. The diagonal of
the figure, which is the PLV of an electrode with itself, is set to
zero rather than one to improve scaling.

The Average PLV difference between control and schizophrenic
subjects ................................. 40

vi

Figure 5.4

Figure 5.5

Figure 5.6

Figure 5.7

Figure 5.8

Figure 5.9

Figure 5.10

Figure 5.11

Figure 5.12

Figure 5.13

Figure 5.14

Figure 5.15

Figure 5.16

Average PLV with respect to electrode F2 for the four control
subjects ................................. 41

Average PLV with respect to electrode F2 for the four
schizophrenic subjects. ........................ 42

Difference of the average PLV with respect to electrode Fz be-
tween the control and schizophrenic subjects. ........... 43

Topographic map of the average PLV with respect to Fz for four
control subjects. ........................... 44

Topographic map of the average PLV with respect to Fz for four
schizophrenic subjects. ....................... 45

Topographic map of the difference of the average PLV with re-
spect to Fz between control and schizophrenic subjects ....... 46

The four most significant cluster centroids of control subject. 1.
The number of IFHs belonging to each cluster: (a) 15, (b) 14, (c)
3, (d) 3. ................................ 49

The four most significant. cluster centroids of control subject 2.
The number of IFHS belonging to each cluster: (a) 12, (b) 6, (c)
5, (d) 4. ................................ 50

The four most significant cluster centroids of control subject 3.
The number of IFHS belonging to each cluster: (a) 21, (b) 6, (c)
5, (d) 4. ................................ 51

The four most. significant cluster centroids of control subject 4.
The number of IFHS belonging to each cluster: (a) 60, (b) 7, (c)
3, (d) 3. ................................ 52

The four most significant cluster centroids of schizophrenic sub-
ject 1. The number of IFHS belonging to each cluster: (a) 27, (b)
21, (c) 5, (d) 3. ............................ 53

The four most significant cluster centroids of schizophrenic sub-
ject 2. The number of IFHS belonging to each cluster: (a) 32, (b)
8, (c) 2, (d) 1 .............................. 54

The four most significant cluster centroids of schizophrenic sub-
ject 3. The number of IFHS belonging to each cluster: (a) 5, (b)
3, (c) 3, (d) 2 .............................. 55

vii

Figure 5.17 The four most significant cluster centroids of schizophrenic sub-
ject 4. The number of IF Hs belonging to each cluster: (a) 40, (b)
9, (c) 5, (d) 5 .............................. 56

viii

KEY TO SYMBOLS AND ABBREVIATIONS

CPT: Continuous Performance Task

CW: Choi-Williams

CWT: Continuous Wavelet. Transform

EEG: Electroencephalograrn

EP: Evoked Potential

FT: Fourier Transform

IF: Instantaneous Frequency

IFH: Instantaneous Frequency Histogram

PLV: Phase Locking Value

RID: Reduced Interference Distribution

SPS: Single Trial Phase Synchrony

STFT: Short-Time Fourier Transform

TFD: Time-Frequency Distribution

ix

CHAPTER 1

INTRODUCTION

1 .1 Problem Statement

Cognitive acts require the integration of numerous functional areas widely dis-
tributed over the brain and in constant interaction with each other. These inter-
actions among separate areas of the brain are referred to as functional integration
[1]. Recent research [2] suggests that the phase synchronization of signals trans-
mitted along reciprocal connections between different areas of the brain allow for
functional integration to occur. Numerous clinical trials have demonstrated that
the phase synchronization of brain signals is related to functional integration [2].
In addition, results from Rodriguez et a1. [3] Show that phase synchrony is a bet-
ter indicator of functional integration than amplitude dependent measures such as
energy.

Furthermore, with the advance of neuroimaging technology, it is now possible to
identify the oscillations of neuronal networks at high temporal and spatial resolu-
tions, using multichannel recordings. Simultaneous recording of multiple oscillations
between different cortical regions offers insight into how distributed neuronal oscil-
lations interact with each other. These interactions, which can be used to study
large-scale functional integration, are transient, time-varying, and frequency spe-
cific.

Therefore, there is a. need for dynamic, high resolution phase synchrony measures
to quantify the degree of functional integration that is occurring within the brain.
Measures of phase synchrony must. account for the nonstationary nature of brain
signals. While some measures that use the Hilbert transform, short-time Fourier

transform, or the continuous wavelet. transform account for the nonstationarity of

the signals, they suffer from limited time and frequency resolution and scaling.
In addition, current measures are limited because they can only quantify pairwise
synchrony between neuronal oscillations, and cannot directly assess the large-scale
interaction between groups of neuronal signals. For these reasons, there is a need for
multivariate, high time-frequency resolution phase distributions based on Cohen’s

class for quantifying phase synchrony.

1.2 Contributions of Thesis

In this thesis, two new methods to measure phase synchrony are proposed. The
first method uses a high resolution time-frequency phase distribution based on Co-
hen’s class to measure time-varying phase estimation and phase synchrony between
a pair of signals. The reduced interference Rihaczek distribution is proposed be-
cause of its abilities to measure complex energy and reduce cross—terms. The high
resolution time-frequency distributions based on Cohen’s class outperform other
time-frequency distributions such as wavelets. Wavelets, on the other hand, use a
non-uniform time—frequency resolution where the frequency resolution is high at low
frequencies and low at high frequencies.

The second proposed method for measuring phase synchrony quantifies multivari-
ate phase synchrony across groups of neuronal oscillations. Unlike previous multi-
variate methods [4], the proposed method uses distributions from Cohen’s class with
uniformly high time—frequency resolution. The proposed method is based on using
the relationship between phase synchrony and frequency locking. Since phase and
frequency are directly related to each other, two signals are synchronous whenever
their instantaneous frequency is approximately the same. Using this relationship,
an instantaneous frequency (IF) estimation method in the time—frequency domain is
proposed to quantify the amount of synchronization between groups of signals. The

IF estimates for different. signals are combined using an instantaneous frequency

histogram (IFH) to indicate the amount of synchronization. The proposed method
is then extended for multiple electrodes and multiple trial EEG recordings.

The proposed methods extend the state of the art in phase synchrony measures.
Improved phase synchrony measures will help neuroscientists and psychologists un-
derstand functional integration within the brain, which is a vital aspect of proper

brain function.

1.3 Organization of Thesis

This thesis is organized as follows. Chapter 2 reviews the background on EEG
signals and time-frequency distributions. Chapter 3 presents and compares several
techniques to measure the phase synchrony between a pair of signals. Chapter 4 pro—
poses a new method to compute the phase synchrony over multiple electrodes and
trials. Chapter 5 presents the results of applying the proposed methods to EEG sig-
nals collected during a study involving control and schizophrenic subjects. Finally,
Chapter 6 discusses the contributions of this thesis along with some suggestions for

future work.

CHAPTER 2

BACKGROUND

2.1 Electroencephalogram

The electroencephalogram (EEG) is a neuroimaging tool that measures the elec-
trical potential over different areas of the brain. EEGs are non-invasive since the
measurements are taken from a set of electrodes placed upon the scalp of the subject
[5]. Each electrode measures the voltage signal generated by its respective area of
the brain. The placement of the electrodes follows a predefined layout where each
electrode is given a specific name [6].

EEG provides excellent temporal resolution compared to other neuroimaging
methods at the expense of spatial localization. Some of the challenges of extracting
useful information from EEG signals include low signal to noise ratio and volume
conduction. EEG signals contain a large amount of “noise” because of the numerous
activities that occur within the brain. Since the EEG measures voltages generated
by the brain from an electrode on the scalp, the electric currents generated by the
different areas of the brain can mix together during conduction through the skull
and scalp resulting in the electrodes measuring a modified voltage signal from the
original brain-generated signal. While most researchers choose to ignore volume
conduction, others have developed techniques to reverse its effects on EEG signals.
For example, the inverse EEG solution estimates the cortical activity from the EEG
signals [2].

The EEG is commonly used in psychology and psychiatry to evaluate different
psychopathologies. In these studies, a visual or audio stimulus is usually presented to
the subject in order to measure the evoked potential (EP). In order to obtain reliable

measurements of the evoked potential, the same stimulus is repeated multiple times

and the average EEG waveform, called the event related potential, is commonly
analyzed.

EEG studies typically focus on specific time frames and frequency bands of
interest. Some common time components of interest include P50, N100, P200, and
P300. The ’P’ indicates a positive deflection of the EEG during the time frame,
and the ’N’ indicates a negative deflection of the EEG. The number following the
direction of the response corresponds to the latency in milliseconds of the deflection
after the stimulus is applied to the subject. For example, the P300 time component
includes the positive deflection of the EEG around 300 ms. In order to observe
this deflection, the time frame from 200 to 600 ms after the stimulus should be
analyzed. Furthermore, the frequency content of an EEG signal is roughly divided
into 5 different frequency bands. Table 2.1 displays the approximate boundary

frequencies of the different bands [6].

Table 2.1. The frequency bands of an EEG signal.

 

 

 

 

 

 

 

Band Name Approx. Frequencies (Hz)
Delta (6) < 4

Theta (0) 4-8

Alpha (0) 8-13
Beta ((3) 13-30

Gamma ('7') 30-55

 

 

 

2.2 Functional Integration

The brain adheres to two fundamental principles of functional organization, func-
tional integration and functional segregation, where the integration within and
among specialized areas is mediated by effective connectivity [1]. Functional seg-

regation is the functional specialization of a group of neurons. A single group of

neurons cannot perform cognitive processing independently. Instead, cognitive pro-
cessing requires that the functionally and anatomically segregated areas of the brain
must combine their functionality [7]. Recognizing a familiar face, for example, re-
quires the functionality of both memory and vision. Functional integration refers
to the interactions among the segregated neuronal areas and how these interactions
depend upon cognition [1]. Functional integration is not static but is a dynamic and
context dependent process.

The most likely mechanism for functional integration uses the strong reciprocal
connections between different areas of the brain. It is believed that the phase syn-
chronization of signals passed along these reciprocal connections allow for functional
integration to occur in the brain [2]. The technical aspects of phase synchrony will
be further discussed in Chapters 3 and 4. Numerous experiments on humans and
animals using invasive and non-invasive techniques have demonstrated that phase
synchronization plays a role in functional integration [8]. Rodriguez et a1. [3] were
one of the first to non—invasively demonstrate that phase synchrony occurs in the hu-
man brain during a cognitive task. They measured the phase synchrony of the EEG
signals in the gamma band generated by a human subject during a facial perception
task using ‘Mooney’ faces. The phase synchrony of the signals was found to be
significantly greater when a face was perceived than when a face was not perceived.
The type of phase synchronization recorded by EEGs is considered large—scale syn-
chronization because of the long distances (> 1cm) between the electrodes. Besides
the work performed by Rodriguez, Von Stein et a1. [9] demonstrated the presence

of long-range synchrony in the beta band of EEG signals from humans.

2.3 Time-Frequency Distributions

Most real life signals including EEG signals have time-varying frequency charac—

teristics. Therefore, conventional signal processing methods based on time domain

modeling or frequency domain analysis are inadequate for the study of these signals.
For example, the Fourier transform (FT) takes a signal from the time domain to the
frequency domain. However, the transform cannot distinguish the frequency content
of a signal at different time instances because it assumes the frequency characteris-
tics of the signal are constant over the duration of the signal. In order to properly
analyze the nonstationary EEG signals, a different signal analysis method must be
utilized.

Time-frequency distributions (TFDS) represent a signal over both the time and
frequency domains. This allows TFDS to handle signals with time-varying frequency
characteristics. Figure 2.1 demonstrates the advantage of using a TFD over a FT
for a complex exponential chirp signal defined as the following: 3(t) = exp(j7r(0.05+
0.3t/400)t). The normalized frequency of the chirp signal increases linearly from 0.05
to 0.35. Given Figure 2.1(b), one is not able to determine that the instantaneous
frequency of the signal is changing with time. All information about the time-
varying frequency characteristics of the signal is lost by the FT. However, with the
Wigner distribution in Figure 2.1(c), one is able to clearly determine that the signal’s
instantaneous frequency is linearly increasing in time. The Wigner distribution will
be defined in Section 2.3.1.

The time frequency analysis methods of wavelets and TFDS have been used
in the study of EEG signals in applications involving neurology and the study of
psychopathologies. The work of Samar [10], Raz [11], and Basar [12] used wavelets
to analyze event related potentials and resulted in a better understanding of the
oscillatory characteristics of EEG. Several authors have used TFDS [13, 14, 15] or
wavelets [16, 17] to detect and analyze epileptic seizures. Polikar has used wavelets
to detect Alzheimer’s disease [18, 19]. Other psychopathologies that have been
analyzed using wavelets include schizophrenia [20] and obsessive compulsive disorder

[21].

2.3.1 Cohen’s Class of TFDS and the Wigner Distribution

Cohen's class of distributions are bilinear TFDS that are expressed as 1[22]:

C(t,w) = if: /// d)(0,'r)s(u + %)s*(u — %)€j(6u—6t—Tw) du d0 dT (2.1)

where the function q‘)(6, T) is the kernel function and s is the signal. The kernel
completely determines the properties of its corresponding TF D.

One of the most popular TFDS is the Wigner distribution with a kernel function
of q§(6, T) = 1. The Wigner distribution is a real-valued distribution that describes
the energy of the signal over time and frequency, simultaneously.

The major differences between Cohen’s class of TFDS compared to other time-
frequency representations such as the wavelet transform are the nonlinearity of the
distribution, energy preservation, and the uniform resolution over time and fre-
quency. The wavelet transform provides a representation over time and scale where
the frequency resolution is high at low frequencies and low at high frequencies. Al-
though this property makes wavelet transform attractive in detecting high frequency
transients in a given signal, it inherently imposes a non-uniform time-frequency tiling
on the analyzed signal and thus results in biased energy representations. Cohen’s
class of bilinear TFDS on the other hand assumes constant resolution over the en-
tire time—frequency plane. Figure 2.2 compares the time-frequency resolution of the

wavelet transform and Cohen’s class of TFDS.
2.3.2 Reduced Interference Distributions

One of the disadvantages of the Wigner distribution is the existence of cross-terms
for multicomponent signals. Figure 2.3(a) displays the Wigner distribution for the

sum of two Gabor logons. A Gabor logon is defined as a time and frequency shifted

 

1 All integrals are from —00 to 00 unless otherwise stated.

_ _ 2
Gaussian window, g(t) = (Tl )(1/4)exp( (t t ) exp(jw t). For this example,
0 TI 20 O

the sum of two Gabor logons signal, 3(t) 2 g(t) + 92(t), is used with t0,1 = 50,
w0,1 = 0.7, 1‘02 == 150, and “’02 = 0.3. As can be seen from the figure, apart from
the individual energy distributions of the two logons, there are interference terms
in the middle of the two auto-terms corresponding to the cross-terms. Attenuation
of the cross-terms is possible through proper kernel selection and design. Distribu-
tions that reduce the presence of cross-terms are referred to as reduced interference
distributions (RIDs).

To aid the kernel design process Equation 2.1 can be redefined into the following

form:

C(t,w) = -4—71r§f/(MB,T)A(6,T)€_j(6t+Tw)de0 (2.2)

where A(6, T) is the ambiguity function of the signal defined as:
A(6, T) = /s(u + %)s*(u — %)ej0udu (2.3)

In the ambiguity domain the auto—terms of every component are located near the
origin, while the cross-terms are located away from the origin. Therefore, in order
to remove the cross—terms, a reduced interference kernel must be chosen with a
passband around the origin. As a result, reduced interference kernels meet the
following condition:

96(6, T) << 1 for 67' >> 0 (2.4)

Figure 2.4 shows the ambiguity domain for the same sum of two Gabor logons. The
desired auto—terms are located at the origin at the center of the figure, and the
cross-terms are located off the axes.

The Choi-Williams kernel (CW) with q$(6,T) = exp(—(6T)2/0) is an example

kernel function that removes cross-terms by filtering in the ambiguity domain. The

resulting CW distribution can be written as the following:

 

2 .
cue) = // exp(-(6;) )A(6, 7—)e—J(9t+w)drd9 (2.5)
W—J

CW kernel

The value of a can be adjusted to achieve a desired trade-off between resolution
and the amount of cross-terms retained. Figure 2.3(b) displays the CW distribution

with a = 4 for the sum of two Gabor logons.
2.3.3 Rihaczek Distribution

Most of the members of Cohen’s class are real valued energy distributions such as
the Wigner and Choi—Williams distributions. These distributions are successful at
describing the energy of the signal over both time and frequency. However, they
do not carry any information about the phase of the signal. For this reason, they
cannot be directly used for describing the phase information in an individual signal
and estimating the phase synchrony between two signals.

In 1968, Rihaczek introduced the complex energy distribution and gave a plau-
sibility argument based on physical grounds [23]. The Rihaczek distribution has a
kernel function of (M6, T) = exp( j6T / 2). Using the Rihaczek kernel reduces Equation

2.1 to the following:

 

C(t,w) = s(t)S*(w)e—jwt (2.6)

1
{2?
which is known as the Rihaczek distribution. This distribution measures the com-
plex energy of a signal around time t and frequency w. The complex energy density
function provides a fuller appreciation of the properties of phase-modulated signals
that is not available with other time—frequency distributions. The Rihaczek distribu-
tion is a bilinear, time and frequency shift. invariant, complex-valued time-frequency
distribution belonging to Cohen’s class. The Rihaczek distribution provides both a

time—varying energy spectrum as well as a phase spectrum, and thus is useful for

10

estimating the phase synchrony between any two signals.
2.3.4 Complex Continuous Wavelet Transforms

While complex continuous wavelet transforms (CWTs) are not. part of Cohen’s class
of TFDS, they can be used as a benchmark during comparisons with the Rihaczek
distribution. A popular complex CWT is the Morlet wavelet transform which is
defined as the inner product of the signal, 3(a), with the complex wavelet function
[24]:

VV(t.w) = /s(u)\Il*(u)du (2.7)

where \I'(u) is a Gaussian window shifted in time and frequency, and scaled propor-

tional to the frequency:

u_ 2
Wu) = (E expwu — t» expr—L—Ug—t)» (2.8)

with a being inversely proportional to w. The time-frequency scaling of the Mor-
let wavelet. is nonuniform as previously demonstrated in Figure 2.2 because the
variance, 0, changes with frequency. Also, the Morlet wavelet. is similar to the
short-time Fourier transform (STFT) where the window is defined as h(u — t) =
(fi)exp(jwt)exp(—(u — t)2/02), which scales with respect to frequency, rather
than as a Gaussian window, h.('u — t) = 1/(27r02)exp(—(u — t)2/(202)), which is

constant with respect. to frequency.

11

 

 

 

 

 

1_ /\ ‘
2:03.] l][[[[[[ [‘ ‘ l[[
e l: [“l[[ '[
0 5 [g] [] [l [U M [] [
_10 i so 160 7150 200
60 'Tige'
§ 40 . .

 

 

 

0.2 0.4 0.6 0.8 1.0
Normalized Frequency

 

 

 

(C)
5‘ 1.0 - B
C
§o.8-
(D
'u: 0.6-
'0
g 0.4»
2 o

 

50 100 150 200
Time

Figure 2.1. a) The real component of a complex exponential chirp signal that in—
creases in normalized frequency from 0.05 to 0.35, b) the magnitude of the FT for

the chirp signal, and c) the Wigner TF D of the same chirp signal.

 

 

 

 

 

 

 

 

 

 

Frequency
Frequency

 

 

 

 

 

 

 

 

 

 

 

 

 

Time Time

(a) (b)

Figure 2.2. a) The non—uniform time-frequency resolution of the wavelet transform
and b) the uniform time-frequency resolution of Cohen’s class of TFDS.

13

Normalized Frequency

Normalized Frequency

.0
N

.4
O

P
on

.0
o:

.0
A

.0
N

_\.
O

.o
oo

.0
c:

.o
4:.

Wigner distribution

 

I

I

 

I H

 

 

50

100 1 50 200
Time
(a)

Choi-Williams distribution

 

 

I I

 

 

50

100 1 50 200
Time

(b)

Figure 2.3. a) Wigner distribution and b) Choi-Williams distribution for the sum

of two Gabor logons

14

100

50

-100

Ambiguity Domain

 

 

 

   

 

 

 

 

 

Auto-terms\ _ ; ;
E ’ j. -: 15
’ i - 10
i 3 - 5
E E : i
1 I 1 O
-50 0 50 100

Theta

Figure 2.4. The ambiguity domain for the sum of two Gabor logons.

15

CHAPTER 3

BIVARIATE PHASE SYNCHRONY

As previously mentioned in Section 2.2, cognitive acts require the integration of nu—
merous functional areas widely distributed over the brain and in constant interaction
with each other. This integration is typically quantified by measuring the amount of
phase synchrony between different regions of the brain. Bivariate phase synchrony
measures the relation between the temporal structures of two signals regardless of
the signal amplitude. Two signals are considered to be synchronous if their rhythms

coincide.

3.1 Background on Phase Synchrony Measures

The most common measures used to quantify dependence between two signals is
time domain cross-correlation [8] and spectral coherence [24]. The time-domain

cross-correlation between two signals, :r:(t) and g(t), is defined as the following:

:r(t) *y(t) 2: /m*(u)y(t + u) du (3.1)

A high cross-correlation between two signals implies that the signals are similar.
Cross-correlation is like covariance except that the signals are not shifted by their
means. Also, the computation of cross-correlation is comparable to the computation
of convolution except neither of the signals is reversed during cross-correlation [25].
The spectral coherence between the signals, T(t) and y(t), at a frequency of interest
f is expressed as:

ley(f)|

(f) = (3-2)
9 [5331?(f) ' Syy(f)l1/2

 

16

where S$y( f) is the cross—spectral density between the signals a: and y computed as
the Fourier transform of the time-domain cross-correlation (Equation 3.1). However,
these measures have limitations for two reasons. First, they assume stationarity
of the underlying signals whereas most real life signals, including EEGS, are not.
Second, coherence is a measure of spectral covariance and does not separate the
effects of amplitude and phase from each other. The ideal phase synchrony measure
should be able to separate the phase and amplitude effects from each other and
take the nonstationary nature of brain activity into account [26]. The phase and
amplitude terms of the input signals in Equations 3.1 and 3.2 can be separated into
different components, but this change will not correct for the nonstationarity of the
input signals.

In order to address these limitations, two different measures for quantifying phase
synchrony have been proposed. In both methods, the amount of synchrony between
two signals is usually quantified by first estimating the instantaneous phase of the
individual signals around the frequency of interest. The first method employs the
Hilbert transform to get the signal into an analytic form and estimates instantaneous
phase directly from its analytic form [27]. In order to be able to estimate the
instantaneous phase of a signal from its analytic form, one has to make sure that
it is a narrow-band signal. For this reason, the Hilbert transform based phase
synchrony measure first bandpass filters the signal around a frequency of interest
and then uses the Hilbert transform to get the instantaneous phase. This is an
indirect way of obtaining the frequency dependent phase estimates and is not exact.
The second approach, on the other hand, computes a time-varying complex energy
spectrum using either the CWT with a complex Morlet wavelet [24, 28] or the STFT
[29] to get the phase of the signal.

For both phase synchrony measurements, the goal is to obtain an expression for

the signal in terms of its instantaneous amplitude, a(t), and instantaneous phase,

17

@(t), at the frequency of interest as follows:

3(t,w) = a(t)€Xp(j(wt + 95(0)) (3-3)

This formulation can be repeated for different frequencies and the relationships
between the temporal organization of two signals, 51(t) and 82(t), can be observed

by their instantaneous phase difference.

em = |n¢1(t) —- m¢2(t)l (3.4)

where 9791(25) and Q5205) are the phases of 31(t) and 32(t), respectively. In addition,
71 and m are integers that indicate the ratios of possible frequency locking. Most
studies are only concerned about the case when n = m = 1. In neuroscience, the
focus is on the case when the phase difference is approximately constant over a
limited time window. This is defined as a period of phase locking between two
events. Phase locking is an indicator of the dynamic phase relationship between two
oscillatory neuronal sources independently of their amplitude.

It has been observed that the two methods are similar in their results with the
time-varying spectrum based methods giving sharper phase synchrony estimates
over time and frequency, especially at the low frequency range [27]. Although
the wavelet and STFT based phase synchrony estimates address the issue of non-
stationarity, they suffer from a number of drawbacks. In the case of the wavelet
transform, a representation over time and scale is obtained where the frequency res-
olution is high at low frequencies and low at high frequencies as was demonstrated by
Figure 2.2. Despite this property making the wavelet transform attractive in detect-
ing high frequency transients in a given signal, it inherently imposes a. non-uniform

time-frequency tiling on the analyzed signal and thus results in biased energy rep-

18

resentations and corresponding phase estimates. In the case of STF T, there is a
tradeoff between time and frequency resolution due to the window. For these rea-
sons, there is a need for high time-frequency resolution phase distributions that can

better track the dynamic changes in phase synchrony.

3.2 Reduced Interference Rihaczek Distribution (RID-Rihaczek)

As mentioned in Section 2.3.3, the Rihaczek distribution is a complex-valued TFD
from Cohen’s class that provides phase information about the signal. Unfortunately,
the Rihaczek distribution, like the Wigner distribution, suffers from the existence
of interference terms for multicomponent signals. For any signal in the form, s(t) =

31(t) + 32(t), the Rihaczek distribution is:

C(W) = (81(t)S’{ (Me—1“" + 32(t)5’2*(w)e-J'wt

I
m (3.5)

+ .91(t)3§(w)e—JWt + .92(t)Sf(w)e_Jwt)

where the last two terms in the above expression are the cross-terms. These cross-
terms are located at the same time locations and occupy the same frequency bands
as the original signals. A reduced interference Rihaczek distribution (RID-Rihaczek)
must be defined to attenuate the cross-terms.

In Section 2.3.2, the CW kernel was used to filter the cross—terms in the ambiguity
domain. After applying the CW kernel to the Rihaczek distribution, the resulting
RID-Rihaczek can be written as [30]:

(9T)2 . _ '(0t+Tw)
#0 .W
CW kernel thaczek kernel

where A(6,T) is the ambiguity function of the signal as defined by Equation 2.3.

This new distribution still satisfies the marginals and preserves the energy. Figure

19

3.1 illustrates the original and the reduced interference Rihaczek distributions using

the same sum of two Gabor logons signal described in Section 2.3.2.

3.3 Time-Varying Phase Difference

The time-varying phase estimate of a signal based on the Rihaczek distribution can

be defined as the following:

_ ‘ C(t,w)
W”) ‘ “3 [Tami] ’

= arg [€j¢(t)e—j0(w)e—jwt] , (3.7)

= ¢)(t) — 9(a)) — wt

where q§(t) and 6(a)) refer to the phase in the time and the frequency domains,
respectively.

Once the time-varying phase spectrum is defined, the phase difference between
two signals, 31(t) and 32(t), with Rihaczek distributions of Cl(t,w) and 02(t,w),

respectively, can be computed as:

 

P C Lu) 0*(t,W)
(1)12(t,w) = arg 1( ) 2 ] a

[IC1(taW)l |02(taW)l
: arg "ej((¢1(t)—91(w)-wt)e-j(¢2(i)—92(w)-wt)]

= mg “ej(<<z>1<t)—¢2<t>)—<01<w)—62(w))]

’ (3.8)

 

3

= ($10) —- (152(0) — (91M - 9201))

where cbl(t) and 952(25) correspond to the phase of the time domain signals, whereas
01(w) and 02(w) correspond to the phase of the Fourier transform of the two signals,
respectively.

It can be shown that for a real-valued signal, the phase between a signal 81(t)

20

and its shifted version 31(t — to) is given by

 

‘1’12(ti’“’):‘irg |s ensue )l Isl(t—to)l|51(w)l

(
51(t) s *(t—to )e M0 (3.9)
_|:‘1(t)| |81(t-t0)| ’

"31(t)sf(w )wc th8T(t—t0)51(w)e Jwtoert]
1 7

 

= arg[

 

= —wt0

which is a linear function of frequency as expected.

Example: T‘ime- Varying Phase Tracking: In this example, two complex exponen-
tial signals are considered with a time-varying phase difference, 2:1(t) = exp(jw1t)
and 122(t) = exp(jw1(t — at2)). Figure 3.2(a) illustrates the real component of the
two signals. First, the Rihaczek distribution of the two complex exponentials is
computed and then the phase difference between them is computed using Equation
3.8. Figure 3.2(b) compares the theoretical phase difference which is a second order
polynomial and the estimated one at the frequency of interest, “’1' As can be seen
from the figure, the proposed method is successful at estimating the time-varying

phase difference with high accuracy. Images in this thesis are presented in color.

3.4 Measuring Bivariate Phase Synchrony

In most applications, the time-varying phase spectrum is not particularly useful
for measuring the synchrony between signals. Similar to the definitions given in
[24, 28, 31], phase synchrony measures based on the Rihaczek distribution time-
varying phase estimate are defined. There are two major phase synchrony measures,
single trial phase synchrony and phase locking value. Both of these measures are

bivariate and therefore can only measure the phase synchrony between two signals.

21

3.4.1 Single Trial Phase Synchrony

Single trial phase synchrony (SPS) measures the phase locking between two elec-

trodes averaged over time and is defined as:

1 t+6/2 .
_/t 6 2 GXI)(](I)12(T,w))dT (3.10)

SPS(t,w) = (5

 

where 6 is the length of the time window used for smoothing the phase difference
estimates. In previous work based on the wavelet transform, 6 is chosen as a function
of frequency since the time-bandwidth product is not a constant for the wavelet.
transform. However, in the proposed approach, the phase synchrony is based on
the time-frequency distribution which has a constant time-bandwidth product over
the whole time-frequency plane. Therefore, 6 will be a constant determined by the
signal.

SPS will always be between 0 and 1, and can be used to detect neural synchrony
between different parts of the brain. When the phase difference between the two

signals is constant over time, SPS will be equal to 1.
3.4.2 Phase Locking Value

The phase locking value (PLV) is defined between two electrodes but is averaged
over all trials. It measures the stability of phase differences across trials and is

defined as:
N

PLV(t,w) = % Z exp(j<I>(1§)(t,w)) (3.11)
k=1

. . A“. . . . .
where N is the number of trials and @g2)(t,w) 1s the time-varying phase estnnate

between two electrodes during the kth trial. Like SPS, PLV ranges from 0 to 1.

3.4.3 Comparison Between the Wavelet and Rihaczek based Phase Syn-

chrony Measures

In this example, we compare the performance of the PLV phase synchrony measure
based on the complex Morlet wavelet transform and the Rihaczek distribution for
signals with time-varying frequency content. Figure 3.3 compares the performance
of wavelet based phase synchrony measure with Rihaczek based phase synchrony
measure. With the Rihaczek based phase synchrony measure, the PLV is equal to 1
around the instantaneous frequency and smoothly tapers off as the frequency moves
away from the instantaneous frequency. With the wavelet based phase synchrony
measure, the PLV values do not taper off gradually but rather sharply and there is a
larger bandwidth around the instantaneous frequency with a phase locking value of 1.
It is also important to note that in the wavelet based phase locking value estimates,
the bandwidth around the instantaneous frequency increases as the instantaneous
frequency increases. This is due to the fact that at high frequencies the wavelet

transform has high time resolution at the expense of low frequency resolution.

23

Magnitude of Rihaczek distribution

 

 

 

 

 

 

 

 

1.0 I . i
5‘ 0.8 - -
C ..
g Q G
E» 0.6 — —
LI.
'0
a)
g 0.4- -
E Q Q
2 0.2 — 4
O l l l
50 100 150 200
Time
(a)
Magnitude of RID-Rihaczek distribution
1.0 T l l
g 0.8 - . -
E» 0.6 - —
u.
13
a)
a 0.4- -
(U
E m
0
2 0.2 - _
0 l l l
50 100 150 200
Time
(b)

Figure 3.1. a) Magnitude of Rihaczek distribution and b) Magnitude of reduced
interference Rihaczek distribution

24

 

 

 

 

 

 

 

 

1 I , ’ 'l\ I
\ .
.’ ' I \
l \
0.5 » - 3 - \
a) l ' I \
'O I - .
:3 ' l \
a 0 ' I I
E ‘ . I
< I " ‘I
'05 ' ‘. I ' ‘
\ ‘.
' I \ .
_1 L 1’ l‘ 1 ’8 u,
0 0 2 0.4 0 6 0 8 1
Time
(b)
8 I I I 1

 

Theoretical Phase Diff.
6 . . - . - . Est. from Rihaczek Dist.

 

 

 

 

Phase Difference

 

 

 

 

Time

Figure 3.2. a) Real component of the two complex exponentials and b) Comparison
of the theoretical and estimated time-varying phase difference of the two complex

exponentials.

Phase Locking Value Over 100 Trials

 

 

0.6 0.6
05 0.5
E E
g 04 g 0.4
s - s
I: u.
“D ‘U
.3 .3 0.3
E 0.3 7°
E E
2° 2°
0.2
0.2
0.1
0.1
0
20 40 60 80 100 120 20 40 60 80 100 120
Time Samples Time Samples
(a) (b)

Figure 3.3. Average PLV for two chirps with a constant phase difference: a) wavelet
distribution, b) Rihaczek distribution.

26

CHAPTER 4

MULTIVARIATE PHASE SYN CHRON Y

Chapter 3 exclusively dealt with measures of phase synchrony between two signals
at a time. Since most EEG signals involve multiple channels of data, the quantifica—
tion of the relationships between different electrodes requires the computation of all
possible pairwise synchronies. For example, 32 channel EEG will require the com-
putation of 512 different. bivariate phase synchrony combinations to determine the
ensemble synchronization of the signals. In this chapter, a method will be presented
to measure the phase synchrony between more than two signals. The proposed mul-
tivariate method is deterministic unlike the statistical methods described in Chapter
3.

The amount of research exploring multivariate phase synchrony has been rel-
atively limited. Allefeld and Kurths developed a method to measure multivariate
phase synchrony using statistical cluster analysis [32, 33]. Another technique devel-
oped by Rudrauf et al. estimates the instantaneous frequency of the signals with
wavelets and finds common frequencies among the signals [4]. In this chapter, we
will propose a multivariate phase synchrony measure based on the estimation of
the instantaneous frequency of signals using high resolution time-frequency distri-

butions.

4.1 Instantaneous Frequency

Two monocomponent. analytic signals, 2:1 (t) = (11(t)e¢1(t) and T2(t) = a2(t)e¢2(t) ,
are phase synchronous if the phase difference between the two signals is constant.

To allow for a small amount of noise in the phase of synchrtmous signals, the phase

27

difference can be approximately constant:
A9512“) 2 mqbl(t) — nq‘)2(t) z constant (4.1)

where m and n. are any two integers, and A961 2 is the phase difference between the

two signals. The derivative of Equation (4.1) can be shown as:

qubLQU) : m’d¢1(t) _ nd¢2(t)

dt dt dt

 

= mw1(t) — nw2(t) z 0 (4.2)
deb.- . .
where cal-(t) = #(t) > 0 and 'Uii(t) IS the 111stantaneous frequency (IF) of :rZ-(t)
in radians/ second. Therefore, when n = m, two monocomponent signals are phase
synchronous if they have approximately the same IF. This relation can be extended
to three or more monocomponent signals. If a set of signals all have approximately
the same IF, then all of the signals in the set are phase synchronous together [4].
The IF of a monocomponent signal can be computed by taking the time deriva-
tive of the phase of the analytic signal. Unfortunately, most real life signals, includ-
ing brain signals, are not necessarily monocomponent. They tend to be multicom—
ponent and are approximately equivalent to separable components in the following
general form: s(t) = Z k“ k(t)ej $18“). While the same concept of phase synchrony
applies to multicomponent signals, computing their lFs requires the use of more

sophisticated methods [4].

4.2 Estimating IF of Multicomponent Signals

Ideally, the IF of a multicomponent signal will be the union of the IF s of the separate
components. Numerous techniques exist to estimate the IF of a multicomponent
signal. Previous methods include time-frequency moments, adaptive recursive least
squares, and adaptive least mean square [34. 35]. Since most biological signals are

non-stationary it is not possible to estimate their IF from the signal in the time

28

or frequency domain. For this reason, a time frequency (TF) peak algorithm for
estimating IFS of multicomponent non-stationary signals is proposed.

For the purposes of IF estimation it is important. to choose ¢(6, T) such that
it corresponds to a reduced-interference distribution in order to remove the cross-
terms. The Choi-Williams distribution, previously mentioned in Section 2.3.2, is
used to estimate the IF [22]. The proposed TF peak method can be summarized as

follows:

1. Using the Choi-Williams kernel, compute the TFD of the signal, 3(t), to get
C(t, f).

2. Find the local peaks of C (t, f) using the following:

,1 “(17W :0} /\ {—1,—82gft’f) <0}

0 otherwise

B(t, f) = (4.3)

where B (t, f) is a binary-valued image with ones at the local peaks of C (t, f)

and zeros everywhere else.

3. Assign all nonzero time-frequency points into connected components. A con-
nected component D contains nonzero time-frequency points such that there
is an 8—neighborhood path containing only points in D. Two time-frequency
points are connected if B(t1,f1) = 1 , B(t2,f2) = 1 , [152 — ill 3 1, and

U2 —f1| S 1.

4. Remove any connected component, Dz" from B(t, f) if [Di] < c where [Di]
is the support. of Dir The removal of connected components with a small
support reduces the effect. of noise in s(t). The threshold, 6, is dependent on

the application and the time support of the signal.

29

5. Compute the average energy of each component, D -, as follows:

t,-=—1— 2 can (4.4)

Remove any connected component, Di, from B (t, f) if 5,; < x\. The removal
of low energy connected components reduces the effect of noise in S(t). The
threshold, A, is application dependent and should depend on the maximum

average energy component.

The remaining connected components of B(t, f) form the IF estimate, F (t, f), of
3(1) The resulting F (t, f) is a binary image with ones indicating the time and
frequency location of the IF.

Consider the following example of estimating the IF for a multicomponent signal.
The signal $(t) = sin(27r(3 + 3t)t) + sin(247rt) contains two different components,
a chirp and a sinusoid. The expected IF of the first. component is f1 (t) = 6t + 3,
and the expected IF of the second component is f2(t) = 12. Figure 4.1 displays the

expected and estimated IFS computed using the TF peak finding method.

4.3 Instantaneous Frequency Histograms

Given the fact that phase synchronous signals have similar IFS, the phase synchrony
of multiple signals or multivariate phase synchrony can be represented using an
instantaneous frequency histogram (IF H). If the collection of signals contains Q
signals, computing the IF for each signal results in Q different IF estimates, Fifi, f)
for 2' = 1, 2, . . . , Q. Summing the IFS from the collection of signals results in the
IFH, IFH(t, f) 2 23:1 [fl-(t, f). Since the IF estimate is a binary image, each value
of the IFH is a discrete, finite value in the range from 0 to Q. For a single trial
EEG consisting of multiple electrodes, the IFH represents the multivariate phase

synchrony of the trial.

30

 

15: .

 

 

1712;- n._..
1:,

5‘ 9-

C -
3 _II
g 6- __..:-r- -
I: -=-""

 

 

O . 1 1 l I 1‘
0.2 0.4 0.6 0.8 1.0
Time (sec)

 

Figure 4.1. Example of estimating IF on a multicomponent signal. Black indicates
the expected IF, and gray indicates the estimated IF.

4.4 Multiple Trial IFHS

Most EEG studies involve multiple trials of the same stimulus. Therefore, an EEG
recording with N trials will result in N IFH surfaces across the electrodes. Each of
the IFHS must be interpreted and analyzed to find general patterns of multivariate
phase synchrony among the trials. Using data reduction techniques on the IFHS will
determine representative patterns of multivariate phase synchrony across trials.
Typically, principal component. analysis (PCA) is the primary tool for the di-
mension reduction of data sets. However, because of the discrete, finite values of the
IFH, PCA is not the best option. For this reason, the K-means clustering approach
is proposed for reducing the N IF H surfaces across trials into a few, distinct IFH

components. The centroids can be restricted to the finite values of the IFH. For

31

K-means clustering, each time-frequency point of the IF HS is treated as a separate
dimension in the clustering vector space, and each of the N IFHS is treated as a
point in the clustering vector Space. For indices 2' = 1, 2, . . . , N and j = 1,2,. . . , K,

K—meanS clustering can be summarized as follows:

1. Choose the initial centroids, Ill j = median({IFHz- : IF Hz’ 6 Rj}) where Rj is

a collection of randomly selected IF HS from the set of N IFHS.

2. Assign each IF H to a cluster by finding its nearest centroid Ci 2
argminj “IFH,- — AIJ-[lg , where C,- is the cluster index that IFH,- is assigned

to.

3. Calculate the new centroids based upon the new cluster assignments using
[if j = median({IFHz- : C,- = j }) The median is used rather than the mean to

keep the centroids in the finite field.

4. If the centroids, I” j, have changed since the last iteration, go back to step 2.

If the centroids, .Mj, have not changed Since the last iteration, stop.

The results of the K-means algorithm depend on the randomly chosen initial
centroids. Multiple instances of the K-means clustering should be performed with
different initial centroids, and the resulting relationships between the IF HS can be
hierarchically clustered together to get a more stable clustering. Multiple K—means

clusterings can be combined into a connection matrix using the following technique:

1. Define A to be a N X N connection matrix that is initially a zero matrix.
2. Perform the K-means algorithm on the N IFHS to get K different clusters.

3. For each t = 2,3,. ..,N and S = 1,2, . . .,t — 1, if IFHS and IFHt are in the
same cluster such that C5 2 Ct , then increment A(S, t). The matrix A will

be upper triangular because S is always less than t.

32

4. Repeat steps 2 and 3 P — 1 more times where P is the desired number of times
that K-means clustering is performed. Each value of the connection matrix,
A(s, t), is an integer ranging from 0 to P which indicates the number of times

that IFHS and IFHt were in the same cluster.

Hierarchial clustering can now be performed on the connection matrix, A, by the

following:

1. Assign each IFH to its own individual cluster by letting the cluster index

Ci=ifori=1,2,...,N.

2. Find the location of the maximum value in A using (S_ma:r,t_ma:c) =
argmaxSJ A(S,t). If IFHS_ma;r and IFHLmam are already in the same
cluster such that. Cs_m,a;1; = Ct_ma:rv then no action is necessary. How-
ever, if IFH5_mag; and IFHLmaa: are not in the same cluster such that
Cs_7n.a1' # CLmazs then the clusters that contain IFHanax and IFHLmar

must be merged together. Merging can be performed by the following:
Ci = Csjnalf If Ci = Ct_ma;r (4.5)

for 2' = 1,2,. . . , N. The merge reduces the number of clusters by one.
3. Set A(S_ma:r, t_ma1‘) = 0.
4. Return to step 2 if the number of different clusters is greater than K.

5. Compute the centroids by using Alj = median({IFHi : Ci. = j}) for i =
1,2,...,Nandj=1,2,...,K.

Once the final clustering is complete, the significance of each centroid is deter-

mined by the number of IFHS in that cluster.

33

4.5 Quantifying IFHS

In order to quantify the phase synchrony from IFHS the correlation average is com—

puted [4]:
Z (1FH2(t.f) — IFH<t. f))
(thew
nc(nc — Uri-Tap

 

Pavg = (4.6)

where ac is the number of Signals, IV is the time-frequency window of interest, ”T
is the number of time samples in W’, and n F is the number of frequency samples
in W. The pavg is an indicator of phase synchrony per TF point. An IFH with
a higher pavg indicates more signals are phase synchronous. In order for pavg to
be one, all time-frequency points in the window of interest need to be equal to ac,
which only occurs when every Signal is phase synchronous at every time-frequency
point in the window of interest. The window of interest, W, allows for the measure

to be localized to a specific region of the TF plane.

4.6 IFH of Synthetic Signals

In this example, two sinusoidal signals are considered that are shown in Figure 4.2,
1131(t) = sin(27rf1 (t) + 27r/3) and x2(t) = cos(27rf2(t)) where f1(t) = 3u(t) + 3u(t —
2) —6u(t—3) and f2(t) = 3u(t)—3u(t—1)+6u(t—2)—6u(t—3) such that u(t) is the
unit step function. The sinusoids have the same IF of 3 Hz from 0 to 1 seconds and
6 Hz from 2 to 3 seconds. The resulting IFH in Figure 4.3 shows that the proposed

method is effective at tracking the time-varying phase synchrony.

34

 

 

 

 

 

 

 

 

 

0.5 -
9:. o — -
><
-O.5 - -
-1 1 1 I
0 0.5 1 1.5 2 2.5 3
t (sec)
(8)
3:3“ _
1.5 2 2.5 3
t (sec)
(b)

Figure 4.2. Synthetic Signals (a) r1(t) and (b) 12(t).

 

 

 

 

 

9 - i l r l I e 2
if, 6 _ "_ 1.5
5 n
C I IIt":I 1
g. 3 1“ L“- - !ll[ _
g a... 0.5

0 ~ __.._......_...~ ' ~ -

1 I I l l o
0.5 1.0 1.5 2.0 2.5 3.0
t(sec)

Figure 4.3. The IF H of a pair of synthetic signals with common IF from 0 to 1
seconds and from 2 to 3 seconds.

35

CHAPTER 5

RESULTS

In this chapter, the proposed bivariate and multivariate methods for measuring
phase synchrony are applied to EEG signals, and the results are analyzed. Images

in this thesis are presented in color.

5.1 Gamma Band Synchrony in Schizophrenic Subjects

One popular theory to the cause of schizophrenia is the reduced coordination of
neural activities in the brain. This theory of schizophrenia is referred as the discon-
nection hypothesis by Friston [36]. Numerous postmortem studies have supported
the disconnect hypothesis in schizophrenic subjects [37, 38]. Further research has
Shown that this theory of reduced coordination causes a lack of functional integration
in the brain [39].

AS previously mentioned in Section 2.2, functional integration in the brain is
indicated by the amount of phase synchrony between neuronal populations. Fur-
thermore, it has been found that oscillations in the 7 band (30-55 Hz) are responsible
for the coordination of neural activity [40, 41]. A schizophrenic subject with low
levels of functional integration Should have low levels of phase synchrony in the 7
band. This hypothesis was shown to be true through experimentation by Spencer
et al. [20] and Uhlhaas et a1. [42] using a measure similar to PLV applied to the
EEG data of control and schizophrenic subjects. The PLV measures used in these
studies were based on the complex Morlet wavelet transform.

In this chapter, we test the hypothesis that schizophrenic patients have lower
7 band synchrony for four schizophrenic and four non—psychiatric control subjects
using the measures defined in Chapters 3 and 4. We have a total of 329 EEG trials

for the four schizophrenia. subjects and a total of 344 trials for the four control

36

subjects.

5.2 Bivariate Phase Synchrony in Schizophrenic Subjects

The 7 band synchrony over 27 different electrodes is compared for the schizophrenic
and control subjects during a cognitive task. The PLV is computed using Equation
3.11 over the P300 window (200—600 ms after the stimulus) and the 7 band over all
trials. Figure 5.1 and Figure 5.2 for the control and schizophrenic subjects Show the
PLV for all electrode pairs averaged over subjects, trials, and the time-frequency
window, respectively. The electrode pairs with high PLV for both the control and
schizophrenic subjects include C3-C3P, C4-C4P, P3-Pz, and P4—Pz. All of these
electrode pairs are located in the central and parietal regions of the brain. The
PLV of these electrode pairs are large likely because of the short physical distances
between the electrodes.

To give a better insight into the results, the difference of the average PLV between
the control and schizophrenic subjects is computed and Shown in Figure 5.3. The
difference between the two subject groups is highest for the electrode pairs C4-P8
and P4-P8. The control subjects have higher PLV for the electrode pairs in the
central and parietal regions of the brain than the schizophrenic subjects. Spencer
et. al also mentioned an increase of the phase synchrony in the parietal electrodes for
the control subjects compared to the schizophrenic subjects in [20]. The increased
PLV for the electrodes in the parietal region are likely caused by the stimulus being
visual. However, for a few of the electrode pairs the schizophrenic subjects showed
increased synchrony, specifically for the pairs: T8-FT8, F7—FT7, and T7-FT7, which
are all located in the frontal and temporal regions of the brain.

The PLV over the time and frequency region with respect. to electrode F 2 aver-
aged over subjects, trials, and the remaining 26 electrodes is shown for control and

schizophrenic subjects in Figure 5.4 and Figure 5.5, respectively. The difference of

37

Average PLV of Control Subjects for each Electrode Pair

electrodes

 

 

FP1
FP2

thNNFN mmvm mvll@9nw wee
uuuomoophmuo daggfigaa EEK

electrodes

Figure 5.1. Average PLV between electrode pairs in the 7 band of the P300 time
window for the four control subjects. The diagonal of the figure, which is the PLV
of an electrode with itself, is set to zero rather than one to improve scaling.

these two figures can be seen in Figure 5.6. The PLV of the schizophrenic subjects
peaks at about 200 milliseconds after the stimulus in the frequency range 35 to 60
Hz. However, the PLV of the control subjects has higher peaks later in the P300
window. Also, the PLV for the control subjects appears highest in the 35 to 45 Hz
frequency range. The higher PLV values for the control subjects become apparent
when looking at the end of the P300 window in the diflerence figure (Figure 5.6)

indicating increased 7 band activity in the P300 window.

38

Average PLV of Schizophrenic Subjects for each Electrode Pair

 

. i. 0.6

electrodes

 

 

   

electrodes

Figure 5.2. Average PLV between electrode pairs in the 7 band of the P300 time
window for the four schizophrenic subjects. The diagonal of the figure, which is the
PLV of an electrode with itself, is set to zero rather than one to improve scaling.

A topo—map of the average PLV with respect to electrode Fz averaged over sub-
jects, trials, and the time-frequency region (P300 time window, '7 frequency band)
is shown for control and schizophrenic subjects in Figures 5.7 and 5.8, respectively.
The highest values of both topo—maps are found at electrodes F3, F4, and Oz. The
difference of these two maps can be seen in Figure 5.9, which makes the results more
apparent. AS previously mentioned during the discussion of Figure 5.3, the control

subjects have higher PLV in the central and parietal regions of the brain, but the

39

Difference of Average PLV for each Electrode Pair

0.05

electrodes

-0.05

 

 

  

electrodes

Figure 5.3. The Average PLV difference between control and schizophrenic subjects.

schizophrenic subjects Show increased synchrony in the temporal and frontal regions
of the brain. The difference topo—map matches this paradigm as well.

The statistical significance of the PLV measure can be shown by comparing it to
the PLV measure using trial-shifted data [28]. The PLV for this surrogate data set
is computed like the original PLV, but the phase difference is between two signals

from a random permutation of the trials rather than two signals from the same trial.

40

Average PLV for Control Subjects

frequency (Hz)

 

 

0.1 0.2 0.3 0.4 0.5 0.6 0.7
time after stimulus (sec)

Figure 5.4. Average PLV with respect to electrode F2 for the four control subjects.

This computation can be seen in the following equation :

erm(k
PLVsurrogate( (t, w) =H 12% 2 exp(] (¢1( (t, w) ——¢p (5.1)

where the constant H is the number of different random permutations. The
PLVsurrogateftvw) is computed 200 times with a H of 40, which results in 200
different phase synchrony surfaces. The maximum of each one of these surfaces is
taken to get 200 surrogate data thresholds. These 200 thresholds are compared to
each time—frequency point’s PLV value from the non-surrogate data set. The PLV

measure is considered statistically significant if it is on average higher than 95% of

41

Average PLV for Schizophrenic Subjects

60

A A 01 01
0 0| O 01

frequency (Hz)

w
01

 

 

30

25

0.1 0.2 0.3 0.4 0.5 0.6 0.7
time after stimulus (sec)

  

Figure 5.5. Average PLV with respect to electrode Fz for the four schizophrenic
subjects.

the surrogate thresholds. Using this measure of statistical significance, it was found
that the PLV measure is significant for both the control and schizophrenic data sets
at all time-frequency points in the P300 time window and 7 frequency band.
While the method described above shows that the PLV measure is significant,
it does not show whether there is significant difference between the control and
schizophrenic subjects based on the phase synchrony measure. The Student’s t—test
can be used to compare the significance of the PLV measure between control and
schizophrenic subjects. The PLV was computed with respect to electrode F2 for
each of the 26 remaining electrodes and each of the subjects averaged over the time

range 500 to 700 ms after the stimulus and the 35-45 Hz frequency range. This

42

Difference of Average PLV

60

A 18 CI 01
O 0| 0 0'!

frequency (Hz)

00
U1

 

30

25

0.1 0.2 0.3 0.4 0.5 0.6 0.7
time after stimulus (sec)

 

Figure 5.6. Difference of the average PLV with respect to electrode Fz between the
control and schizophrenic subjects.

results in the control and schizophrenic subjects each having 26 different PLVS, one
PLV for each electrode pairing. The t-statistic is computed using the samples, and
the significance level is found by looking up the respective t-statistic in the t-table.
Using this measure of statistical significance, it was found that the PLV for the
control subjects is Significantly higher than the PLV for the schizophrenic subjects
at the a = 0.1 significance level. This corresponds to a late response time in the

lower 7 band.

43

Average PLV for Control Subjects

A

 

 

 

Figure 5.7. Topographic map of the average PLV with respect to F2 for four control
subjects.

5.3 Multivariate Phase Synchrony in Schizophrenic Subjects

The PLV method used in the previous section is bivariate and, therefore, can mea-
sure the phase synchrony between only two electrodes. On the other hand, a multi—
variate method can measure the phase synchrony between more than two electrodes
and will give a better indication of ensemble phase synchronization among the elec—
trodes. The multivariate phase synchrony measure proposed in Chapter 4 is applied
to 6 of the 27 electrodes. The six chosen electrodes (C3, C4, P3, P4, P7, and P8 )

are all from the central and parietal regions of the brain and were chosen based on

44

Average PLV for Schziophrenic Subjects

 

 

 

Figure 5.8. Topographic map of the average PLV with respect to Fz for four
schizophrenic subjects.

their increased phase synchrony indicated by the bivariate results of Figure 5.3.
The IFHS for the multivariate method are computed with 6 = 10 and /\ = 0.05
over the P300 window and the 7 frequency band for the six chosen electrodes in
each of the trials for all subjects. Multiple K-means clustering with K = 40 and
P = 500 is then performed on the IFHS from the trials of each subject. This results
in 8 different clusterings one for each subject, 4 control subjects and 4 schizophrenic
subjects. Figures 5.9 through 5.13 display the four most significant cluster centroids
for each control subject, and Figures 5.14 through 5.17 display the four most signif-

icant cluster centroids for each schizophrenic subject. The significance of a cluster

45

Difference of Average PLV

A

0.1

. 0.05

-0.05

 

 

 

Figure 5.9. Topographic map of the difference of the average PLV with respect to
F2 between control and schizophrenic subjects.

is determined by the number of IFHS in the cluster. The magnitude of the cluster
centroids indicates the number of electrodes that were phase synchronous at that
time—frequency point.

The top centroid of control subject 1, Figure 5.10(a), and control subject 2,
Figure 5.11(a), contain most of their phase synchrony in the 40—50 Hz band with
control subject 2 having four phase synchronous electrodes at 45 Hz. On the other
hand, the top centroid of control subject 3, Figure 5.12(a), and control subject
4, Figure 5.13(a). contain most of their phase synchrony in the 50—55 Hz band.

However, the second centroids of control subject 3, Figure 5.12(b), and control

46

subject 4, Figure 5.13(b), contain most. of their phase synchrony in the 40-50 Hz
band.

The top centroids of schizophrenic subject 2, Figure 5.15(a), and schizophrenic
subject 4, Figure 5.17(a), contain relatively little phase synchrony in the 7 frequency
band. Schizophrenic subject 1 has some phase synchrony between two electrodes
in the 50-55 Hz frequency band of the top centroid, Figure 5.14(a). Schizophrenic
subject 3 is a notable exception among the schizophrenic subjects. This subject
contains phase synchrony between three electrodes at the 45-55 Hz frequency band
of its top centroid, Figure 5.16(a).

Tables 5.1 and 5.2 display the correlation averages (Pavg), computed using Equa-
tion 5.2, in five different frequency bands for the four most significant centroids of
the control and schizophrenic subjects, respectively. Table 5.3 shows the weighted
average of the correlation averages, W, from each subject’s four most significant
centroids where the centroid’s Significance is used as the weight. The following

equation shows how Pavg was calculated for each subject:

N “1195,52;

Z:

 

(5.2)

1

C

2: NW)
i=1

(it)

where C is the number of clusters averaged over, paég is the correlation average for
the k-th cluster, and N (k) is the number of IFHS in the k-th cluster. An average of
the W is then taken over control and schizophrenic subjects, and the difference
is taken between the control and schizophrenic subjects which can be seen at the
bottom of Table 5.3.

From the difference row in Table 5.3 it can be seen that. the control subjects

have a higher correlation average than the schizophrenic subjects in the 40-50 Hz

47

frequency band. The higher correlation averages indicate that the control subjects
Show increased phase synchrony over the schizophrenic subjects in the six chosen
electrodes in the P300 time window of the 40-50 Hz frequency band. However,
the schizophrenic subjects do have a higher correlation average for the 50-55 Hz
frequency band.

This result is similar to the results from the bivariate measure as seen in Figure
5.6. With the bivariate measure, the control subjects exhibited increased synchrony
over the schizophrenic subjects in the 35-45 Hz frequency band rather than in the
40-50 Hz frequency band for the multivariate measure. However, the measures did
not use the same electrode pairings. The bivariate measure in Figure 5.6 paired
electrode Fz with the 26 other electrodes while the multivariate measure used the
electrodes: C3, C4, P3, P4, P7, and P8. Despite this difference, the results are

relatively similar.

48

   
     

$3 1?
SE 5
a a
C C
a) d)
3 3
0' O'
2 2
LL u.
30 30
200 400 600 200
Time (ms)
(C)

55 55
1e 50 Tu‘ 50
5 a
5 45 a 45
5 5
g 40 (:7 40
93 E
u- 35 L 35

30 30

200 400 600 200 400 600

Time (ms) Time (ms)

Figure 5.10. The four most significant cluster centroids of control subject 1. The
number of IFHS belonging to each cluster: (a) 15. (b) 14, (c) 3, (d) 3.

49

   

 

55 6
1? ’N‘ 50
5 5 4
a a 45
s s
a a 40 2
9.’ 9
u. LL 35
30 30 o
200 400 600 200 400
Time (ms) Time d(ms)
(C)
55 6 55 6
E 50 E: 50
5‘ 45 5 45
5 s
40 40
S 2 S 2
LL 35 ] LL 35
30 0 30 0
200 400 600 200 400
Time (ms) Time (ms)

Figure 5.11. The four most significant cluster centroids of control subject 2. The
number of IFHS belonging to each cluster: (a) 12, (b) 6, (c) 5. (d) 4.

50

01
U1

  
    

 

if 50 ’5?
3.5, 5
5‘ 45 5‘
8 a“:
3 40 3
E E
U- 35 LL
0
200
55
ET 50 ’57 50
55 E
5 45 5‘ 45
5 g
8 4o] 2", 40
u“. 35’ U': 35
30 30
200 400 200 400 600

 

Time (ms) Time (ms)

Figure 5.12. The four most significant cluster centroids of control subject 3. The
number of IFHS belonging to each cluster: (a) 21, (b) 6, (c) 5, (d) 4.

51

   

   

1? 131‘
5 :5 4
5‘ 5‘
C C
d) d)
3 3
E E 2
u. u.
30 30 o
200 400 600 200 400 600
Time (ms) Time (ms)
(C) (d)
6 55
1? T»? 50
5 4 5
5 3 45
s s
a 2 a 40
93 2
u. U- 35
i
30 o 30
200 400 600 200 400 600
Time (ms) Time (ms)

Figure 5.13. The four most significant cluster centroids of control subject 4. The
number of IF Hs belonging to each cluster: (a) 60, (b) 7, (c) 3, (d) 3.

52

   

   

55 55 6
’N‘ 50 ’N‘ 50
3:, 5 4
5 45 5‘ 45
5 5
g 40 g 40 2
2 9
LL 35 U- 35

30 30 o

200 400 600 200 400
Time (ms) Time d(ms)
(C)

55 55 6
’N‘ 50 1? 50
E 5 4
a 45 3 45
s s
g, 40 g 40 2
E 2
L 35 LL 35

30 30 o

200 400 600 200 400
Time (ms) Time (ms)

Figure 5.14. The four most significant cluster centroids of schizophrenic subject 1.
The number of IFHS belonging to each cluster: (a) 27, (b) 21, (c) 5, (d) 3.

53

   

     

55 55
’N‘ 50 1:: 50
e a
5 45 f; 45
5 8
g 40 g, 40
9 92
u. 35 U- 35

30 30

200 400 600 200 400
Time (ms) Time d(ms)
(C)

55 6 55
1? 50 RT 50
5 4 5
a 45 5 45
s s
g 40 2 a 40
“.3 9
LL 35 LL 35

30 o 30

200 400 600 200 400
Time (ms) Time (ms)

Figure 5.15. The four most significant cluster centroids of schizophrenic subject 2.
The number of IFHS belonging to each cluster: (a) 32, (b) 8, (c) 2, (d) 1.

54

   
 

     

 

55 55
’N‘ 50 if 50
35, E,
a 45 5" 45
5 5
g 40 g 40
2 9
U- 35 IL 35

30 30

200 400 600 200 400
Time (ms) Time d(ms)
(C)

55
1: 50 1? 50
23, E,
g- 45 a 45
s s
g, 40 g 40
9 2
u. 35 LL 35

30

200 400 600 400

 

Time (ms) Time (ms)

Figure 5.16. The four most significant cluster centroids of schizophrenic subject 3.
The number of IFHS belonging to each cluster: (a) 5, (b) 3, (c) 3, (d) 2.

55

   

   

55 6 55
7: 50 RT 50
3:, 4 5
3 45 g 45
5 5
g 40 2 a 40
2 2
LI. 35 U- 35

30 o 30

200 400 600 200 400 600
Time (ms) Time (ms)
(C) (d)

55 55
1: 50 T; 50
E. E
5 45 a 45
5 5
g 40 g 40
2 93
U- 35 U- 35

30 30

200 400 600 200 400 600

   

Time (ms) Time (ms)

Figure 5.17. The four most significant cluster centroids of schizophrenic subject 4.
The number of IFHS belonging to each cluster: (a) 40, (b) 9, (c) 5, (d) 5.

56

Table 5.1. The correlation averages in five different. frequency bands for the four
most significant centroids of each control subject.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Subject Significance Pavg by Frequency Band
30-35 Hz 35-40 Hz 40-45 Hz 45-50 Hz 50—55 Hz
Control 1 15/35 0.0031 0.0002 0.0028 0.0028 0.0000
14/35 0.0003 0.0010 0.0106 0.0000 0.0000
3/35 0.0034 0.0103 0.0201 0.0036 0.0007
3/35 0.0095 0.0093 0.0021 0.0147 0.0002
Control 2 12/ 27 0.0036 0.0000 0.0039 0.0212 0.0000
6/27 0.0096 0.0072 0.0175 0.0034 0.0028
5/27 0.0000 0.0052 0.0021 0.0199 0.0330
4/ 27 0.0033 0.0225 0.0201 0.0000 0.0000
Control 3 21 / 36 0.0000 0.0000 0.0000 0.0003 0.0100
6 / 36 0.0020 0.0039 0.0025 0.0198 0.0000
5 / 36 0.0002 0.0036 0.0095 0.0007 0.0000
4/36 0.0150 0.0008 0.0070 0.0013 0.0000
Control 4 60/73 0.0000 0.0000 0.0000 0.0000 0.0031
7/73 0.0000 0.0062 0.0162 0.0015 0.0000
3/73 0.0062 0.0060 0.0176 0.0078 0.0015
3/ 73 0.0077 0.0134 0.0049 0.0000 0.0000

 

57

 

Table 5.2. The correlation averages in five different frequency bands for the four
most significant centroids of each schizophrenic subject.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Subject Signif. 190129 by Frequency Band
30-35 Hz 35-40 Hz 40—45 Hz 45-50 Hz 50-55 Hz
Schizophrenic 1 27/56 0.0000 0.0000 0.0000 0.0000 0.0018
21/ 56 0.0000 0.0000 0.0000 0.0000 0.0062
5 / 56 0.0069 0.0000 0.0002 0.0031 0.0265
3 / 56 0.0074 0.0047 0.0095 0.0000 0.0003
Schizophrenic 2 32/ 43 0.0000 0.0000 0.0007 0.0005 0.0000
8/43 0.0041 0.0015 0.0051 0.0039 0.0147
2/43 0.0059 0.0263 0.0008 0.0131 0.0211
1/43 0.0312 0.0284 0.0242 0.0337 0.0335
Schizophrenic 3 5/ 13 0.0029 0.0051 0.0038 0.0036 0.0157
3/13 0.0016 0.0015 0.0160 0.0000 0.0000
3/13 0.0078 0.0067 0.0078 0.0010 0.0175
2/13 0.0105 0.0391 0.0075 0.0034 0.0013
Schizophrenic 4 40/ 59 0.0010 0.0000 0.0000 0.0000 0.0007
9/59 0.0000 0.0002 0.0008 0.0018 0.0088
5/59 0.0011 0.0013 0.0020 0.0000 0.0069
5 / 59 0.0038 0.0095 0.0051 0.0052 0.0000

 

 

 

Table 5.3. The weighted average of the correlation averages for the four most signif-
icant centroids of each control and schizophrenic subject in five different frequency
bands. The final row contains the difference between the average control subject’s
correlation average and the average schizophrenic subject’s correlation average.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Subject Weighted Avg Pang by Frequency Band
30-35 Hz 35—40 Hz 40-45 Hz 45-50 Hz 50-55 Hz
Control 1 0.0026 0.0022 0.0073 0.0028 0.0001
Control 2 0.0042 0.0059 0.0090 0.0139 0.0067
Control 3 0.0020 0.0012 0.0025 0.0037 0.0058
Control 4 0.0006 0.0014 0.0025 0.0005 0.0026
Schizophrenic 1 0.0010 0.0003 0.0005 0.0003 0.0056
Schizophrenic 2 0.0018 0.0022 0.0021 0.0025 0.0045
Schizophrenic 3 0.0049 0.0099 0.0081 0.0021 0.0103
Schizophrenic 4 0.0011 0.0009 0.0007 0.0007 0.0024
Control Avg. 0.0024 0.0027 0.0053 0.0052 0.0038
Schizophrenic Avg. 0.0022 0.0033 0.0029 0.0014 0.0057
Difference 0.0002 -0.0006 0.0024 0.0038 —0.0019

 

58

 

CHAPTER 6

CONCLUSIONS AND FUTURE WORK

In this thesis, two different measures of phase synchrony, bivariate and multivari-
ate measures, were introduced. The bivariate phase synchrony measure computes
the phase locking value (PLV) using the reduced interference Rihaczek distribution
between pairs of signals and gives a measure of the variability of phase differences
across trials. The multivariate phase synchrony measure, on the other hand, is a
deterministic measure of phase synchrony that is directly related to the instanta-
neous frequency. Multiple K-means clusterings are performed on the instantaneous
frequency histograms in the time—frequency plane to estimate the phase synchrony

across multiple electrodes.

6. 1 Contributions

One of the major contributions of this thesis is the introduction of the reduced-
interference Rihaczek distribution. Time-frequency distributions like the Wigner
and Rihaczek distributions suffer from cross-terms when applied to multi-component
signals. The reduced interference Rihaczek distribution removes the cross-terms by
applying a Choi-Williams kernel in the ambiguity domain. In addition, the reduced
interference Rihaczek distribution, like the Rihaczek distribution, is complex and
retains the phase information unlike classical time-frequency methods such as the
Morlet wavelet and Wigner distribution.

While previous bivariate measures computed the PLV using a different time-
frequency distribution such as the Morlet wavelet, this thesis proposed using the
reduced interference Rihaczek distribution to compute the PLV. The Morlet wavelet
suffers from biased phase estimates because of its nonuniform frequency resolution.

The reduced interference Rihaczek distribution, like all time-frequency distributions

 

in Cohen’s class, has a uniform resolution and as a result has an unbiased phase
estimate. The introduction of Cohen’s class of time-frequency distributions to the
computation of PLV is a major contribution of this thesis.

Another contribution is the proposed multivariate phase synchrony measure,
which can be applied to multiple trials, unlike previous multivariate measures. Pre-
vious multivariate measures could only be applied to a single trial, and as a result
every instantaneous frequency histogram (IFH) would have to be analyzed indepen-
dently to get an overview of the phase synchrony across the trials. By applying
K-means clustering to the IFHS, the major frequency locking patterns can be ex-
tracted and analyzed to give a general overview of the phase synchrony across the

trials and electrodes.

6.2 Future Work

6.2.1 Multivariate Method Improvements

One disadvantage of the proposed multivariate method is the need for parameter
calibration. The IF estimation requires two different parameters, the time support
threshold (6) and the average energy threshold (A). These parameters determine
which components are filtered out of the IF estimate. The time support threshold
depends on the length of the signal, and the average energy threshold depends on
the energy of the signal. Multiple K-means requires two additional parameters, the
number of clusters (K ) and number of times K-means clustering is performed (P).
A larger P value gives a more accurate clustering. Both of these parameters depend
on the number of trials and time-frequency points. An automated technique to
optimize these parameters would significantly improve the efficiency of this method.

Furthermore, another disadvantage of the proposed multivariate method is the
inability to determine which electrodes are phase synchronous. Only the number

of electrodes that are phase synchronous are known from the IFH. This issue could

60

be solved by multiplying each electrode’s IF estimate by a unique number that is a
power of two. The resulting IFH can then be decoded to determine the original IF
estimates from each electrode.

In addition, the multivariate method cannot be easily extended to multiple sub-
jects. The cluster centroids cannot be averaged together because of their discrete
finite values, and using K-means clustering on the IFHS from multiple subjects

proved to be ineffective.
6.2.2 Volume Conduction

Volume conduction, which was first defined in Section 2.1, is the mixing of brain
signals traveling from the brain to the electrodes on the scalp. A reduction of volume
conduction would assure that only true phase synchrony is measured. The effects
of volume conduction can be seen in the topo-maps of Figures 5.7 and 5.8. The
electrodes around the reference electrode (Fz) have higher PLVs than the electrodes
further away from the electrode.

Lachaux suggests two different techniques to reduce volume conduction, scalp
current density and inverse deblurring [31]. The scalp current density is computed
by multiplying the Laplacian of the EEG by the negative of the scalp conductivity
[43]. Measures of phase synchrony can then be applied directly to the scalp current
density. Inverse deblurring estimates the potentials of the cortical sources in a finite
element model by optimizing the difference between the measured EEG data and the
estimated potentials of cortical sources conducted through the skull and scalp [44].
Measures of phase synchrony can then be applied to the estimated potentials of the
cortical sources. Another popular technique to reduce volume conduction is baseline
correction, which subtracts the phase synchrony measure during the pre—stimulus
from the phase synchrony measure during the P300 time window. The measure from
the pre—stimulus should be mostly volume conduction, and by subtracting it. from

the measure during the P300 time window, the phase synchrony caused by volume

61

conduction will be removed.
6.2.3 Directional Measures

Neither of the measures presented in this thesis are directional. From directionality,
the causality or the direction of the coupling could be determined. Pereda et a1.
[45] discuss several directionality measures such as transfer entropy, partial directed
coherence, and Granger causality. However, none of these are measures of phase
synchrony because they take the amplitude of the signals into account. The proposed
phase synchrony methods could be extended to measure the directionality of the
phase synchrony either by using the directionality index proposed by Rosenblum
and Pikovsky [46] or by computing the phase synchrony measures using time-shifted

signals.
6.2.4 Functional Networks

The proposed phase synchrony measures can be extended to determine the func-
tional network patterns in the brain across space, time, and frequency. The mea—
sures will be the first step in defining adjacency matrices which can be used in graph
theoretic analysis of functional networks to determine the organization patterns of
the brain such as small world network models. Furthermore, if a directional mea-
sure is used, a directional network could be created. Once a network is constructed,
methods in graph theory could be applied for analysis. Stam et al. converted syn-
chronization matrices similar to the ones seen in Figures 5.1 and 5.2 into graphs
and applied methods from graph theory such as the clustering coefficient and path

length[47].

62

BIBLIOGRAPHY

63

BIBLIOGRAPHY

[1] R.S.J. Frackowiak, K. Friston, C. Frith, R. Dolan, C.J. Price, S. Zeki, J. Ash-
burner, and W. Penny, Human Brain Function, Academic Press, 2003.

[2] 0. David, D. Cosmelli, J.P. Lachaux, S. Baillet, L. Garnero, and J. Martinerie,
“A theoretical and experimental introduction to the non-invasive study of large-
scale neural phase synchronization in human beings,” International Journal of
Computational Cognition, vol. 1, pp. 53—77, 2002.

[3] E. Rodriguez, N. George, J.P. Lachaux, J. Martinerie, B. Renault, and F.J.
Varela, “Perception’s shadow: long-distance synchronization of human brain
activity,” Nature, vol. 397, no. 6718, pp. 430—433, 1999.

[4] D. Rudrauf, A. Douiri, C. Kovach, J.P. Lachaux, D. Cosmelli, M. Chavez,
C. Adam, B. Renault, J. Martinerie, and M. Le Van Quyen, “Frequency flows
and the time-frequency dynamics of multivariate phase synchronization in brain
signals,” Neuroimage, vol. 31, no. 1, pp. 209—227, 2006.

[5] E. Basar, EEG-Brain Dynamics: Relation Between EEG and Brain Evoked
Potentials, Elsevier-North—Holland Biomedical Press, 1980.

[6] B.J. Fisch, F isch and Spehlmann’s EEEG Primer: Basic Principles of Digital
and Analog EEC, Elsevier Science Health Science, 1999.

[7] G. Tononi, G.M. Edelman, and O. Sporns, “Complexity and coherency: inte-
grating information in the brain,” Trends in Cognitive Sciences, vol. 2, no. 12,
pp. 474—484, 1998.

[8] F. Varela, J .P. Lachaux, E. Rodriguez, and J. Martinerie, “The brainweb: phase
synchronization and large-scale integration,” Nature Reviews Neuroscience, vol.
2, no. 4, pp. 229—239, 2001.

[9] A. Von Stein, P. Rappelsberger, J. Sarnthein, and H. Petsche, “Synchronization
between temporal and parietal cortex during multimodal object processing in
man,” Cerebral Cortex, vol. 9, no. 2, pp. 137—150, 1999.

[10] V. Samar, K. Swartz, and M. Raghuveer, “Multiresolution analysis of event—
related potentials by wavelet decomposition,” Brain and Cognition, vol. 27, no.
3, pp. 398—438, 1995.

[11] J. Raz, L. Dickerson, and B. Turetsky, “A wavelet packet model of evoked
potentials,” Brain and Language, vol. 66, no. 1, pp. 61—88, 1999.

[12] R. Quiroga, O.W. Sakowitz, E. Basar, and M. Schiirmann, “Wavelet trans-

64

[13]

[14]

[15]

[15]

[17]

[18]

[19]

[20]

[21]

[22]

form in the analysis of the frequency composition of evoked potentials,” Brain
Research Protocols, vol. 8, no. 1, pp. 16—24, 2001.

P. Zarjam, G. Azemi, M. Mesbah, and B. Boashash, “Detection of newborns’
EEG seizure using time-frequency divergence measures,” in Proceedings of

the IEEE International Conference of Acoustics, Speech and Signal Processing,
2004, vol. 5, pp. 429—432.

M. Sun, M.L. Scheuer, S. Qian, S.B. Baumann, P.D. Adelson, and R.J.
Sclabassi, “Time-frequency analysis of high-frequency activity at the start
of epileptic seizures,” in Proceedings of the IEEE International Conference on
Engineering in Medicine and Biology, 1997, vol. 3, pp. 1184—1187.

M. Carlos Guerrero, A.M. Trigueros, and J .I. Franco, “Time-frequency EEG
analysis in epilepsy: What is more suitable?,” in Proceedings of the IEEE

International Symposium on Signal Processing and Information Technology,
2005, pp. 202—207.

P. Zarjam and M. Mesbah, “Discrete wavelet transform based seizure detection
in newborns EEG signals,” in Proceedings of the Seventh International Sympo-
sium on Signal Processing and Its Applications, 2003, vol. 2, pp. 459—462.

P. Zarjam, M. Mesbah, and B. Boashash, “Detection of newborn EEG seizure
using optimal features based on discrete wavelet transform,” in Proceedings of

the IEEE International Conference on Acoustics, Speech and Signal Processing,
2003, vol. 2, pp. 265—268.

R. Polikar, M.H. Greer, L. Udpa, and F. Keinert, “Multiresolution wavelet
analysis of ERPs for the detection of Alzheimer’s disease,” in Proceddings of

the IEEE 19th International Conference of Engineering in Medicine and Biology
Society, 1997, pp. 1301—1304.

G. Jacques, J.L. Frymiare, J. Kounios, C. Clark, and R. Polikar, “Multires-
olution analysis for early diagnosis of Alzheimer’s disease,” in Proceddings of
26th Annual International Conference of IEEE Engineering in Medicine and
Biology Society, 2004, vol. 1, pp. 251—254.

K.M. Spencer, P.G. Nestor, M.A. Niznikiewicz, D.F. Salisbury, M.E. Shenton,
and R.W. McCarley, “Abnormal neural synchrony in schizophrenia,” Journal
of Neuroscience, vol. 23, no. 19, pp. 7407—7411, 2003.

N. Hazarika, J.Z. Chen, A.C. Tsoi, and A. Sergejew, “Classification of EEG
signals using the wavelet transform,” Signal Processing, vol. 59, no. 1, pp.
61-i-72, 1997.

L. Cohen, Time-frequency analysis: theory and applications, Prentice Hall,

[23]

[24]

[25]

[26]

[27]

[28]

[29]

[30]

[31]

[32]

[33]

[34]

1995.

A. Rihaczek, “Signal energy distribution in time and frequency,” IEEE Trans-
actions on Information Theory, vol. 14, no. 3, pp. 369—374, 1968.

J.P. Lachaux, A. Lutz, D. Rudrauf, D. Cosmelli, M. Le Van Quyen, J. Mar-
tinerie, and F. Varela, “Estimating the time-course of coherence between single-
trial brain signals: an introduction to wavelet coherence,” Neurophysiologie
Clinique/ Clinical Neurophysiology, vol. 32, no. 3, pp. 157—174, 2002.

A. Ambardar, Analog and digital signal processing, PWS, 2002.

S. Haykin, R.J. Racine, Y. Xu, and CA. Chapman, “Monitoring neuronal os-
cillations and signal transmission between cortical regions using time-frequency
analysis of electroencephalographic activity,” Proceedings of the IEEE, vol. 84,
no. 9, pp. 1295—1301, 1996.

M. Le Van Quyen, J. Foucher, J.P. Lachaux, E. Rodriguez, A. Lutz, J. Mar-
tinerie, and F.J. Varela, “Comparison of Hilbert transform and wavelet methods

for the analysis of neuronal synchrony,” Journal of Neuroscience Methods, vol.
111, no. 2, pp. 83—98, 2001.

J.P. Lachaux, E. Rodriguez, M. Le Van Quyen, A. Lutz, J. Martinerie, and
F.J. Varela, “Studying single-trials of phase synchronous activity in the brain,”
International Journal of Bifurcation and Chaos, vol. 10, pp. 2429—2439, 2000.

D. Li and R. Jung, “Quantifying coevolution of nonstationary biomedical sig-
nals using time-varying phase spectra,” Annals of Biomedical Engineering, vol.
28, no. 9, pp. 1101—1115, 2000.

S. Aviyente, W.S. Evans, E.M. Bernat, and S. Sponheim, “A time-varying
phase coherence measure for quantifying functional integration in the brain,”
in Proceddings of the IEEE International Conference on Acoustics, Speech and
Signal Processing, 2007, vol. 4, pp. 1169—1172.

J .P. Lachaux, E. Rodriguez, J. Martinerie, and F.J. Varela, “Measuring phase
synchrony in brain signals,” Human Brain Mapping, vol. 8, pp. 194—208, 1999.

C. Allefeld and J. Kurths, “An approach to multivariate phase synchronization
analysis and its application to event-related potentials,” International Journal
of Bifurcation and Chaos, vol. 14, no. 2, pp. 417—426, 2004.

C. Allefeld and J. Kurths, “Multivariate phase synchronization analysis of EEG
data,” IE1 CE Transactions on Fundamentals of Electronics, Communications,
and Computer Sciences, vol. 86, no. 9, pp. 2218—2221, 2003.

B. Boashash, “Estimating and interpreting the instantaneous frequency of a.

66

 

[38]

[39]

[40]

[41]

[42]

[43]

[44]

[4

l

signal—part 2: algorithms and applications,” Proceedings of the IEEE, vol. 80,
no. 4, pp. 540—568, 1992.

L. Rankine, M. Mesbah, and B. Boashash, “IF estimation for multicompo-
nent signals using image processing techniques in the time—frequency domain,”
Signal Processing, vol. 87, no. 6, pp. 1234—1250, 2007.

K.J. Friston, “Schizophrenia and the disconnection hypothesis,” Acta Psychi-
atrica Scandinavica, vol. 99, pp. 68—79, 1999.

L.D. Selemon and PS. Goldman-Rakic, “The reduced neuropil hypothesis: a
circuit based model of schizophrenia,” Biological Psychiatry, vol. 45, no. 1, pp.
17—25, 1999.

FM. Benes, “Emerging principles of altered neural circuitry in schizophrenia,”
Brain Research Reviews, vol. 31, no. 2—3, pp. 251—269, 2000.

G. Tononi and GM. Edelman, “Schizophrenia and the mechanisms of conscious
integration,” Brain Research Reviews, vol. 31, no. 2-3, pp. 391—400, 2000.

P. Fries, S. Neuenschwander, A.K. Engel, R. Goebel, and W. Singer, “Rapid
feature selective neuronal synchronization through correlated latency shifting,”
Nature Neuroscience, vol. 4, pp. 194—200, 2001.

J.S. Kwon, B.F. O’Donnell, G.V. Wallenstein, R.W. Greene, Y. Hirayasu,
P.G. Nestor, M.E. Hasselmo, G.F. Potts, M.E. Shenton, and R.W. McCarley,
“Gamma frequency-range abnormalities to auditory stimulation in schizophre-
nia,” Archives of General Psychiatry, vol. 56, no. 11, pp. 1001—1005, 1999.

F.J. Uhlhaas, D.E.J. Linden, W. Singer, C. Haenschel, M. Lindner, K. Maurer,
and E. Rodriguez, “Dysfunctional long-range coordination of neural activity

during gestalt perception in schizophrenia,” Journal of Neuroscience, vol. 26,
no. 31, pp. 8168—8175, 2006.

F. Perrin, O. Bertrand, and J. Pernier, “Scalp Current Density Mapping:
Value and Estimation from Potential Data,” IEEE Transactions on Biomedical
Engineering, pp. 283—288, 1987.

J. Le and A. Gevins, “Method to reduce blur distortion from EEG’s using a
realistic headmodel,” IEEE Transactions on Biomedical Engineering, vol. 40,
no. 6, pp. 517—528, 1993.

E. Pereda, R.Q. Quiroga, and J. Bhattacharya, “Nonlinear multivariate analysis
of neurophysiological signals,” Progress in Neurobiology, vol. 77, no. 1-2, pp.
1-37, 2005.

67

[46] MG. Rosenblum and AS. Pikovsky, “Detecting direction of coupling in inter—
acting oscillators,” Physical Review E, vol. 64, no. 4, pp. 45202, 2001.

[47] C.J. Stam, B.F. Jones, G. Nolte, M. Breakspear, and P. Scheltens, “Small-
world networks and functional connectivity in alzheimer’s disease,” Cerebral
Cortex, vol. 17, no. 1, pp. 92—99, 2007.

68