CC BY
45
Прикладные задачи
м^^и^^^м^м1 нелинейной теории колебании и волн
УДК 530.182 https://doi.org/10.18500/0869-6632-2019-27-6-8-24
Влияние нелинейности на оценки связанности осцилляторов методом частной направленной когерентности
Д. И. Астахова1'3, М.В. Сысоева2'3, И. В. Сысоев1'3
Саратовский национальный исследовательский государственный университет им. Н.Г. Чернышевского
Россия, 410012 Саратов, Астраханская, 83 2Саратовский государственный технический университет им. Гагарина Ю.А. Россия, 410054 Саратов, Политехническая, 77 3Институт радиотехники и электроники им. В.А. Котельникова РАН, Саратовский филиал Россия, 410019 Саратов, Зелёная, 38 E-mail: astakhovadi@mail.ru, bobrichek@mail.ru, ivssci@gmail.com Автор для переписки Дарья Ивановна Астахова, astakhovadi@mail.ru Поступила в редакцию 18.07.2019, принята к публикации 13.10.2019, опубликована 2.12.2019
Цель настоящего исследования - определить, может ли метод частной направленной когерентности правильно определять направленное взаимодействие между нелинейными системами при нелинейных связях между этими системами и в случае, когда измеренные сигналы порождены объектами высокой размерности (ансамблями). Также определить зависимость результатов оценки связанности методом частной направленной когерентности от параметров: длины реализации, частоты дискретизации, размерности модели и от архитектуры связей в системах. Методы. Ансамбли из четырёх связанных осцилляторов различных типов при различном внесении связей между ними и высокоразмерная динамическая математическая модель эпилепсии использованы в качестве тестовых систем. Для определения значимости результатов использовались суррогатные временные ряды, построенные путём перестановки реализаций. Результаты. Показано, что в ансамблях маломерных осцилляторов архитектура связи может быть правильно выявлена для линейных и нелинейных систем, связанных как линейной, так и нелинейной связью. Для сложных составных сигналов, когда каждый измеряемый временной ряд представляет собою сумму сигналов многих отдельных осцилляторов, методика оказывается недостаточно специфична, выявляя несуществующие связи, и недостаточно чувствительна, пропуская имеющиеся. Заключение. Сформулированы критерии применения метода частной направленной когерентности к различным сигналам. При достаточной длине ряда, частоте выборки и размерности мера не показывает опосредованных связей в отличие от попарных методов причинности по Грейнджеру и энтропии переноса и хорошо работает для зашумлённых временных рядов. Метод позволяет изучать связанность в ансамбле из произвольного числа осцилляторов и можно определить, на каких частотах происходит взаимодействие. По сравнению с методом причинности по Грейнджеру, для которого работоспособность декларируется уже при 4-16 характерных периодах, метод частной направленной когерентности даёт приемлемые результаты для рядов от 80 характерных периодов.
Ключевые слова: частная направленная когерентность, связанность, нелинейность, нелинейные системы, нелинейная связь, многомерные системы.
Образец цитирования: Астахова Д.И., Сысоева М.В., Сысоев И.В. Влияние нелинейности на оценки связанности осцилляторов методом частной направленной когерентности//Известия вузов. ПНД. 2019. T. 27, № 6. С. 8-24. https://doi.org/10.18500/0869-6632-2019-27-6-8-24
Финансовая поддержка. Работа выполнена при поддержке РФФИ, грант № 17-02-00307.
8
© Астахова Д.И., Сысоева М.В., Сысоев И.В., 2019
https://doi.org/10.18500/0869-6632-2019-27-6-8-24
Effect of nonlinearity on coupling estimations between oscillators using partial directed coherence approach
D. I. Astakhova1,3, M. V. Sysoeva2,3, I. V. Sysoev1,3
1 Saratov State University 83, Astrakhanskaya, Saratov 410012, Russia 2 Yuri Gagarin Saratov State Technical University 77, Politekhnicheskaya, Saratov 410054, Russia 3V.A. Kotelnikov Institute of Radioengineering and Electronics of RAS, Saratov Branch 38, Zelenaya, Saratov 410019, Russia E-mail: astakhovadi@mail.ru, bobrichek@mail.ru, ivssci@gmail.com Correspondence should be addressed to Daria Ivanovna Astakhova, astakhovadi@mail.ru Received 18.07.2019, accepted 13.10.2019, published 2.12.2019
The purpose of this work is to determine the ability of the partial directed coherence method to identify directed interactions between nonlinear systems correctly in presence of nonlinear couplings between systems, as well as in the case of measured signals generated by objects of high dimension. The other purpose is to determine the dependence of the coupling estimation results on the parameters: series length, sampling rate, model dimension and coupling architecture. Methods. Ensembles composed of four differently coupled oscillators and dynamical mesoscale model of epilepsy are considered as test systems. Surrogate time series constructed by permutation of realization are used to determine the significance of the results. Results. Coupling architecture in ensembles of small-dimensional oscillators can be correctly identified for linear and nonlinear systems in both cases of linear and nonlinear coupling. For complex composite signals, when each measured time series is the sum of signals from many individual oscillators, the technique is not specific enough, revealing non-existent connections, and it is not sensitive enough, missing the existing ones. Outcomes. The criteria for applying the partial directed coherence method to different signals are formulated. The measure does not show indirect couplings at sufficient series length, sampling rate and model dimension in contrast to the pairwise methods like Granger causality or transfer entropy. The measure works well for noisy time series. The method allows to study connectivity in an ensemble of arbitrary number of oscillators. The method allows to determine at what frequencies the interaction occurs. The partial directed coherence method gives acceptable results for series of length of 80 and more characteristic periods in comparison with the Granger causality method, for which the efficiency is declared already at 4-16 characteristic periods.
Key words: partial directed coherence, coupling, nonlinearity, nonlinear systems, nonlinear coupling, high-dimensional system.
Reference: Astakhova D.I., Sysoeva M.V., Sysoev I.V. Effect of nonlinearity on coupling estimations between oscillators using partial directed coherence approach. Izvestiya VUZ. Applied Nonlinear Dynamics, 2019, vol. 27, no. 6, pp. 8-24. https://doi.org/10.18500/0869-6632-2019-27-6-8-24
Acknowledgements. This work was supported by Russian Foundation for Basic Research, grant no. 17-02-00307.
Введение
Взаимодействие между отделами головного мозга в норме и при патологии активно изучается в настоящее время. Основным источником информации об активности мозга являются электроэнцефалограммы (ЭЭГ) и магнитоэнцефалограммы (МЭГ). Для исследования связанности структур мозга применялись самые различные меры, как линейные, так и нелинейные [1,2]. Значительное внимание уделялось сопоставлениям линейных и нелинейных мер [1,3-6], а также исследованиям синхронизации в мозге [7]. Накоплен значительный фактический материал по детектированию связанности между структурами мозга во многих специфических случаях, например, для изучения математических способностей [8], при абсансной эпилепсии [9,10], височной эпилепсии средней доли [11], при паркинсонизме [12]. При этом важная и до сих пор далёкая от окончательного решения проблема - интерпретация результатов оценок связанности [13].
В данной работе исследуются новые возможности и ограничения метода частной направленной когерентности (partial directed coherence, PDC), предложенного в [14]. Это частотно разрешённый подход к описанию взаимосвязей (направления информационного потока) между многомерными временными рядами, основанный на декомпозиции многомерных частных коге-рентностей, вычисленных по многомерным авторегрессионным моделям. Метод был введён в 1999 году группой Луиса Баккала [15] первоначально для двух сигналов и расширен на множество сигналов двумя годами позже [14]. Было заявлено, что PDC показывает прямое влияние одной подсистемы на другую, не показывая косвенных связей [16]. Данный метод основан на построении линейных прогностических моделей, неизвестные коэффициенты которых оцениваются методом наименьших квадратов. Затем коэффициенты, зависящие от времени, преобразуются в коэффициенты, зависящие от частоты. Поскольку сами модели, как и Фурье преобразование, линейны, эффективность для нелинейных систем неочевидна.
Хотя данная мера не столь популярна, как, например, причинность по Грейнджеру [17], что может быть обусловлено относительною сложностью процедуры расчёта и интерпретации результатов, она имеет два несомненных преимущества и используется рядом учёных из различных научных коллективов [5, 18-20]. Во-первых, PDC позволяет при оценке связанности получить разрешение результатов по частоте (альтернативным подходом являются различные вариации частотно разрешённого метода причинности по Грейнджеру), во-вторых, в модель изначально закладывается возможная зависимость каждого из исследуемых сигналов от всех прочих, что подразумевает разделение прямых и опосредованных взаимодействий по построению методики. В случае метода причинности по Грейнджеру используют различные варианты условной причинности [21,22].
К настоящему времени появились модификации данного метода (см., например, [16,20]), позволяющие повысить эффективность метода для определённого типа сигналов, однако в данной работе будет рассмотрена оригинальная методика, предложенная в [14]. Основные модификации, сделанные в работах [16,20], касаются нормировок меры для того, чтобы иметь возможность сравнивать друг с другом абсолютные значения PDC, полученные для исследования сигналов, имеющих различную амплитуду. Однако в случае, если необходимо сделать вывод о наличии или отсутствии связанности (а не о силе воздействия), опирающийся на статистическое тестирование с помощью суррогатных временных рядов, предложенные в [16,20] изменения не должны из общих соображений играть существенную роль.
В сравнительно недавних работах PDC была применена к результатам магниторезонанс-ной томографии (МРТ) [23], сигналам сердечно-сосудистой системы [19], и хорошо показала себя в этих исследованиях. Полученные в этих работах результаты были не только правдоподобны (не противоречили ранее проведённым исследованиям других учёных), но и содержали ряд новых результатов, лучше объясняющих наблюдаемые данные.
Несмотря на большое число тестирований изложенного в [14] подхода на различных модельных примерах [24], вопрос о том, насколько метод способен выявлять связанность в нелинейных системах по сравнительно коротким временным рядам (от 80 характерных периодов колебаний) остаётся открытым. Кроме этого, непонятно, как метод будет реагировать на нелинейную связь и/или параметрическую связь между подсистемами.
Цель настоящего исследования - определить, будет ли частная направленная когерентность правильно определять направленное взаимодействие нелинейных по природе систем и выявлять нелинейную связь, а также проверить зависимость результатов оценки связанности от параметров метода и данных: длины реализации, частоты дискретизации, размерности модели и от архитектуры связей в системах.
1. Методика
х1(п) Р х1(п - г) w1 (п)
= £ Аг +
XN (п)_ г= 1 xN (п — г) WN (п)_
1.1. Частная направленная когерентность. Метод заключается в следующем: для начала представим сигнал в виде
(1)
где |хга}^=1 - исходный векторный временной ряд, состоящий из N скалярных рядов (полученных от N осцилляторов, возможная связанность между которыми исследуется) и измеренный с шагом выборки Д£, р - число предыдущих моментов времени (размерность модели), учтённых в модели и, как следствие, число матриц Аг, Аг - матрицы коэффициентов авторегрессионных моделей, состоящие из значений а^(г), оцениваемых в типичном случае методом наименьших квадратов, г - сдвиг во времени. Элементы матрицы Аг отражают линейное влияние сдвинутого во времени на гДЪ ^'-го ряда на г-й ряд в текущий момент времени п, где г, ] - номера исследуемых осцилляторов (подсистем).
Следует отметить, что авторегрессионные модели, содержащие коэффициенты а^(г), на практике часто удобнее строить не в форме (1), а в форме
N
хг(п) = £ £ а%,3 (г)х (п - г) + -Шг(п).
г=13=р+1
(2)
Уравнение (2) даёт возможность оценить г-е строки всех матриц Аг сразу. Чтобы получить все матрицы Аг целиком, реконструкцию моделей (2) нужно провести для всех компонентов сигнала
г = !,..., N.
С помощью 2-преобразования перейдём из временной области в частотную, получив вместо набора матриц Аг одну матрицу А(/), зависящую от частоты. Для этого воспользуемся формулой
А(/) = £ Аг ¿-г |г=е-«я/.
Г=1
Далее удобнее работать с матрицей А(/), получаемой из матрицы А(/) по формуле
А(/) = Е — А(/),
(3)
(4)
где Е - единичная матрица.
Используя введённые обозначения, можно посчитать частную направленную когерент ность РОС следующим образом:
Чз(/)
РОСгз (/) =
оЯ *
а^ (/)£-1 аз (/)
(5)
где / - частота, Н* - операция эрмитова сопряжения (взятие комплексного сопряжения элементов плюс транспонирование), (/) являются элементами матрицы А(/) (скалярные значения, то есть просто числа), аз - столбцы матрицы А(/) (векторы, или в терминах программирования -одномерные массивы).
1.2. Параметры метода и архитектуры связей. В классическом методе частной направленной когерентности фактически есть один изменяющийся параметр: количество анализируемых предыдущих значений р. Этот параметр ещё можно назвать размерностью модели. Также
а
А
© ©
b
d
ФЧ2) © ©
Рис. 1. Архитектура связей. a - первый осциллятор воздействует на второй, первый и третий взаимно влияют друг на друга, четвёртый осциллятор колеблется изолированно; b - две пары двунаправленно связанных осцилляторов: первый со вторым и третий с четвёртым; c - кольцо из трёх осцилляторов: первый воздействует на второй, второй -на третий, а третий - на первый, четвёртый осциллятор вновь изолированный; d - однонаправленная связь между первым и вторым осцилляторами, третий и четвертый осцилляторы колеблются изолированно
Fig. 1. Coupling architecture. a - the first oscillator leads the second, the first and the third mutually affect each other, the fourth oscillator oscillates in isolation; b - two pairs of bidirectionally coupled oscillators: the first with the second and the third with the fourth; c - the ring of three oscillators: the first drives the second, the second drives the third, and the third drives the first, the fourth oscillator is again isolated; d - unidirectional coupling between the first and second oscillators, the third and fourth oscillators oscillate in isolation
С
есть ещё три параметра анализируемых тестовых временных рядов: частота дискретизации, длина анализируемого участка и архитектура связей.
Чтобы воспроизвести основные возможные варианты связанности подсистем (однонаправленную связь, двунаправленную связь, отсутствие взаимодействия, опосредованное воздействие), для каждого типа систем были построены маломерные ансамбли, состоящие из 4 осцилляторов. Для каждого типа систем осцилляторы связывали четырьмя способами, которые схематически представлены на рис. 1.
a) Первый воздействует на второй, первый и третий взаимно влияют друг на друга, четвёртый осциллятор колеблется изолированно (рис. 1, a). При такой архитектуре связей, кроме выявления однонаправленного и двунаправленного взаимодействия, можно проверить, как метод реагирует на опосредованную связь с третьего на второй осциллятор, а также реакцию метода на изолированный осциллятор.
b) Две пары двунаправленно связанных осцилляторов: первый со вторым и третий с четвёртым (рис. 1, Ь). При такой архитектуре связей проверяем реакцию метода на двунаправленное взаимодействие, и как метод реагирует на отсутствующее взаимодействие между парами связанных осцилляторов.
c) Кольцо из трёх осцилляторов: первый воздействует на второй, второй - на третий, а третий - на первый, четвёртый осциллятор вновь изолированный (рис. 1, с). В такой системе каждый осциллятор в кольце влияет на каждый (частично опосредованно), но непосредственных двунаправленных взаимодействий нет.
ё) Однонаправленная связь между первым и вторым осцилляторами, третий и четвертый осцилляторы колеблются изолированно (рис. 1, ё). Здесь в первую очередь проверяется, как метод реагирует на одностороннее воздействие без опосредованной связи. А также интересна реакция метода на отсутствующее взаимодействие между третьим и четвёртым осцилляторами.
1.3. Статистический анализ полученных результатов. Для определения работоспособности меры введём критерии оценивания [25].
1. Число ложноположительных заключений не должно превышать порог уровня значимости (хорошая специфичность метода);
2. Положительные связи должны быть обнаружены верно, если таковые имеются (хорошая чувствительность метода).
Первый критерий является обязательным, так как на отсутствии детектирования опосредованных и ложных связях строится вся методика.
Для определения статистической значимости полученных значений РОС были построены суррогатные временные ряды путём перестановки реализаций. Для этого рассматривался ансамбль из 10 реализаций, полученных при идентичных параметрах и архитектуре связей, но различных начальных условиях. Генерировались все возможные перестановки из 10 по 4 (число осцилляторов), полученные наборы из четырёх чисел использовались как номера реализаций, которые брались для составления суррогатных. Таким образом, каждый суррогатный векторный временной ряд представлял собой набор скалярных реализаций для различных осцилляторов, соответствующих разным исходным рядам. Поскольку общее число реализаций, которые можно составить таким образом, очень велико (для перестановок из 10 по 4 составляет 5040), было решено ограничиться двадцатью, что даёт возможность определить значимость на общепризнанном уровне 0.05. Эти 20 суррогатных реализаций были выбраны из всего ансамбля так, чтобы каждый номер исходного ряда примерно одинаковое число раз встречался в соответствующей перестановке на каждой позиции.
1.4. Тестовые системы. В численном эксперименте использовались два типа ансамблей. Во-первых, маломерные ансамбли, состоящие из четырёх однотипных осцилляторов; при этом доступными измерению считались временные ряды каждого отдельного осциллятора. Во-вторых, ансамбли большой размерности (динамическая мезомасштабная модель эпилепсии [26]), состоящие из 500 модельных нейронов ФитцХью-Нагумо, для которых в качестве измеряемого сигнала рассматривались и анализировались только суммарные сигналы больших групп нейронов, как это обычно бывает при измерении сигналов локальных потенциалов мозга в эксперименте (всего три сигнала, соответствующие коре, вентропостериальному медиальному и ретикулярному ядрам таламуса), а временные ряды отдельных осцилляторов не были доступны.
Для удобства сопоставления с экспериментальными данными локальных потенциалов мозга и сигналами динамической мезомасштабной модели эпилепсии (которая была получена ранее в [26]), где измерение осуществлялось с частотами 512, 1024, 2048 и 4096 Гц, а основной период колебаний составлял примерно 1/8 с, перенормируем время: будем считать, что 1 единица безразмерного времени (100 шагов интегрирования) соответствует 0.02441406 с. Таким образом, эффективная частота выборки будет составлять 212 Гц. Такой подход позволит, с одной стороны, использовать значения параметров для далее рассмотренных модельных систем из областей или вблизи широко распространённых в литературе значений, а с другой - даст возможность переносить полученные в работе результаты на изучение, например, абсансных приступов, для которых характерна частота порядка 8 Гц [27].
1.4.1. Линейный осциллятор (линейная система с линейной связью). Так как мера является линейной, логично, что мы начнём тестировать ее на примере системы, состоящей из линейных осцилляторов, связанных линейной связью:
Частота гармонических колебаний для четырёх осцилляторов немного варьировалась: mi = = 0.95, m2 = 0.99, m3 = 0.98, m4 = 0.96. При воздействии i-го осциллятора на j-й все значения ki,j для связей, присутствующих в соответствующей архитектуре, задавались следующим образом (для отсутствующих связей они были равны нулю): k1,j = 0.65, k2,j = 0.35, k3,j = 0.45, k4,j = 0.87. Коэффициенты линейной диссипации выбирались из соображений возможной устойчивой генерации колебаний в присутствии связей: у1 = 0.01, у2 = 0.08, у3 = 0.03, у4 = 0.07. Источником энергии этих колебаний выступал динамический белый нормальный шум ^i(t) с нулевым средним и среднеквадратичным отклонением an = 2.5. При значительно меньших значе-
d Xi dxi
(6)
j=i
К
О 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100
. 3->2 . 4->2
О §
0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 , 2->3 , 4->3
О §
о §
0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 1->4 . 2->4 . 3->4
О 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 /, Hz /, Hz /, Hz /, Hz
Рис. 2. Зависимость значений PDC от частоты для линейного осциллятора с архитектурой связей, показанной на рис. 1, a: 1 ^ 2,1 о 3. Частота дискретизации fd = 512 Гц, временной ряд длиной 50 с, количество анализируемых предыдущих значений p = 100 точек. Чёрным цветом обозначена частная направленная когерентность, посчитанная по действительным временным рядам, серым цветом обозначен суррогатный уровень. На диагонали показаны усреднённые спектры мощности в логарифмическом масштабе
Fig. 2. Dependence of the PDC values on the frequency for linear oscillator with the coupling architecture shown in Fig. 1, a: 1 ^ 2,1 о 3. Sampling rate was equal to fd = 512 Hz, time series length was 50 s, the number of analyzed previous values was p =100 points. Black color shows the partial directed coherence calculated from the real time series, gray color indicates the surrogate level. The diagonal shows the averaged power spectra on a logarithmic scale
1
1
1
ниях уг (например, в два раза), чем рассмотренные, добротность колебаний отдельных элементов ансамбля оказывалась столь велика, что из-за связей амплитуда сигнала неограниченно росла. При значительно больших величинах у г роль шума в динамике существенно возрастала, а частота колебаний слишком сильно отличалась от ю, что неудобно при дальнейшем рассмотрении.
Уравнения интегрировались методом Эйлера-Маруямы с шагом Н = 0.01, поскольку в систему вводился гауссовский шум. Методы Рунге-Кутты высоких порядков (например, часто используемый метод 4-го порядка) в таком случае не имеют практического смысла, так как требуют дополнительных предположений о природе шума [28]. Шаг интегрирования подбирался эмпирически: величина Н уменьшалась до тех пор, пока временные реализации, полученные при двух различающихся значениях Н (в нашем случае Н = 0.01 и Н = 0.005) статистически не отличались (точное наложение траекторий в данном случае не может быть достигнуто из-за наличия шума).
Связанные линейные осцилляторы находились в режиме колебаний с несколькими модами на частотах 3...8 Гц (двумя или тремя) с сильно разнесёнными парциальными частотами, что обусловлено большим коэффициентом связи и высокою добротностью [29]. Ведомые осцилляторы при этом демонстрировали суперпозицию колебаний на собственной частоте и вынужденных колебаний. Уединённый осциллятор колебался на собственной частоте (рис. 2). Источником энергии всех колебаний выступал шум.
1.4.2. Система ФитцХью-Нагумо (нелинейная система с линейной связью). Усложним задачу и сконструируем систему, состоящую из нелинейных осцилляторов, связанных линейной связью. Например, систему из четырех связанных осцилляторов ФитцХью-Нагумо:
Жг(а» - Жг)(Жг - 1) - Уг + 1а,г + ^Ху,
вгХг - угУг (7)
dt
dyi dt
О 20 40 60 80 100 °0 20 40 60 80 100 °0 20 40 60 80 100 °0 20 40 60 80 100
3->2
4->2
Ü §
0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 , 2->3 , 4->3
Ü §
о
!_
О:
Qh- , -1---- ,--- -,- -г-—-г—— о
0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 . 2->4 . 3->4
Ü §
1
0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100 /, Hz /, Hz /, Hz /, Hz
1
1
1
1
1
1
Рис. 3. Зависимость значений PDC от частоты для системы ФитцХью-Нагумо с архитектурой связей, показанной на рис. 1 b: 1 о 2, 3 о 4. Частота дискретизации fd = 512 Гц, временной ряд длиной 50 с, количество анализируемых предыдущих значений p = 100 точек. Чёрным цветом обозначена частная направленная когерентность, посчитанная по действительным временным рядам, серым цветом обозначен суррогатный уровень. На диагонали показаны усреднённые спектры мощности в логарифмическом масштабе
Fig. 3. Dependence of the PDC values on the frequency for the FitzHugh-Nagumo system with the coupling architecture shown in Fig. 1, b: 1о2, 3о4. Sampling rate was equal to fd=512 Hz, time series length was 50 s, the number of analyzed previous values was p=100 points. Black color shows the partial directional coherence calculated from the real time series, gray color indicates the surrogate level. The diagonal shows the averaged power spectra on a logarithmic scale
Параметры модели использовались следующие: коэффициенты связи k1,j = 0.35, k2,j = 0.38, k3,j = 0.39, k4,j = 0.37, внешний стимул Ia,i = 0.85, коэффициент ai = 0.8, коэффициент ß менялся от 0.15 до 0.17, а коэффициент у - от 0.06 до 0.068; ^(t) - динамический белый нормальный шум с нулевым средним и среднеквадратичным отклонением an = 2.5. Параметры были в целом взяты из работы [30], где они подбирались таким образом, чтобы воспроизводить характеристики лимбических разрядов.
Все рассмотренные системы ФитцХью-Нагумо находились в режиме нелинейных колебаний (под режимом близких к линейным колебаний мы понимаем режим, соответствующий низкоамплитудной близкой к гармонической генерации), мощность была распределена в спектре в значительном диапазоне частот, при этом основная мощность приходилась на диапазон 3...5 Гц (рис. 3).
1.4.3. Генератор с жестким возбуждением и потенциалом Тоды (нелинейная система с нелинейной связью). Ещё усложним задачу и введём систему, состоящую из нелинейных осцилляторов, связанных нелинейной связью (8). Данные осцилляторы можно назвать осцилляторами c жёстким возбуждением и потенциалом Тоды (далее сокращённо «ж/в-Тоды»).
$ - + £ ki,jx? - x4 j df + » (1 - e-Xi) = ^i(t), (8)
где i - номер текущего осциллятора, j - номер воздействующего осциллятора, коэффициент связи ki,j = 0.65, r1 = -0.05, r2 = -0.08, r3 = -0.06, r4 = -0.07, ^i(i) - нормальный белый шум с нулевым средним и среднеквадратичным отклонением as = 2.5, частота гармонических колебаний »i = 1 (одинаковая для всех четырёх осцилляторов), шаг интегрирования методом Эйлера-Маруямы равен 0.01 [30].
Ü §
о
^Vn^^--——
100 °0 1
20
40 60 4->2
80 100
0 20
40 60
80
100 О 1
20
40 60 2->3
80 100 °0 20 40 60 80
О
О 20
40 60 1->4
80
ü §
О
100 °0 1
100 °0 1
20
40 60 4->3
80 100
20
40 60 2->4
80
О 20
40 60 /, Hz
100 О 1
20
40 60 3->4
. 0.
80 100 0 20 40 60 80 100
80 100 0 20
40 60 /, Hz
80 100 0 20
40 60 /, Hz
80 100 0 20
40 60 /, Hz
80 100
Рис. 4. Зависимость значений PDC от частоты для генератора с жёстким возбуждением и потенциалом Тоды с архитектурой связей, показанной на рис. 1, с: 1 ^ 2, 2 ^ 3, 3 ^ 1. Частота дискретизации fd = 512 Гц, временной ряд длиной 10 с, количество анализируемых предыдущих значений p = 100 точек. Чёрным цветом обозначена частная направленная когерентность, посчитанная по действительным временным рядам, серым цветом обозначен суррогатный уровень. На диагонали показаны усреднённые спектры мощности в логарифмическом масштабе
Fig. 4. Dependence of the PDC values on the frequency for the Van der Pol like-Toda oscillators with the coupling architecture shown in Fig. 1, с: 1^2, 2^3, 3^1. Sampling rate was equal to fd=512 Hz, time series length was 10 s, the number of analyzed previous values was p=100 points. Black color shows the partial directed coherence calculated from the real time series, gray color indicates the surrogate level. The diagonal shows the averaged power spectra on a logarithmic scale
1
1
1
Связанные генераторы находились в режиме сильно нелинейных релаксационных колебаний с основной частотой ~ 8 Гц и большим числом высших гармоник. Такой режим рассматривался как простейшая модель абсансов в [31]. Уединённый осциллятор находился в режиме подпороговых индуцированных шумом колебаний на собственной частоте ~ 8 Гц, которые в [31] сопоставлялись фоновой динамике (см. рис. 4).
1.4.4. Динамическая мезомасштабная модель абсансной эпилепсии (многомерная нелинейная система с большим числом степеней свободы). В качестве экспериментальных моделей абсансной эпилепсии наиболее широко используются генетические модели грызунов. Они позволяют изучать нейробиологические основы данного заболевания с сопровождающими его пик-волновыми комплексами. Однако для тестирования связанности экспериментальные сигналы плохо подходят, поскольку задать архитектуру связей в них напрямую нельзя. Поэтому в данной работе вместо реальных экспериментальных данных использовали динамическую мезо-масштабную модель (далее «ДММ») пик-волновых разрядов (9), схематически представленную на рис. 6, a (см. ниже), с одной стороны, отражающую сложный состав колебаний, с другой -имеющую априорно известную архитектуру связей. Подробно данная система была разработана в [32] и улучшена в [26]. Модель состоит из системы связанных ОДУ, описывающих модельные нейроны типа ФитцХью-Нагумо для трансмембранного потенциала, и включает 3 типа клеток -1 из коры и 2 из таламуса (каждый тип моделируется в виде сосредоточенной системы). Кроме того, учитываются специфические свойства гамма-аминомасляной кислоты (ГАМК) и используется сигмоидальная передаточная функция в связях.
xi(ai - xi)(xi - 1) - yi + ^ ki,jh(xj(t - t)),
j=i
ßixi - Yi yi, (9)
где функция h(x) = 1 + th(x)/2, а = 0.8, ß = 0.008, y = 0.0033. Запаздывание t подбиралось
dt
dyi dt
VPM->RTN
Cortex-s-RTN
О 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100
RTN->VPM
Cortex-s-VPM
О 20 40 60 80 100 0 20 40 60 80 100 0 20 40 60 80 100
RTN-^Cortex
VPM-s-Cortex
ü §
0
20
40 60 /, Hz
80
100 °0
20
40 60 /, Hz
80
100 0
20
40 60 /, Hz
80
100
Рис. 5. Зависимость значений PDC от частоты для динамической мезомасштабной модели. Частота дискретизации fd = 425 Гц, временной ряд длиной 10 с, количество анализируемых предыдущих значений p = 100 точек. Чёрным цветом обозначена частная направленная когерентность, посчитанная по действительным временным рядам, серым цветом обозначен суррогатный уровень. На диагонали показаны усреднённые спектры мощности в логарифмическом масштабе
Fig. 5. Dependence of the PDC values on the frequency for the dynamic mesoscale model. Sampling rate was equal to fd = 425 Hz, time series length was 10 s, the number of analyzed previous values was p = 100 points. Black color shows the partial directional coherence calculated from the real time series, gray color indicates the surrogate level. The diagonal shows the averaged power spectra on a logarithmic scale
1
1
1
1
1
для каждой матрицы связей k индивидуально, в диапазоне от 9 до 13 единиц времени и было одинаковым для всех связей внутри одной матрицы. Уравнения модели численно интегрировались методом Эйлера-Маруямы с шагом 0.5.
В итоге получались сильно нелинейные колебания с основной частотой ~ 8 Гц и большим числом высших гармоник (рис. 5).
2. Результаты
Численный эксперимент: изменение архитектуры связей, частоты дискретизации, размерности модели и длины ряда было организовано одинаково для всех исследуемых типов осцилляторов, кроме динамической мезомасштабной модели.
Всего было рассмотрено четыре типа подсистем: линейный осциллятор, система ФитцХью-Нагумо, генератор «ж/в-Тоды» и ДММ. Для первых трёх систем генерировалось четыре связанных осциллятора и исследовались все четыре типа архитектур связей, представленных на рис. 1, по 10 рядов по 100 секунд каждый. Для ДММ архитектура связей оставалась неизменной, как показано ниже на рис. 6, a (выбиралась из физиологических соображений); выбирались реализации длиной не менее 10 с.
Также было четыре изменяющихся параметра: частота дискретизации, длина анализируемого участка, количество анализируемых предыдущих значений и архитектура связей, которые выбирались вручную. Частота дискретизации fd бралась 256 Гц, 512 Гц, 1024 Гц. Количество анализируемых предыдущих значений p бралось от 50 до 250 с шагом 50. Длина временного ряда составляла 10 с, 25 с, 50 с, 75 с и 100 с.
2.1. Линейный осциллятор. На рис. 2 представлен типичный график, построенный по временным реализациям линейного осциллятора.
Для линейного осциллятора мера при всех наборах параметров правильно выявляет заложенную архитектуру связей: правильно определяется одно- и двунаправленное воздействие, не выявляются опосредованные связи, PDC для изолированного узла ведёт себя адекватно.
Для линейного осциллятора получились идеальные значения, поэтому были расширены границы тестирования: выставлена минимальная длина ряда 3 с (24 характерных периода колебаний), а минимальный сдвиг - 10 точек. Как выяснилось, для этой системы метод способен работать и при таких значениях длины ряда, что не наблюдалось для других осцилляторов. Результат следует признать в целом ожидаемым, так как мера является линейной.
2.2. Система ФитцХью-Нагумо. Для осциллятора ФитцХью-Нагумо мера при всех наборах параметров правильно выявляет заложенную архитектуру связей: правильно определяется одно- и двунаправленное воздействие, РОС для изолированного узла ведёт себя адекватно (рис. 3). При увеличении длины ряда чувствительность метода (способность выявлять правильные взаимодействия) увеличивается. При увеличении количества анализируемых предыдущих значений усиливается влияние шумов измерения, график становится «изрезанным», но архитектура связей всё равно определяется верно. Соответственно, можно сделать вывод, что метод частной направленной когерентности хорошо работает не только для линейных систем, но и для нелинейных, для чего необходимо увеличивать длину реализации.
2.3. Генератор с жёстким возбуждением и потенциалом Тоды. При малой длине ряда (10 с, см. рис. 4) чувствительность метода оказывается недостаточна при кольцевой архитектуре связей, и существующие связи не выявляются: хотя соответствующие пики в спектре имеются, суррогаты также демонстрируют на этих же частотах повышенные значения РОС, полностью перекрывая значения, полученные для данных. С увеличением длины ряда уровень суррогатов падает, а уровень РОС не меняется. Для остальных архитектур связей (а, Ь, ё) метод при всех параметрах правильно определяет направление взаимодействий.
Для «ж/в-Тоды» есть особенность: очень высокий уровень РОС для изолированного узла, хотя суррогатный уровень тоже высокий, поэтому влияние со стороны изолированного осциллятора всё же оказывается незначимым. Это происходит из-за того, что мы используем оригинальный метод подсчёта РОС, который чувствителен к отношению дисперсий рассматриваемых сигналов. Поскольку у изолированного осциллятора «ж/в-Тоды» низкая амплитуда колебаний в сравнении с остальными осцилляторами системы, наблюдается высокий уровень РОС со стороны этого узла на другие и, наоборот, очень низкий уровень РОС на этот узел. Данный эффект подробно рассмотрен в работе [25].
2.4. Динамическая мезомасштабная модель абсансной эпилепсии. Сравнив полученные методом частной направленной когерентности связи рис. 6, Ь со схемой, представленной на рис. 6, а, можно заметить, что метод показал влияние НТК (ретикулярное таламическое ядро) на кору, которое в модель изначально заложено не было, при этом влияние с НТК на УРЫ (вентропостериальное медиальное таламическое ядро) метод в нужном частотном диапазоне выявить не смог.
Рис. 6. Архитектура связи для динамической мезомасштабной модели: a - по построению, b - реконструированная. Чёрные стрелки показывают реально существующие связи, серая стрелка показывает ложно выявленную связь, штриховая стрелка - пропущенную методом связь
Fig. 6. Coupling architecture for dynamic mesoscale model: a - as constructed, b - as reconstructed. Black arrows show the actual couplings, gray arrow shows the falsely identified coupling, and the bar arrow shows the missing by method PDC coupling
a
Интересной особенностью работы метода для ДММ является то, что связь выявляется как значимая в широком частотном диапазоне. Это может быть вызвано наличием высших гармоник в спектре (см. рис. 5) или иными, пока неясными причинами.
Также можно констатировать, что для сильно связанных систем, где лишь малое число связей отсутствуют, выявление прямых и опосредованных связей достаточно затруднено. Аналогичные выводы были сделаны для метода причинности по Грейнджеру в работе [22]. Также стало ясно, что непосредственное применение рассматриваемой методики даже к модельным данным высокой сложности, в которых каждый сигнал представляет собой не реализацию отдельного сосредоточенного элемента, а суммарный сигнал большого числа осцилляторов, связанных как друг с другом, так и с осцилляторами, генерирующими другие сигналы, грозит серьёзными ошибками. Следовательно, о применении к нейрофизиологическим или климатическим данным такого типа также следует говорить с осторожностью.
3. Обсуждение и заключение
Была исследована зависимость результатов оценки частной направленной когерентности от длины используемого временного ряда, числа анализируемых предыдущих значений модели, частоты дискретизации для ансамблей линейных осцилляторов и осцилляторов второго порядка при разной архитектуре связей и различном типе нелинейных функций, а также для мезомас-штабной модели таламо-кортикальной системы головного мозга во время абсансной эпилепсии.
Результаты проведённого анализа показывают, что с помощью РОС архитектура связи может быть правильно выявлена для линейных и нелинейных систем, связанных как линейной, так и нелинейной связью. Но работает это только в случае систем, представляющих собой сосредоточенные элементы. Когда метод был применён к системе, представляющей собой сумму большого количества отдельных осцилляторов, то метод выявил несуществующую связь (плохая специфичность метода), при этом пропустил одно реальное воздействие (плохая чувствительность метода).
Что касается зависимости работоспособности меры от параметров метода, можно заключить следующее.
1. Зависимость от длины реализации. С увеличением длины анализируемого ряда значения когерентности на частотах, на которых нет взаимодействия, а также суррогатный уровень падают, а значения когерентности на релевантных частотах остаются практически неизменными. В этом смысле мера ведёт себя полностью аналогично классической ненаправленной функции когерентности. Поэтому фактически присутствующие связи выявляются всё лучше. Таким образом, можно констатировать, что увеличение длины ряда положительно сказывается на работоспособности меры.
2. Зависимость от частоты дискретизации. При тех значениях частоты дискретизации, которые мы рассматривали, метод везде показывал достаточную чувствительность и специфичность. Это может быть обусловлено тем, что всегда рассматривались достаточно длинные ряды и имелось как минимум 32 точки на характерный период. Также следует иметь в виду, что большая частота выборки требует большого числа коэффициентов модели, чтобы захватить основные медленные масштабы колебаний.
3. Зависимость от размерности модели. С увеличением размерности модели работоспособность метода сначала улучшается (повышается чувствительность), а потом - падает, так как специфичность утрачивается - появляются ложные связи на произвольных, по всему спектру, частотах. Этот эффект, по-видимому, обусловлен двумя факторами: во-первых, недостаточным усреднением при оценке коэффициентов (чем больше размерность, тем больше коэффициентов и хуже статистические свойства их оценок при той же длине ряда), во-вторых, численными ошибками алгоритма, которые начинают играть существенную роль
при количестве коэффициентов ~ 1GGG. Возможно, чтобы компенсировать второй эффект, следует использовать специализированные библиотеки алгоритмов, реализующие расчёт с точностью представления данных выше стандартных б4-80 бит для чисел с плавающею запятой.
4. Зависимость от архитектуры связи. Эффективность исследуемого подхода оказалась сильно зависима от архитектуры связей. Лучше всего метод работает при взаимодействии в парах. Хуже всего - в кольце, как было показано на примере со связанными генераторами с жёстким возбуждением и потенциалом Тоды (8). Изолированные узлы показывают высокое (хотя, как правило, незначимое по суррогатам) влияние на другие узлы, низкое на них и умеренное между ними, что обусловлено чувствительностью метода к дисперсии сигнала, как показано в [1б].
В результате данной работы можно сформулировать следующие преимущества метода:
• при достаточной длине ряде, частоте выборки и размерности мера не показывает опосредованных связей в отличие от попарных методов причинности по Грейнджеру и энтропии переноса;
• можно изучать связанность в ансамбле из произвольного числа осцилляторов;
• хорошо работает для зашумлённых временных рядов [2S];
• можно определить, на каких частотах происходит взаимодействие; это может быть актуально для анализа связанности в популяциях нейронов, где различные частоты генерации спайков и бёрстов соответствуют разных биологическим процессам.
Также в процессе работы с методом были выявлены недостатки:
• по сравнению с методом причинности по Грейнджеру, для которого работоспособность декларируется уже при 4-1б характерных периодах [4,30], метод частной направленной когерентности даёт приемлемые результаты для рядов от 8O колебаний;
• для сложных составных сигналов, когда каждый измеряемый временной ряд представляет собою сумму сигналов многих отдельных осцилляторов, методика оказывается недостаточно специфична, выявляя несуществующие связи, и недостаточно чувствительна, пропуская имеющиеся.
Библиографический список
1. Gourevitch B., Le Bouquin-Jeannes R., Faucon G. Linear and nonlinear causality between signals: Methods, examples and neurophysiological applications // Biological Cybernetics. 200б. Vol. 9S(4). P. 349-3б9.
2. Sakkalis V. Review of advanced techniques for the estimation of brain connectivity measured with EEG/MEG // Computers in Biology and Medicine. 2011. Vol. 41, iss. 12. P. lllO-lln.
3. Pereda E., QuianQuiroga R., Bhattacharya J.Nonlinear multivariate analysis of neurophysiological signals // Progress in Neurobiology. 200S. Vol. TT, iss. 1-2. P. 1-3T.
4. Sysoeva M.V., Sitnikova E., Sysoev I.V., Bezruchko B.P., van Luijtelaar G. Application of adaptive nonlinear Granger causality: Disclosing network changes before and after absence seizure onset in a genetic rat model // Journal of Neuroscience Methods. 2014. Vol. 22б. P. 33-41.
5. Курганский А.В. Некоторые вопросы исследования корково-корковых функциональных связей с помощью векторной авторегрессионной модели многоканальной ЭЭГ // Журнал высшей нервной деятельности им. И.П. Павлова. 2010. Т. б0(б). С. T40-TS9.
6. Курганский А.В. Количественные меры кортико-кортикального взаимодействия: ^временное состояние // Физиология человека. 2013. Т. 39(4). С. 112-122.
T. Sakkalis V., Doru Giurcaneanu C., Xanthopoulos P., Zervakis M.E., Tsiaras V., Yang Y., Kara-konstantaki E., Micheloyannis S. Assessment of linear and nonlinear synchronization measures
for analyzing EEG in a mild epileptic paradigm // IEEE Transactions on Information Technology in Biomedicine. 2009. Vol. 13, iss. 4. P. 433-441.
8. Micheloyannis S., Sakkalis V., Vourkas M., Stam C.J., Simos P.G. Neural networks involved in mathematical thinking: Evidence from linear and non-linear analysis of electroencephalographs activity // Neuroscience Letters. 2005. Vol. 373, iss. 3. P. 212-217.
9. Luttjohann A., van Luijtelaar G. The dynamics of cortico-thalamo-cortical interactions at the transition from pre-ictal to ictal LFPs in absence epilepsy // Neurobiology of Disease. 2012. Vol. 47. P. 47-60.
10. Sysoeva M.V., Luttjohann A., van Luijtelaar G. and Sysoev I.V. Dynamics of directional coupling underlying spike-wave discharges // Neuroscience. 2016. Vol. 314. P. 75-89.
11. Lehnertz K., Andrzejak R.G., Arnhold J., Kreuz T., Mormann F., Rieke C., Widman G., and Elger C.E. Nonlinear EEG analysis in epilepsy: Its possible use for interictal focus localization, seizure anticipation, and prevention // Journal of Clinical Neurophysiology. 2001. Vol. 18, iss. 3. P. 209-222.
12. Tass P., Smirnov D., Karavaev A., Barnikol U., Barnikol T., Adamchic I., Hauptmann C., Pawelcyzk N., Maarouf M., Sturm V., Freund H.-J., and Bezruchko B. The causal relationship between subcortical local field potential oscillations and Parkinsonian resting tremor // Journal of Neural Engineering. 2010. Vol. 7. 016009.
13. Smirnov Dmitry A. Quantifying causal couplings via dynamical effects: A unifying perspective // Phys. Rev. E. 2014. Vol. 90. 062921.
14. Baccala L.A., Takahashi D.Y. Partial directed coherence: A new concept in neural structure determination // Biological Cybernetics. 2001. Vol.84. P.272-273.
15. Sameshima K., Baccala L.A. Using partial directed coherence to describe neuronal ensemble interactions // Journal of Neuroscience Methods. 1999. Vol. 94. P. 93-103.
16. Takahashi D.Y., Baccala L.A., Sameshima K. Connectivity inference between neural structures via partial directed coherence // Journal of Applied Statistics. 2007. Vol. 34(10). P. 1255-1269.
17. Granger C.W.J. Investigating causal relations by econometric models and cross-spectral methods // Econometrica. 1969. 37(3). P. 424-438.
18. Baccala L.A., Takahashi D.Y., Sameshima K. Directed transfer function: Unified asymptotic theory and some of its implications // IEEE Transactions on Biomedical Engineering. 2016. Vol. 63(12). P. 2450-2460.
19. Milde T., Schwab K., Walther M., EiseltM., Schelenz C., Voss A., Witte H. Time-variant partial directed coherence in analysis of the cardiovascular system: A methodological study // Physiological measurement. 2011. Vol. 32. P. 1787-1805.
20. Schelter B., Timmer J., Eichler M. Assessing the strength of directed influences among neural signals using renormalized partial directed coherence // Journal of Neuroscience Methods. 2009. Vol.179. P.121-130.
21. Chen Y., Rangarajan G., Feng J., Ding M. Analyzing multiple nonlinear time series with extended Granger causality // Phys. Lett. A. 2004. 324(1). P. 26-35.
22. Корнилов М.В., Сысоев И.В. Реконструкция архитектуры связей в цепочке из трех однона-правленно связанных систем методом причинности по Грейнджеру // Письма в ЖТФ. 2018. Т. 44, № 10. С. 86-95.
23. Sato J.R., Takahashi D.Y., Arcuri S.M., Sameshima K., Morettin P.A., Baccala L.A. Frequency domain connectivity identification: An application of partial directed coherence in fMRI // Human Brain Mapping. 2009. Vol. 30. P. 452-461.
24. Sommerlade L., Eichler M., Jachan M., Henschel K., Timmer J., Schelter B. Estimating causal dependencies in networks of nonlinear stochastic dynamical systems // Physical Review E. 2009. Vol. 80, iss. 5.
25. Smirnov D., Schelter B., Winterhalder M., Timmer J. Revealing direction of coupling between neuronal oscillators from time series: Phase dynamics modeling versus partial directed coherence // Chaos. 2007. Vol. 17. 013111(11).
26. Medvedeva T.M., Sysoeva M.V., Sysoev, I.V. Coupling analysis between thalamus and cortex in mesoscale model of spike-wave discharges from time series of summarized activity of model neurons // Proceedings of 2nd School on Dynamics of Complex Networks and their Application in Intellectual Robotics. 2018. 8589208. P. 137-138.
27. Вольнова А.Б., Ленков Д.Н. Абсансная эпилепсия: механизмы гиперсинхронизации нейронных ансамблей. Медицинский академический журнал. 2012. Т. 12, № 1. С. 7-19.
28. Тихонов В.И., Миронов М.А. Марковские процессы. М.: Советское радио, 1977. 488 с.
29. Рабинович М.И., Трубецков Д.И. Введение в теорию колебаний и волн. НИЦ «Регулярная и хаотическая динамика». 2000. 560 с.
30. Сысоева М. В., Медведева Т. М. Оптимизация параметров метода причинности по Грейн-джеру для исследования лимбической эпилепсии // Известия вузов. ПНД. 2018. Т. 26, № 5. С. 39-62.
31. Сысоева М.В., Кузнецова Г.Д., Сысоев И.В. Модел^ование стгналов элекфоэнцефало-грамм ^ью ^и абcанcной эпилепсти в ^иложении к анализу cвязанноcти между отделами мозга// Биофизика. 2016. Т. 61(4). С. 782-792.
32. Medvedeva T.M., SysoevaM.V., van Luijtelaar G., Sysoev I.V. Modeling spike-wave discharges by a complex network of neuronal oscillators // Neural Networks. 2018. Vol. 98. P. 271-282.
References
1. Gourevitch B., Le Bouquin-Jeannes R., Faucon G. Linear and nonlinear causality between signals: Methods, examples and neurophysiological applications. Biological Cybernetics, 2006, vol. 95(4), pp. 349-369.
2. Sakkalis, V. Review of advanced techniques for the estimation of brain connectivity measured with EEG/MEG. Computers in Biology and Medicine, 2011, vol. 41, iss. 12, pp. 1110-1117.
3. Pereda E., QuianQuiroga R., Bhattacharya J. Nonlinear multivariate analysis of neurophysiological signals. Progress in Neurobiology, 2005, vol. 77, iss. 1-2, pp. 1-37.
4. Sysoeva M.V., Sitnikova E., Sysoev I.V., Bezruchko B.P., van Luijtelaar G.. Application of adaptive nonlinear Granger causality: Disclosing network changes before and after absence seizure onset in a genetic rat model. Journal of Neuroscience Methods, 2014, vol. 226, pp. 33-41.
5. Kurgansky A.V. Study of cortico-cortical functional connectivity with vector autoregressive model of multichannel EEG. Zh Vyssh Nerv Deyat+, 2010. vol. 60(6), pp. 740-759 (in Russian).
6. Kurgansky A.V. Quantitative measures of cortical functional connectivity: A state-of-the-art brief survey. Human Physiology, 2013, vol. 39, no. 4, pp. 432-440.
7. Sakkalis V., Doru Giurcaneanu C., Xanthopoulos P., Zervakis M.E., Tsiaras V., Yang Y., Karakon-stantaki E., Micheloyannis S. Assessment of linear and nonlinear synchronization measures for analyzing EEG in a mild epileptic paradigm. IEEE Transactions on Information Technology in Biomedicine, 2009, vol. 13, iss. 4, pp. 433-441.
8. Micheloyannis S., Sakkalis V., Vourkas M., Stam C.J., Simos P.G. Neural networks involved in mathematical thinking: evidence from linear and non-linear analysis of electroencephalographic activity. Neuroscience Letters, 2005, vol. 373, iss. 3, pp. 212-217.
9. Luttjohann A., van Luijtelaar G. The dynamics of cortico-thalamo-cortical interactions at the transition from pre-ictal to ictal LFPs in absence epilepsy. Neurobiology of Disease, 2012, vol. 47, pp. 47-60.
10. Sysoeva M.V., Luttjohann A., van Luijtelaar G. and Sysoev I.V. Dynamics of directional coupling underlying spike-wave discharges. Neuroscience, 2016, vol. 314, pp. 75-89.
11. Lehnertz K., Andrzejak R.G., Arnhold J., Kreuz T., Mormann F., Rieke C., Widman G., and Elger C.E. Nonlinear EEG Analysis in Epilepsy: Its Possible Use for Interictal Focus Localization, Seizure Anticipation, and Prevention. Journal of Clinical Neurophysiology, 2001, vol. 18, iss. 3, pp. 209-222.
12. Tass P., Smirnov D., Karavaev A., Barnikol U., Barnikol T., Adamchic I., Hauptmann C., Pawelcyzk N., Maarouf M., Sturm V., Freund H.-J., and Bezruchko B. The causal relationship between subcortical local field potential oscillations and parkinsonian resting tremor. Journal of Neural Engineering, 2010, vol. 7, 016009.
13. Smirnov Dmitry A. Quantifying causal couplings via dynamical effects: A unifying perspective. Phys. Rev. E, 2014, vol. 90, 062921.
14. Baccala L.A., Takahashi D.Y. Partial directed coherence: a new concept in neural structure determination. Biological Cybernetics, 2001, vol. 84, pp. 272-273.
15. Sameshima K., Baccala L.A. Using partial directed coherence to describe neuronal ensemble interactions. Journal of Neuroscience Methods, 1999, vol. 94, pp. 93—103.
16. Takahashi D.Y., Baccala L.A., Sameshima K. Connectivity Inference between Neural Structures via Partial Directed Coherence. Journal of Applied Statistics, 2007, vol. 34(10), pp. 1255-1269.
17. Granger C.W.J. Investigating causal relations by econometric models and cross-spectral methods. Econometrica, 1969, vol. 37(3), pp. 424-438.
18. Baccala L.A., Takahashi D.Y., Sameshima K. Directed Transfer Function: Unified Asymptotic Theory and Some of its Implications. IEEE Transactions on Biomedical Engineering, 2016, vol. 63(12), pp. 2450-2460.
19. Milde T., Schwab K., Walther M., Eiselt M., Schelenz C., Voss A., Witte H. Time-variant partial directed coherence in analysis of the cardiovascular system: A methodological study. Physiological measurement, 2011, vol. 32, pp. 1787-1805.
20. Schelter B., Timmer J., Eichler M. Assessing the strength of directed influences among neural signals using renormalized partial directed coherence. Journal of Neuroscience Methods, 2009, vol. 179, pp. 121-130.
21. Chen Y., Rangarajan G., Feng J., Ding M. Analyzing multiple nonlinear time series with extended Granger causality. Phys. Lett. A., 2004, 324(1), pp. 26-35.
22. Kornilov M.V., and Sysoev I.V. Recovering the architecture of links in a chain of three unidirec-tionally coupled systems using the Granger-causality test. Technical Physics Letters, 2018, vol. 44, no. 5, pp. 445-449.
23. Sato J.R., Takahashi D.Y., Arcuri S.M., Sameshima K., Morettin P.A., Baccala L.A. Frequency domain connectivity identification: An application of partial directed coherence in fMRI. Human Brain Mapping, 2009, vol. 30, pp. 452-461.
24. Sommerlade L., Eichler M., Jachan M., Henschel K., Timmer J., Schelter B. Estimating causal dependencies in networks of nonlinear stochastic dynamical systems. Physical Review E, 2009, vol. 80, iss. 5.
25. Smirnov D., Schelter B., Winterhalder M., Timmer J. Revealing direction of coupling between neuronal oscillators from time series: Phase dynamics modeling versus partial directed coherence. Chaos, 2007, vol. 17, p. 013111(11).
26. Medvedeva T.M., Sysoeva M.V., Sysoev, I.V. Coupling analysis between thalamus and cortex in mesoscale model of spike-wave discharges from time series of summarized activity of model neurons. Proceedings of 2nd School on Dynamics of Complex Networks and their Application in Intellectual Robotics, 2018, 8589208, pp. 137-138.
27. Volnova A.B., Lenkov D.N. Absence epilepsy: Mechanisms of hypersyncronization of neuronal networks. Medical academic journal, 2012, vol. 12, no. 1, pp. 7-19.
28. Tikhonov V.I., Mironov M.A. Markovian processes. Moscow: Soviet Radio, 1977, 488 p.
29. Rabinovich M.I., Trubetskov D.I. Introduction to the Theory of Oscillations and Waves. «Regular and chaotic dynamics», 2000, 560 p. (in Russian).
30. Sysoeva M.V., Medvedeva T.M. Optimization of Granger causation method parameters for the study of limbic epilepsy. Izvestia VUZ. Applied Nonlinear Dynamics, 2018, vol. 26, iss. 5, pp. 39-62 (in Russian).
31. Sysoeva M.V. , Kuznetsova G.D., Sysoev I.V. Modelling EEG signals from rats when analysing absence epilepsy in application to analysis of coupling between brain areas Biophysics, 2016, vol. 61(4), pp. 661-669.
32. Medvedeva T.M., SysoevaM.V., van Luijtelaar G., Sysoev I.V. Modeling spike-wave discharges by a complex network of neuronal oscillators. Neural Networks, 2018, vol. 98, pp. 271-282.
Астахова Дарья Ивановна - родилась в городе Балаково Саратовской области (1997). Окончила с золотой медалью МАОУ СОШ № 27 с углубленным изучением отдельных предметов (2015). Окончила бакалавриат с красным дипломом факультета нано- и биомедицинских технологий СГУ имени Н.Г. Чернышевского по направлению «Биотехнические системы и технологии» (2019). На данный момент продолжает обучение в магистратуре по той же специальности.
Россия, 410012 Саратов, Астраханская, 83
Саратовский государственный университет имени Н.Г. Чернышевского Россия, 410019 Саратов, Зелёная, 38
Институт радиотехники и электроники им. В.А. Котельникова РАН, Саратовский филиал E-mail: astakhovadi@mail.ru
Сысоева Марина Вячеславовна - родилась в Саратове (1987), окончила с отличием бакалавриат (2009) и магистратуру (2011) факультета нано- и биомедицинских технологий СГУ имени Н.Г. Чернышевского по направлению «Биомедицинская инженерия». Защитила диссертацию на соискание учёной степени кандидата физико-математических наук на тему «Особенности реализации метода причинности по Грейнджеру для исследования электроэнцефалограмм при абсансной эпилепсии» по специальностям «Биофизика» и «Радиофизика» (2015, СГУ). С 2015 года работает на кафедре «Радиоэлектроника и телекоммуникации» Саратовского государственного технического университета им. Гагарина Ю.А. в должности доцента. Научные интересы - анализ временных рядов, нейронаука, математическое моделирование. Опубликовала более 20 статей в российских и зарубежных журналах.
Россия, 410054 Саратов, Политехническая, 77
Саратовский государственный технический университет имени Гагарина Ю.А. Россия, 410019 Саратов, Зелёная, 38
Институт радиотехники и электроники им. В.А. Котельникова РАН, Саратовский филиал E-mail: bobrichek@mail.ru
Сысоев Илья Вячеславович - родился в Саратове (1983), окончил факультет нелинейных процессов СГУ (2004), защитил диссертации на соискание учёной степени кандидата физико-математических наук (2007) и доктора физико-математических наук (2019). Доцент кафедры динамического моделирования и биомедицинской инженерии, ответственный секретарь редакционной коллегии журнала «Известия вузов. ПНД». Научные интересы - исследование сигналов биологической природы методами нелинейной динамики, исследование эффективности и модернизация подходов к анализу сигналов. Автор более 40 публикаций.
Россия, 410012 Саратов, Астраханская, 83
Саратовский государственный университет имени Н.Г. Чернышевского Россия, 410019 Саратов, ул. Зеленая, 38
Институт радиотехники и электроники им. В.А. Котельникова РАН, Саратовский филиал E-mail: ivssci@gmail.com
