CC BY
62
Наука и Образование
МГТУ им. Н.Э. Баумана
Сетевое научное издание
Наука и Образование. МГТУ им. Н.Э. Баумана. Электрон. журн. 2014. № 10. С. 114-136.
Б01: 10.7463/1014.0728495
Представлена в редакцию: 27.09.2014
© МГТУ им. Н.Э. Баумана
УДК 519.6
Метод сегментации пространственно-распределенных временных рядов на основе бегущих волн
Трофимов А. Г.1'2'", Колодкин И. В.1, 2 12 Ушаков В. Л. , Величковский Б. М. '
:Россия, Национальный исследовательский ядерный университет
«МИФИ»
2Россия, НИЦ «Курчатовский институт»
ап-ойтоу glist.ru
Предложен метод сегментации многомерных временных рядов, обладающих пространственной организацией, на основе характеристик бегущих волн. Для оценки степени выраженности бегущих волн введены показатели синхронности и когерентности. Расчёт этих показателей проводился на основе результатов анализа кросс-спектров сигналов в пространственно близких точках наблюдения. В результате экспериментальных исследований на реальных данных электроэнцефалографии (ЭЭГ) головного мозга показано, что предложенный метод обеспечивает более качественную сегментацию по сравнению традиционными методами сегментации сигналов ЭЭГ.
Ключевые слова: пространственно-распределенный временной ряд, бегущая волна, сегментация, синхронность, когерентность, микросостояния ЭЭГ, электроэнцефалограмма.
Введение
При наблюдении выходов распределённой динамической системы можно заметить в некоторые моменты времени резкое изменение характеристик наблюдаемых сигналов. Эти моменты времени могут соответствовать смене режима функционирования системы, изменению характеристик внешних воздействий на систему или являться следствием присущей системе нестационарности. Задача обнаружения этих моментов времени и описания наблюдаемого сигнала как последовательности сменяющих друг друга участков, обладающих схожими характеристиками (сегментов), является одной из ключевых задач обработки временных рядов [1].
Выделенные сегменты не только позволят представить временной ряд в более компактной форме, но и лучше понять природу наблюдаемых сигналов. Анализ характеристик выделенных сегментов, их продолжительностей, последовательностей смены одних сегментов другими может предоставить исследователю важную
информацию о наблюдаемом явлении и служить основой для решения других задач обработки временных рядов, в частности, классификации.
Выделение сегментов (сегментация) временного ряда, как правило, предполагает следующую последовательность действий.
1) Предобработка временного ряда. На этом шаге осуществляется фильтрация сигнала от шумов, децимация, устранение артефактов, сокращение размерностей или другие операции, позволяющие на ранних этапах отсеять неинформативные данные, используя априорно известную информацию о природе сигнала.
2) Построение пространства признаков для сегментации. Каждый сегмент временного ряда должен обладать схожими характеристиками. От того, как определены эти характеристики, зависят результат и качество сегментации. Как правило, эти характеристики выбираются исследователем, исходя из целей сегментации, и отражают некоторые определённые аспекты временного ряда. Например, в качестве характеристик временного ряда могут быть выбраны среднее значение или дисперсия, рассчитанные в скользящем окне, спектральные характеристики, рассчитанные в результате частотно-временного анализа [2] или др.
3) Поиск кластеров в выделенном пространстве признаков. Как только определено пространство характерных признаков временного ряда, далее возникает задача определения групп близко расположенных данных в нём. Для решения этой задачи могут быть использованы известные методы кластеризации статических данных, например, метод ^-средних или самообучающиеся нейронные сети [3], а также специальные методы, учитывающие временную организацию данных, например, метод агломеративной кластеризации темпоральных данных [4], темпоральные нейронные сети Кохонена [5] или др.
В настоящей работе рассматриваются многомерные временные ряды, имеющие пространственно-временную структуру, т.е. временные ряды, характеризующиеся в каждый момент времени вектором наблюдений, каждое из которых соответствует некоторой точке пространства. Примерами таких рядов являются сигналы электроэнцефалограмм (ЭЭГ), магнитоэнцефалограмм (МЭГ), сигнал сейсмической активности Земли, результаты наблюдения погоды и т.д. Сегментацию таких рядов целесообразно проводить на основе характеристик, учитывающих распределённую природу наблюдаемых динамических данных. В работе [6] использованы характеристики пространственно распределённых временных рядов, основанные на результатах корреляционного анализа. В [7] рассмотрен метод анализа пространственных распределений значений временных рядов. В [8] приводятся методы моделирования пространственно организованных временных рядов, параметры которых могут быть использованы для их описания. В настоящей работе предложен подход к выделению характеристик пространственно организованных временных рядов, основанный на анализе бегущих волн [9].
Бегущая волна - это волновое движение, при котором поверхность равных фаз (фазовые волновые фронты) перемещается с конечной скоростью, постоянной в случае однородных сред [9]. Бегущие волны наблюдаются во многих физических, химических и биологических процессах, в частности, в реакционно-диффузных системах, процессах изменения численности популяций животных, в распределённой динамике солнечной активности.
Цель данной работы - предложить алгоритм сегментации многомерного временного ряда, обладающего пространственно-временной структурой, и продемонстрировать его работу на реальных данных.
В первом разделе приводится математическая постановка задачи сегментации многомерного временного ряда. Второй раздел посвящён формальному описанию бегущей волны и показателей выраженности бегущей волны между двумя пространственными точками. В третьем разделе рассматривается метод построения пространства признаков для сегментации многомерного временного ряда, основанных на характеристиках бегущих волн. В четвертом разделе приводятся результаты экспериментальных исследований предложенного алгоритма сегментации на реальных данных электроэнцефалографии. В заключении формулируются основные результаты работы и направления дальнейших исследований, связанных с выделением когнитивных микросостояний головного мозга человека и построением его функциональных коннектомов.
1. Постановка задачи
Рассматривается многомерный временной ряд X = ( х(1),..., х(Т)), характеризующийся
в каждый момент времени t вектором наблюдений х(±) =(х(^),..., X (/)), t = 1, Т, где Т -
число временных отсчётов. Каждое наблюдение х ^), ^ = 1, Ь, получено в точке на плоскости с координатами (^ ) в некоторой системе координат, каждая точка
р = (^, щ) , / = 1, Ь, принадлежит односвязной области О. Ряд X получен в результате дискретизации некоторого непрерывного сигнала с постоянной частотой дискретизации Ставится задача выделения микросостояний ^,...,Бк временного ряда Хи отнесения
вектора наблюдений х^) в каждый момент времени t, t = 1, Т, одному из выделенных
микросостояний. Микросостоянием будем считать область в некотором пространстве
признаков, характеризующих временной ряд X в каждый момент времени t, t = 1, Т.
Формально задача выделения микросостояний ^,...,8К - это задача нахождения
оператора Е, сопоставляющего временному ряду X в каждый момент времени t номер к(^ микросостояния:
k(t) = F [X, t], t = 1, T, k(t) e {1,..., K}, (1)
где K - число микросостояний. Отметим, что число микросостояний временного ряда априорно может быть неизвестно.
Введённое понятие микросостояния требует определения пространства признаков, в котором будут выделяться эти микросостояния. Обозначим через ф оператор, сопоставляющий временному ряду X в каждый момент времени t некоторый вектор характерных признаков y(t) = (y (t),..., yM (t)):
y(t) = ф[ X, t], t = \T, y(t)EWM. (2)
В частном случае значение оператора ф в момент времени t может зависеть не от всего временного ряда X, а от его фрагмента во временном окне с центром в t и шириной т:
y(t) = ф[х(t-т),..., x(t+т)], t = т +1, T-т, y(t)EWM. (3)
В другом частном случае вектор y(t) может совпадать с вектором x(t) при всех t, t = 1, T, т.е. т = 0, M = L, ф[x(t)] = x(t). Множество векторов характерных признаков образует многомерный временной ряд, который далее будем обозначать Y = ( y(1),..., y(T)).
В настоящей работе предлагается каждое микросостояние Si характеризовать соответствующим эталоном wt =(wy,...,wMi) в пространстве характерн^1х признаков y,
i = 1, K. Таким образом, задача выделения микросостояний состоит в определении их
эталонных представителей. Как только эталонные представители w1,...,wK найдены,
принадлежность вектора значений временного ряда X тому или иному микросостоянию в момент времени t может быть определена по критерию близости соответствующего вектора характерных признаков y(t) к эталонному представителю микросостояния:
x(t) e Qk ^Р (y(t), Wk ) = min p (y(t), w ), t = IJ, (4)
i=1, K
где p( x, y) - расстояние (в частности, евклидово) между векторами x и y.
Сопоставляя выражения (1), (2) и (4), определим вид рассматриваемых в работе операторов F:
F [ X, t ] = argmin р(ф[ X, t ], w), t = VT. (5)
i=1, K
Таким образом, задача выделения микросостояний (нахождения оператора F) сводится к двум подзадачам:
1) нахождение оператора ф, осуществляющего переход к пространству характерных признаков y;
2) нахождение эталонных представителей w1,...,wK микросостояний в пространстве выделенных признаков.
Для решения каждой из этих подзадач может быть предложено множество способов. Учитывая пространственно-временную организацию данных, в настоящей работе предлагается использовать характеристики бегущих волн, возникающих в области Q.
2. Формальное описание бегущей волны
Согласно [9], возмущение среды z(r, t), создаваемое плоской бегущей волной в момент времени t в точке с координатой r, описывается выражением:
z (r, t) = A(r, t)sin (kr - rat+ф0 ), (6)
где A(r, t) - амплитуда волны в момент времени t в точке с координатой r, k - волновое число, ш - круговая частота, фо - начальная фаза волны.
Предположим, что от точки Pi в направлении точки Pj в области Q распространяется монохроматическая бегущая волна с частотой ш (рис. 1).
а б
Рис. 1. Иллюстрация понятия бегущей волны
Выберем точку P¿ в качестве начала координат, а ось координат зададим в направлении точки Pj. Тогда, учитывая (6), возмущения среды z¿(t) и zj(t) в точках P¿ и Pj соответственно, создаваемые рассматриваемой бегущей волной в момент времени t, будут описываться выражениями:
Z(t) = z (0, t) = Д. (ra, t) sin (-rat+ф0), (7)
zj (t) = z (tj, t) = aj (ra, t) sin (ktj - rat + фо), (8)
где t = P(P,Pj) - расстояние между точками P¡ и Pj, A¿(ш, t) и Д7(ш, t) - амплитуды волны частоты ш в точках P¿ и Pj соответственно.
Отметим, что для однородной стационарной среды волновое число к не зависит ни от времени, ни от пространственной координаты. Разность фазы Дфу волны в один и тот же момент времени в точках P, и Pj постоянна и равна:
ДФУ- = ^ . (9)
Условие (9) является необходимым условием наличия бегущей монохроматической волны между точками Pi и Pj. Если Дф > 0 (т.е. к > 0), то волна распространяется в
направлении от Pi к Pj, если Дф < 0 (т.е. к < 0), то - в направлении от Pj к Pi. Отметим,
что условие (9) не является достаточным, поскольку постоянство разности фаз между точками Pi и Pj может наблюдаться и при распространении двух независимых бегущих волн, проходящих через эти точки в различных направлениях.
На практике бегущая волна между точками Pi и Pj может не являться монохроматической. В этом случае для выделения монохроматических компонентов волны используется спектральное разложение (например, с помощью преобразования Фурье). Кроме того, бегущая волна может проявляться лишь на некотором интервале (или интервалах) времени.
Для выявления наличия бегущей волны между точками Pi и Pj в момент времени t перейдём к спектральному представлению наблюдаемых в этих точках временных рядов
(t), t = 1,T| и (t), t = 1,Tj соответственно. Обозначим через {Хг (ю, t),ra = юх,...,юГ } результат применения к временному ряду |хг (t), t = 1, T j дискретного оконного преобразования Фурье с центром окна в момент времени t, i = 1, L :
T 2 m(т-У) _
X,(®m,t) = £ Xi(t)h(T - t)e T , m = 1,T, (10)
t=1
где h(t) - оконная функция, i в показателе экспоненты - мнимая единица.
f
Мгновенные амплитуды и фазы гармоник на частоте ют = ^¡р (т -1) в момент
времени t, т = 1, T, t = 1, T, рассчитываются на основе результатов преобразования Фурье:
Л, (Ют , t) =\X, (Ют , t)|, i = \L , (11)
Фг (Ют, t) = arg Хг (Ют, t) , i = 1, L . (12)
Учитывая (9), сформулируем необходимое условие наличия бегущей волны частоты &т между точками Pi и Pj на интервале времени ; t2 ]:
ДФу (Ют , t) = Ф, (Ют , t) - Ф; (Ют, t) = к(Ют Уг] , Vt G [t{, t2], (13)
где к(шт) - волновое число на частоте &т, т = 1, T .
Вследствие погрешности измерения исходных временных рядов, неточностей оценки их спектрального состава, возможной нестационарности (квазистационарности) или других причин условие (13) носит теоретических характер. На практике для проверки необходимого условия наличия бегущей волны используем статистические
характеристики выборки разностей фаз в каждый момент времени t е \tx; t2 ]. В настоящей работе предлагается несколько характеристик, чувствительных к условию (13).
Пусть Лф.(rom,t1),...,Лф.(rom,t2) - выборка разностей фаз волны частоты &m в точках
Pi и Pj, i е{1,..., L}, j е{1,..., L}, в моменты времени tx,..., t2. Для оценки степени выраженности условия (13) введём комплексные числа z. (rom, t1),..., z. (rom, t2) с единичным модулем и фазами, равными Лф..(rom,t1),...,Лф..(rom, t2) соответственно: | z.(rom, t) |= 1,
ar§ Z.j (®m ' t) = АФ. (®m > t) , t = tj, t2. Рассчитаем среднее для них:
iz, .t))'
1
t=ti t — t
e2 П t=tj
i Zj , t) .
(14)
Используя представление комплексного числа в экспоненциальном виде, запишем среднее (14) в виде:
zn (fflm > t))^ = Т^7 i eXP [/АФ j (Ют > t)] =
1 ^ tj t=t1
t2 t2 i cos АФУ. (юи,t) + iL sin аф . (юи,t)
t2 t1
(15)
В качестве меры разброса комплексных чисел zij(ют,^),...,zij(ют,¿2) относительно среднего значения выберем величину отклонения комплексного среднего ^ (юот, ^)^2 от единичной окружности (рис. 2). Критерием выполнения условия (13) является равенство единице модуля среднего ^(ют,. Таким образом, определим показатель уг>(юот,Хх,)
степени выраженности условия (13) как модуль комплексного среднего ^(ют,:
Yj К > t1 > t2 ) =
К,t )>
t2 — ^
i2 Y f h ^
ii cos АФУ. (юи,t) + ii sm АФУ. (юи,t)
У V t=t1
Vt=t1
. (16)
Возможные значения показателя (16) лежат в отрезке от 0 до 1.
Показатель (16) характеризует только степень разброса разностей фаз Лф.(rom,tj),...,Лф.(rom,t2) волны на частоте Qm (phase-locking value), при этом никак не
учитывает амплитуды волны Д (ют, ^) и А}- (ют, ^), ^ = ^, ¿2 . Синхронность амплитуд Д (юот, ^) и А. (ют, ¿) не является необходимым условием наличия бегущей волны между
точками Р1 и Pj, однако может говорить о согласованности изменения энергии волны в этих точках. Для оценки степени выраженности бегущей волны частоты шт между точками Р^ и Pj и синхронности изменения её амплитуд в этих точках на интервале времени ; ] может быть использован показатель когерентности [10]:
1
г, (шя>'1,Ч) =
X (шт) X* (ая))
^г К )К ))I¿X, (Шя )X* К ))*=*
(17)
где символ означает комплексное сопряжение.
Рис. 2. Иллюстрация к расчёту показателя у..(шт, ^, '2 )
Отметим, что модуль усреднённого кросс-спектра сигналов ('), ' = 1, Т | и
('), ' = 1, Т| (числитель показателя когерентности) является обобщением показателя (16):
{Хг (Шя ) X* (Шя )) ^
*2 '
"2
X Аг (Шя ' ')А, (Шя > ') COS А% (Шя , ')
V'='1
,Л1/2
2
X Аг (Шя > ') А, (Шя > ') АФ, (Шя > ')
V'='1
Возможные значения когерентности лежат в отрезке от 0 до 1, причём равенство когерентности единице является критерием выполнения условия (13) и согласованности изменения амплитуд бегущей волны в точках Рг и Р, на интервале времени ; '2 ].
2
1
3. Метод построения пространства признаков, основанных на характеристиках бегущих волн
Введённые показатели (16) и (17) выраженности бегущей волны могут быть рассчитаны для любой пары точек Рг и Р, из области О. В настоящей работе предложено рассчитывать эти показатели лишь для соседних точек. Другой способ может быть в расчёте этих показателей относительно некоторой референтной точки.
Обозначим через В = (Ъ,...,Ъ) множество пар соседних точек среди р,...,р . Для
нахождения соседних точек среди р,...,р используем диаграмму Вороного. Этот подход в последнее время начинает применяться в нейрокогнитивных исследованиях сознания и памяти [11]. Диаграмма Вороного представляет собой разбиение области О на Ь областей таких, что каждая область представляет собой геометрическое место точек, более близких к одной из точек р,...,р, чем к любой другой точке из этого множества (рис. 3). При
использовании евклидовой метрики каждая область будет являться выпуклым многоугольником. Для построения диаграммы Вороного может быть использован, например, алгоритм Форчуна [12].
Рис. 3. Диаграмма Вороного
Пару соседних точек определим как пару точек, ячейки Вороного для которых имеют общее ребро.
( А* А* Л
Пусть уг I ют, * - —,* + у I - значение показателя (16), рассчитанное для г-ой пары
соседних точек из множества В, г = 1, Я, А* - фиксированная ширина временного окна. При фиксированном значении частоты шт совокупность значении показателя (16), рассчитанных для каждои пары соседних точек можно рассматривать как вектор характерных признаков у (*1 , описывающих временной ряд X в момент времени * (см.
(2)):
У ('К X =
г г
Ух
А* А* Л
«т > * ~ — , * + — | — Уя
( А* А* ^
Ют , * "Т , * + А
* = 1, Т . (18)
V ч 2 2 у V 2 2 уу Компоненты этого вектора характеризуют степень выраженности бегущей волны частоты шт между всеми парами соседних точек в области О.
Аналогично может быть определён вектор у (* | ) характерных признаков,
составленный из значений показателя (17). Далее в обозначении вектора характерных признаков будем опускать нижний индекс у или Г, если вид этого показателя неважен.
Для сокращения размерности вектора у (* | ) может быть проведена
пространственная группировка пар соседних точек. Разобьём область О на множество областей О1,...,Од (рис. 4). Каждую область Оф q = 1,Q, характеризуем средним значением показателя (16):
_ ( А* А* Л 1 ^ ( А* А* Л , — ^(«т,*—,*+Т]=^ I Тг[,*-у• *+Т) * = Т-
Вектор характерных признаков у(V | ют) , описывающих временной ряд X в момент времени *, составляем из рассчитанных средних значений:
г г
А* А*
Г А* А*^
Аналогично определяем вектор >'(/1 шж).
Микросостояниями ^,..., 8К временного ряда X являются непересекающиеся области в пространстве предложенных характерных признаков. В настоящей работе для их выделения используется кластерный анализ. Пусть у(1| )у (Т | ) - множество векторов, характеризующих временной ряд X в каждый момент времени. В результате их кластеризации будут найдены центы кластеров w1,...,(эталонные представители микросостояний). Принадлежность вектора х(*) значений временного ряда X в момент времени * определяется по критерию близости соответствующего вектора характерных признаков у (* | ) к эталонному представителю микросостояния (см. (4)). Для
кластеризации данных у (1| « ),...,у (Т | ) может быть использован, например, метод ^-средних [13]. Отметим, что микросостояния ^,..., 8К зависят от анализируемой частоты шт бегущей волны.
а б
4. Результаты экспериментальных исследований
4.1. Описание исходных данных
Экспериментальные исследования предложенных алгоритмов сегментации проводились на данных электроэнцефалографии головного мозга (ЭЭГ), записанных в НБИКС-центре НИЦ «Курчатовский институт». Запись проводилась в состоянии покоя испытуемого в течение нескольких минут, число каналов электроэнцефалографа Ь = 62, расположение каналов соответствует стандартной схеме «10-20» (рис. 5), частота дискретизации /ц = 5000 Гц. Для проведения исследований выбран промежуток записи длиной 10 с (Т = 50000).
0.5 -0.4 0.3 0.2 0.1 0 -0.1 -0.2 -0.3 -0.4 -0.5
-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5 Рис. 5. Расположение электродов электроэнцефалографа
Согласно [14], работу головного мозга можно представить как последовательность сменяющих друг друга микросостояний. Выдвигается гипотеза, что микросостояния, связанные с характеристиками бегущих волн ЭЭГ, будут иметь высокую степень выраженности и повторяемости, что позволит говорить от их объективности и использовать в дальнейших нейрокогнитивных исследованиях.
4.2. Исследование показателей выраженности бегущих волн ЭЭГ
Для расчёта спектров сигналов ЭЭГ использовано дискретное оконное преобразование Фурье. В качестве оконной функции выбрана оконная функция Хемминга, ширина окна равнялась 4096 отсчетов (« 0,82 с.). Расчёт спектров для каждого канала ЭЭГ проводился независимо. На рис. 6 представлены примеры спектрограмм для некоторых каналов ЭЭГ.
0 0.5 1 1.5 2 2.5 3 3.5 4 0 0.5 1 1.5 2 2.5 3 3.5 4
Время, с Время, с
а б
Рис. 6. Примеры спектрограмм для каналов F5 (а) и FC3 (б)
На основе результатов спектрального разложения вычислены показатели уДшт, tx, t2)
и Гг:/ (шт, tx, t2) выраженности бегущей волны в соответствии с формулами (16) и (17) для
всех пар соседних электродов и на всех частотах ат от 5 до 30 Гц с шагом 1 Гц. Ширина окна i2 - f изменялась в зависимости от анализируемой частоты ат:
t2 -ti = f . (21)
а
т
Таким образом, на каждой частоте ат расчёт показателей проводился по числу временных отсчётов, соответствующих одному периоду гармоники анализируемой частоты ат.
Значения показателей уДшт, tx, t2) и Г (ют, tx, t2) для одной из пар соседних
электродов приведены на рис. 7. Из рисунка видно, что на некоторых частотах, в частности, на частотах альфа-ритма (8-13 Гц) и тета-ритма (5-8 Гц) присутствуют
интервалы времени, характеризующиеся высокой выраженностью бегущей волны, разделенные моментами десинхронизации. Моменты резкой десинхронизации обсуждаются в работе [14]. Предполагается, что в эти моменты происходит переход от одного микросостояния к другому.
а
б
Рис. 7. Значения показателей синхронности у„ (шт, *1, *2 ) (а) и когерентности Г.. («т, , *2 ) (б) для пары
соседних электродов Б5 и БС3
Картина выраженности бегущей волны у обоих показателей схожа, незначительная разница проявляется в крутизне межсинхронных переходов и значениях показателей на низких частотах.
На рис. 8 приведены диаграммы значений показателей выраженности бегущей волны на частоте 6 Гц для всех пар соседних электродов. Из диаграмм видно, что выделяются группы пар электродов со схожей динамикой выраженности бегущей волны, при этом показатель синхронности (16) более отчётливо отражает закономерность смены степени выраженности волны.
О 0.5 1 1.5 2 2.5 3 3.5 4 0 0.5 1 1.5 2 2.5 3 3.5 4
Время, с Время, с
а б
Рис. 8. Значения показателей синхронности угДшт, ^, *2 ) (а) и когерентности Г^ («,, , *2 ) (б) бегущих
волн на частоте 6 Гц для всех пар соседних электродов.
4.3. Исследование качества сегментации ЭЭГ
Признаки у (7 | ), у (71 )г, в пространстве которых проводится выделение
микросостояний, зависят от анализируемой частоты сош, а признаки у(71 со/н ) , >'(7 | со/н) -
ещё и от разбиения множества пар соседних электродов на группы (см. (18)—(20)). В зависимости от выбора частоты шт, разбиения множества пар электродов на группы и числа кластеров К, представление исходного временного ряда как последовательности номеров микросостояний будет получаться различным. Для оценки качества формируемой последовательности номеров микросостояний введём следующие показатели.
1. Коэффициент детерминации.
Я2 = , (22)
—общ
где —внутр и —общ - внутрикластерная дисперсия и общая дисперсия кластеризуемых
2 и данных. Показатель Я2 характеризует степень выраженности кластеров (микросостояний)
в исследуемом пространстве признаков.
2. Средняя продолжительность микросостояний.
К _
_ Е п А7к
А? = -, (23)
к=1_ К
Е пк
к=1
где Пк - число непрерывных фрагментов временного ряда, соответствующих к-му микросостоянию, А^ - средняя продолжительность к-го микросостояния, к = 1, К . Малые
значения показателя А? означают частую смену микросостояний, что может говорить о чрезмерной детализации сегментов.
3 . Среднеквадратичное отклонение продолжительности микросостояний.
а -
К
А-=1
К
Е пк
к=1
(24)
где дк - среднеквадратичное отклонение продолжительности к-то микросостояния,
к = 1,К. Высокие значения показателя а говорят о сильной вариабельности продолжительностей микросостояний.
На рис. 9, 10 представлены значения показателей (22)-(24), рассчитанные для различных частот шт в пределах от 5 Гц до 30 Гц с шагом 1 Гц и различного числа микросостояний К. На рис. 11 приведены значения показателей (22)-(24), рассчитанные по результатам кластеризации в классических пространствах признаков [15]: в
пространстве, составленном из амплитуд гармоник на частотах от 5 Гц до 30 Гц с шагом 1 Гц, и в пространстве исходных амплитуд сигнала ЭЭГ (размерности Ь = 62). Кластеризация проводилась с помощью метода ^-средних.
Рис. 9. Показатели качества сегментации Л2, А? и о в пространствах признаков, связанных с показателем синхронности. По оси абсцисс отложено число микросостояний К, по оси ординат - анализируемая частота бегущей волны ют. Группировка пар соседних электродов на Q = 4 и Q = 9 групп проводилась в
соответствии со схемой на рис. 4
Рис. 10. Показатели качества сегментации Л2, А? и (Т в пространствах признаков, связанных с показателем когерентности. По оси абсцисс отложено число микросостояний К, по оси ординат - анализируемая частота бегущей волны ют. Группировка пар соседних электродов на Q = 4 и Q = 9 групп проводилась в
соответствии со схемой на рис. 4
Рис. 11. Показатели качества сегментации Я^, А/ и О в пространствах признаков, сформированных из амплитуд гармоник частоты ю по всем каналам (верхняя часть) и из амплитуд исходного сигнала ЭЭГ по всем каналам (нижняя часть). На верхних графиках по оси абсцисс отложено число микросо стояний К, по
оси ординат - частота гармоник ю
Из рисунков видно, что свойства микросостояний в пространствах, связанных с показателями синхронности и показателями когерентности, сильно схожи. С увеличением числа микросостояний наблюдается тенденция к уменьшению средней продолжительности микросостояний при уменьшении среднеквадратичного отклонения продолжительностей. Качество кластеризации (показатель Я2) почти не зависит от анализируемой частоты бегущей волны, при этом с увеличением частоты средняя продолжительность микросостояний уменьшается. Причиной этому может быть то, что на более высоких частотах интервалы времени, на которых бегущие волны обладают схожими характеристиками, становятся уже. Иными словами, бегущие волны на высоких частотах становятся более стохастичными.
В сравнении с классическими пространствами признаков для сегментации, кластеры в построенных пространствах признаков, связанных с характеристиками бегущих волн, в большинстве случаев более выражены при проведении группировки пар соседних электродов. Исключение составляет пространство, составленное из амплитуд гармоник на частотах альфа-ритма 8-12 Гц, где показатель Я для фиксированного числа кластеров
принимает максимальные значения. Кластеры в пространстве, составленном из показателей синхронности/когерентности для всех пар соседних электродов менее чётко выражены, что, вероятно, связано с высокой размерностью этого пространства.
На рис. 12 приведён график последовательности смены номеров микросостояний для одного из лучших (по показателю ^ ) результатов сегментации на К = 8 микросостояний, и визуально изображены их эталонные представители. Характеристики микросостояний
для полученной сегментации: Я2 « 0,86, Л/ ~ 0,15 с, а « 0,12 с).
ос ^
х к
о ь
о о о о
о. ^
5
о. ф
г
о
X
3 4 5 6
время, с а
о;
1 6
о;
О 5 о о о.
о.
а) -1 5 1 о
X
н
-
1
10 95
0.9
0.85
■О 75
2 3
номер группы
б
Рис. 12. График последовательности смены номеров микросостояний для сегментации на К = 8 микросостояний в пространстве показателей синхронности (7 | С)ш ) на частоте ют = 8 Гц при
группировке пар соседних электродов на 2 = 4 группы (а) и эталонные представители микросостояний (б). Группировка проводилась в соответствии со схемой на рис. 4а
Заключение
В работе предложен новый метод сегментации многомерных временных рядов, обладающих пространственной организацией, на основе характеристик бегущих волн. Для оценки степени выраженности бегущих волн использованы показатели синхронности и когерентности. Расчёт этих показателей проводился по результатам анализа кросс-спектров сигналов, соответствующих пространственно близким точкам наблюдения.
В результате обработки реальных экспериментальных данных электроэнцефалографии головного мозга показана эффективность предложенного метода сегментации. Качество сегментации в пространстве характеристик бегущих волн в ряде случаев оказывается выше, чем при использовании традиционных методов. Полученные сегменты легко интерпретируемы, каждый сегмент отражает степень синхронности изменений электрофизиологических процессов в различных областях мозга. Таким образом, анализ характеристик выделенных сегментов и взаимосвязей между ними может
дать важную информацию о поведении исследуемой пространственно распределённой системы с точки зрения пространственной синхронности протекающих в ней процессов.
В развитие работы авторы планируют усовершенствование применяемых методов выделения бегущих волн и методов кластеризации в пространствах признаков, связанных с их характеристиками. Результаты сегментации сигналов ЭЭГ предполагается применить при моделировании ЭЭГ и фМРТ-активностей головного мозга с целью выделения его когнитивных микросостояний и построения функциональных коннектомов.
Исследование выполнено за счет гранта Российского научного фонда (проект № 14-28-00234, Отделение нейрокогнитивных и социогуманираных наук НБИКС-Центра НИЦ «Курчатовский институт»).
Список литературы
1. Abonyi J., Feil B. Cluster Analysis for Data Mining and System Identification. Birkhauser Basel Publ., 2007. 306 p. DOI: 10.1007/978-3-7643-7988-9
2. Percival D.B., Walden A.T. Wavelet Methods for Time Series Analysis. Cambridge University Press, 2006. (Cambridge Series in Statistical and Probabilistic Mathematics).
3. Seiffert U.S., Jain L.C., eds. Self-organizing Neural Networks: Recent Advances and Applications. Physica-Verlag HD Publ., 2002. (Ser. Studies in Fuzziness and Soft Computing; vol. 78). DOI: 10.1007/978-3-7908-1810-9
4. Liao T. Warren. Clustering of time series data—a survey // Pattern Recognition. 2005. Vol. 38, no. 11. P. 1857-1874. DOI: 10.1016/j.patcog.2005.01.025
5. Varstal M., Millan J.R., Heikkonen J. A recurrent self-organizing map for temporal sequence processing // In: Artificial Neural Networks—ICANN'97. Springer Berlin Heidelberg, 1997. P. 421-426. (Ser. Lecture Notes in Computer Science; vol. 1327). DOI: 10.1007/BFb0020191
6. Pusheng Zhang, Yan Huang, Shashi Shekhar, Vipin Kumar. Correlation analysis of spatial time series datasets: A filter-and-refine approach // In: Advances in Knowledge Discovery and Data Mining. Springer Berlin Heidelberg, 2003. P. 532-544. (Ser. Lecture Notes in Computer Science; vol. 2637). DOI: 10.1007/3-540-36175-8 53
7. Park J.Y. The spatial analysis of time series. Working Papers 2005-07. Rice University, Department of Economics, 2005. Available at:
http://econpapers.repec.org/paper/eclriceco/2005-07.htm , accessed 01.09.2014.
8. Kamarianakis Y. Spatial time series modeling: A review of the proposed methodologies. Regional Economics Applications Lab Technical Series, University of Illinois at Cham-paign-Urbana, REAL 03-T-19.
9. Sandstede B. Stability of travelling waves // In: Fiedler B., ed. Handbook of Dynamical Systems. Vol. 2. North-Holland, Amsterdam, 2002. P. 983-1055.
10. Nunez P L., Srinivasan R., Westdorp A.F., Wijesinghe R.S., Tucker D.M., Silberstein R.B., Cadusch P.J. EEG coherency: I: statistics, reference electrode, volume conduction,
Laplacians, cortical imaging, and interpretation at multiple scales // Electroencephalography and Clinical Neurophysiology. 1997. Vol. 103, no. 5. P. 499-515.
11. Hartwig J., Schnitzspahn K.M., Kliegel M., Velichkovsky B.M., Helmert J.R. I see you remembering: What eye movements can reveal about process characteristics of prospective memory // International Journal of Psychophysiology. 2013. Vol. 88, no. 2. P. 193-199. DOI: 10.1016/j.ijpsycho.2013.03.020
12. Fortune S. A Sweepline algorithm for Voronoi diagrams // Algorithmica. 1987. Vol. 2, no. 1-4. P. 153-174. DOI: 10.1007/BF01840357
13. Jain A.K. Data clustering: 50 years beyond K-means // Pattern Recognition Letters. 2010. Vol. 31, no. 8. P. 651-666. DOI: 10.1016/j.patrec.2009.09.011
14. Фингелькурц А.А. Пространственно-временная организация сегментной структуры ЭЭГ человека: дис. ... канд. биол. наук. М., МГУ, 1998. 273 с.
15. Lehmann D., Strik W.K., Henggeler B., Koenig T., Koukkou M. Brain electric microstates and momentary conscious mind states as building blocks of spontaneous thinking: I. Visual imagery and abstract thoughts // International Journal of Psychophysiology. 1998. Vol. 29, no. 1. P. 1-11. DOI: 10.1016/S0167-8760(97)00098-6
p t n i t | Science and Education of the Bauman MSTU,
science &Lducation 2014 °10 pp 114136
r . _ ^m^mwT DOI: 10.7463/1014.0728495
of the Bauman MSTU Re e d
ISSN 1994-0448 © Bauman Moscow State Technical Unversity
Spatial Time Series Segmentation Method Based on Traveling Waves
A.G. Trofimov1'2'*, I.V. Kolodkin1, V.L. Ushakov2, 'atrofimovaiast.ru
B.M. Velichkovsky1,2
1NRNU «MEPhI», Moscow, Russia 2NRC «Kurchatov Institute», Moscow, Russia
Keywords: spatial time series, travelling-wave, segmentation, phase-locking value, coherence, eeg microstates, electroencephalogram.
This paper describes a method for spatial time series segmentation based on the analysis of traveling waves and demonstrates its efficiency using real data.
Phase-locking value and coherence measure are used to evaluate the intensity of the traveling waves. The calculation of these indicators has been based on the cross-spectral analysis of the signals corresponding to the spatially adjacent points of observation.
The described method can be applied to analysis of spatially organized physical, chemical, and biological processes. In particular, the method can be used to neurophysiological, geother-mal and ocean researches and to studies of the distributed solar dynamics.
Experimental researches of the method have been carried out on the real data of human brain electroencephalography (EEG). It is shown that the proposed method outperforms the classical methods of EEG segmentation. The results achieved have a simple interpretation.
The first section provides a mathematical formulation of the multivariate time series segmentation problem. The second section gives a formal description of the traveling wave and indicators of traveling wave intensity between two spatial points. The third section considers the method of the features extraction based on the characteristics of traveling waves for multivariate time series segmentation. The fourth section presents the simulation results of the proposed segmentation algorithm on real EEG data. In conclusion the main results and future lines of research are discussed.
References
1. Abonyi J., Feil B. Cluster Analysis for Data Mining and System Identification. Birkhäuser Basel Publ., 2007. 306 p. DOI: 10.1007/978-3-7643-7988-9
2. Percival D.B., Walden A.T. Wavelet Methods for Time Series Analysis. Cambridge University Press, 2006. (Cambridge Series in Statistical and Probabilistic Mathematics).
3. Seiffert U.S., Jain L.C., eds. Self-organizing Neural Networks: Recent Advances and Applications. Physica-Verlag HD Publ., 2002. (Ser. Studies in Fuzziness and Soft Computing; vol. 78). DOI: 10.1007/978-3-7908-1810-9
4. Liao T. Warren. Clustering of time series data—a survey. Pattern Recognition, 2005, vol. 38, no. 11, pp. 1857-1874. DOI: 10.1016/j.patcog.2005.01.025
5. Varstal M., Millan J.R., Heikkonen J. A recurrent self-organizing map for temporal sequence processing. In: Artificial Neural Networks—ICANN'97. Springer Berlin Heidelberg, 1997, pp. 421-426. (Ser. Lecture Notes in Computer Science; vol. 1327). DOI: 10.1007/BFb0020191
6. Pusheng Zhang, Yan Huang, Shashi Shekhar, Vipin Kumar. Correlation analysis of spatial time series datasets: A filter-and-refine approach. In: Advances in Knowledge Discovery and Data Mining. Springer Berlin Heidelberg, 2003, pp. 532-544. (Ser. Lecture Notes in Computer Science; vol. 2637). DOI: 10.1007/3-540-36175-8 53
7. Park J.Y. The spatial analysis of time series. Working Papers 2005-07. Rice University, Department of Economics, 2005. Available at: http://econpapers.repec.org/paper/eclriceco/2005-07.htm , accessed 01.09.2014.
8. Kamarianakis Y. Spatial time series modeling: A review of the proposed methodologies. Regional Economics Applications Lab Technical Series, University of Illinois at Champaign-Urbana, REAL 03-T-19.
9. Sandstede B. Stability of travelling waves. In: Fiedler B., ed. Handbook of Dynamical Systems. Vol. 2. North-Holland, Amsterdam, 2002, pp. 983-1055.
10. Nunez PL., Srinivasan R., Westdorp A.F., Wijesinghe R.S., Tucker D.M., Silberstein R.B., Cadusch P.J. EEG coherency: I: statistics, reference electrode, volume conduction, Laplacians, cortical imaging, and interpretation at multiple scales. Electroencephalography and Clinical Neurophysiology, 1997, vol. 103, no. 5, pp. 499-515.
11. Hartwig J., Schnitzspahn K.M., Kliegel M., Velichkovsky B.M., Helmert J.R. I see you remembering: What eye movements can reveal about process characteristics of prospective memory. International Journal of Psychophysiology, 2013, vol. 88, no. 2, pp. 193-199. DOI: 10.1016/j.ijpsycho.2013.03.020
12. Fortune S. A Sweepline algorithm for Voronoi diagrams. Algorithmica, 1987, vol. 2, no. 1-4, pp. 153-174. DOI: 10.1007/BF01840357
13. Jain A.K. Data clustering: 50 years beyond K-means. Pattern Recognition Letters, 2010, vol. 31, no. 8, pp. 651-666. DOI: 10.1016/j.patrec.2009.09.011
14. Fingelkurts A.A. Prostranstvenno-vremennaia organizatsiia segmentnoi struktury EEG cheloveka. Kand. diss. [Time-spatial organization of the human EEG segmental structure. Ph.D. Diss.]. Moscow, MSU, 1998. 273 p. (in Russian).
15. Lehmann D., Strik W.K., Henggeler B., Koenig T., Koukkou M. Brain electric microstates and momentary conscious mind states as building blocks of spontaneous thinking: I. Visual imagery and abstract thoughts. International Journal of Psychophysiology, 1998, vol. 29, no. 1, pp. 1-11. DOI: 10.1016/S0167-8760(97)00098-6
