CC BY
15
ОЦЕНКА НАПРАВЛЕНИЯ ВЗАИМОДЕЙСТВИЯ МЕЖДУ МОДЕЛЬНЫМИ СИСТЕМАМИ СВЯЗАННЫХ ОСЦИЛЛЯТОРОВ ПРИ СИЛЬНОЙ СВЯЗИ
В. С. Хорев
С помощью метода моделирования фазовой динамики исследована задача выявления по временным рядам значимого взаимодействия двух колебательных систем, а также преимущественного направления и времени задержки в связи для случая сильной связи систем. Рассмотрены модельные системы связанных осцилляторов с различными видами функций связи.
Ключевые слова: Сильная связь, взаимодействие, направление связи, задержка связи, колебательные системы.
Введение
Задача определения характера связи между двумя колебательными системами по временным рядам, которые могут быть зашумлены и хаотичны, возникает в различных приложениях [1,2]. Особенно сложно выявить направленность достаточно сильной связи, вызывающей синхронизацию систем, по коротким временным рядам - порядка нескольких десятков характерных периодов колебаний [3]. Для достижения успеха в этой сложной ситуации требуется высокая чувствительность метода диагностики связи, которую обеспечивает слежение за фазами колебаний [4], так как фаза является характеристикой процесса, чувствительной к внешним воздействиям на систему [5].
Один из наиболее чувствительных методов оценки связи между двумя осцилляторами был предложен в [6]. Он основан на моделировании фазовой динамики и применим в случае достаточно длинных временных рядов (несколько сотен характерных периодов) и малых шумов. В работе [7] метод обобщен на случай более коротких временных рядов и существенных шумов за счет специальных поправок и эмпирически найденного порогового значения характеристики связи, соответствующего 95%-й доверительной вероятности вывода о наличии связи. Метод нашел
применения в нейрофизиологии [8] и климатологии [9]. Он был развит далее в работе [10], где предложены модифицированные характеристики связи и получено аналитическое выражение для доверительной вероятности, с которой можно сделать вывод о наличии связи. Однако оба упомянутых улучшения метода ориентированы на случай систем со слабой связью, описываемой в моделях членами невысокого порядка в уравнениях фазовой динамики, в связи с чем возникает вопрос о применимости метода в случае сильной связи.
1. Метод моделирования фазовой динамики
Для оценки направления взаимодействия между двумя осцилляторными системами по их временным рядам использован метод, основанный на моделировании фазовой динамики. Основная идея метода - оценить, насколько сильно зависят будущие значения фазы одной системы от текущего значения фазы другой системы. Для этого по исходным временным рядам {жх(¿1),..., ж^м)} и ,..., ж2(¿м)} от
двух систем (где ¿г = ¿Д£ - интервал выборки, г = 1,..., N, N - длина ряда) рассчитываются временные ряды фаз колебаний {ф1(^1),..., ф1(¿м)} и {фг^), •••, ф2(¿м)}. На основе последних строится эмпирическая математическая модель [6], по коэффициентам которой оцениваются характеристики связанности. По временным рядам фаз построим модель динамики в виде
Дф1^) = *1(ф1(*), ф2^)+ ^(¿)),
(1)
Дфг(£) = ¿ЬЫ*), Ф1^)+ ^(¿)),
где Дфк (¿) = фк(£ + т) — фк(¿), к = 1, 2, т - конечный временной интервал; ек (¿) -шумы с нулевым средним; ¿к - тригонометрические многочлены следующего вида:
^\(ф1, Ф2, ai) = + (а^и eos(m^1 - Пф2) + вт,п sin(m^1 - пф2П ,
m,n v ' (2)
(2) _ / (2) (2) \
F2(ф2, Ф1, a2) = ад + ( am'n eos(m^2 - ПФ1) + pm'n sin(m^2 - пф1Л ,
т п ^ '
где ak = (а0к), {а^ге,™}, {Ргтп}) - векторы коэффициентов, k = 1, 2. Диапазон суммирования, то есть набор значений индексов m и n, определяющих, какие слагаемые присутствуют в многочлене, может быть свой для каждого осциллятора [11]. А сами оценки силы воздействия У12 рассчитываются следующим образом:
2П 2П
Vl'2 = / J Р1,2>ф2,Ьа1)2)/9ф2,1)2Йф1Йф2- (3)
о о
Данные оценки выводились для случая слабой связи, когда стационарное вероятностное распределение циклических фаз (ф1 mod 2п, ф2 mod 2п) является равномерным в квадрате [0, 2п) х [0, 2п) , а слагаемые в многочлене Fk при равномерной плотности распределения являются взаимно ортогональными функциями в этой области [11]. Вопрос о работоспособности данной методики в случае сильной связи остаётся открытым и требует дальнейшего исследования.
2. Тестовые модели
В качестве тестовых моделей выступали:
— фазовые осцилляторы;
— однонаправленно связанные осцилляторы ван дер Поля.
Для широкого круга ситуаций фазовая динамика осцилляторов, имеющих ярко выраженный ритм, адекватно описывается стохастическими дифференциальными уравнениями [12], поэтому в качестве наиболее простой и универсальной модели взаимодействующих систем были выбраны однонаправленно связанные фазовые осцилляторы следующего вида:
<1щ/<И = /1 +
(4)
^ф2М = /2 + к 8т(ф1 (£ — А) — ф2(¿)) +
где /1,2 - частоты осцилляторов, ^1,2(г) - белые шумы с нулевым средним, к - коэффициент связи между осцилляторами, А - задержка в связи между системами.
Для качественного описания динамики автоколебательных систем в ряде случаев хорошо подходит модель в виде осциллятора ван дер Поля [13, 14]. Поэтому в качестве второй исследуемой модели рассмотрим однонаправленно связанные осцилляторы ван дер Поля
(5)
где /1,2 - собственные частоты, X - коэффициент нелинейности, ^1,2 - белые шумы с нулевым средним, к - коэффициент связи между осцилляторами, А - задержка в связи между системами.
3. Результаты
3.1. Два однонаправленно связанных фазовых осциллятора. Рассмотрим зависимости оценки направления воздействия и времени запаздывания в связи от вариации различных параметров осцилляторов. При этом ответим на следующие вопросы:
— при каких условиях оценка времени задержки в связи остается несмещенной;
— при каких условиях вероятность ошибочных выводов о направлении связи не превышает 5%?
Кроме оценки количества ошибочных выводов важен вопрос и о том, при каких условиях достаточно высока вероятность правильных выводов о направлении связи. Для ответа на эти вопросы рассчитываем зависимости оценок сил связи от пробного времени задержки по временному ряду. Смещение оценки определяется разницей между значением, заданным в уравнениях модели, и значением, полученным по временному ряду. 95%-й уровень значимости рассчитывался по 100 реализациям суррогатных рядов, полученных из исходных временных рядов фаз осцилляторов путём случайной перестановки участков ряда длительностью 2п. При таком
способе приготовления суррогатных данных происходит разрушение связи между исследуемыми временными рядами при сохранении их основных свойств.
3.1.1. Влияние силы связи между осцилляторами на оценку времени задержки в связи. Используемые в эксперименте временные ряды систем (4) были получены путём численного интегрирования уравнений методом Эйлера с шагом 0.02п при частотах /1,2 ~ 0.1. Каждое из начальных условий выбиралось случайно из отрезка [0..2п] (использовалось равномерное распределение).
На рис. 1 представлены зависимости оценок сил связи от пробного времени задержки. Видно, что максимум У12 - оценки воздействия со стороны первого осциллятора на второй - соответствует реальному времени задержки в связи и превышает уровень значимости, рассчитанный по суррогатным данным. Оценки воздействия со стороны второго осциллятора на первый (У21) не превышают 95%-й уровень значимости. Это свидетельствует об отсутствии значимого взаимодействия.
На рис. 2, а представлена зависимость оценки времени задержки в связи от коэффициента связи между системами. Видно, что оценки времени задержки в связи имеют меньший разброс и практически не смещены, когда коэффициент связи принимает значения выше 0.4. На рис. 2, б представлена зависимость индекса фазовой синхронизации от коэффициента связи между системами. Из графика видно, что с увеличением коэффициента связи уровень синхронизации между осцилляторами растёт по линейному закону. Однако высокий уровень синхронизации не мешает получению правильной оценки времени задержки в связи (см. рис. 2, а).
3.1.2. Влияние расстройки частот на оценку времени задержки в связи. На рис. 2, в представлена зависимость оценки времени задержки в связи от расстройки частот. Видно, что большая расстройка частот приводит к ухудшению оценок времени задержки в связи, что отражено на графике значительными отклонениями от значения реального времени задержки в связи между системами. На рис. 2, г представлена зависимость индекса фазовой синхронизации от расстройки частот. Видно, что с увеличением расстройки значения индекса фазовой синхронизации немного уменьшаются.
3.1.3. Влияние шума на оценку времени задержки в связи. На рис. 2, д и е
представлены зависимости оценок времени задержки в связи от уровня шума при одинаковом коэффициенте связи и разной длине ряда. Из рис. 2, д видно, что график оценки времени задержки имеет пологий вид. Можно сделать вывод, что значение оценки времени задержки в связи слабо зависит от уровня шума в случае достаточно длинного временного ряда при большом уровне связи. Из рис. 2, е, построенного для случая коротких временных рядов, видно, что при слабых уровнях шума оценки времени задержки определяются хуже и имеют тенденцию к смещению и увеличению разброса. Кроме того, увеличивается риск получения ложного вывода о преимущественном направлении воздействия, поскольку для случая очень слабых шумов при сильной синхронизации можно получить значения оценок сил связи, превышающие 95%-й уровень значимости, рассчитанный с помощью суррогатных данных.
3.2. Два однонаправленно связанных осциллятора ван дер Поля. Рассмотрим зависимости оценки направления воздействия и времени запаздывания в связи от шума в случае сильной связи. При этом ответим на следующие вопросы:
Рис. 1. Зависимости оценок сил связи от пробного времени задержки при к = 0.7, Аист = 40, расстройке частот ô/ = 0.02, длине ряда 350 характерных периодов и дисперсии шума а = 0.8. Пунктирной линией показан 95%-й уровень значимости, рассчитанный по суррогатным рядам
Рис. 2. Зависимости от коэффициента силы связи оценки времени задержки (а) и индекса фазовой синхронизации (б) при а = 0.6, А = 40, расстройке частот б/ = 0.02 и длине ряда 350 характерных периодов. Зависимости от расстройки частот оценки времени задержки (в) и индекса фазовой синхронизации (г) при к =0.7, а = 0.3, А = 40, длине ряда 350 характерных периодов. Зависимости оценки времени задержки от уровня шума при к = 0.7, А = 40, б/ = 0.02, длине ряда 350 характерных периодов (д) и 70 (е). Вертикальными ограничителями показан разброс относительно среднего
— возможно ли получение несмещенных оценок запаздывания при большом уровне шума;
— существует ли зависимость смещения оценок задержки в связи от силы
связи;
— какой вклад вносит динамический и измерительный шум в динамику фазы осциллятора ван дер Поля?
Для ответа на эти вопросы рассчитывались зависимости оценок сил связи от пробного времени задержки по временному ряду, а также автокорреляционная функция и распределение свернутой фазы.
3.2.1. Зависимость оценки времени задержки в связи от длины ряда. Используемые в эксперименте временные ряды систем (5) были получены путём численного интегрирования уравнений методом Эйлера с шагом 0.01 при частотах /1,2 ~ 0.1. Участок, содержащий переходный процесс (50 характерных периодов), во временной ряд не записывался. Ошибка определения времени задержки в связи систем рассчитывалась следующим образом:
е = ^^ • 100%, (6)
где А - истинное время задержки в связи, а А - оценка, полученная из анализа временных рядов.
На рис. 3 представлена зависимость ошибки определения времени задержки в связи от длины ряда и уровня шума при большом значении коэффициента связи. Видно, что области, в которых ошибка мала, расположены в правой верхней части рисунка, то есть при высоких значениях уровня шума и длины ряда. Таким образом, увеличение длины ряда позволяет более точно определить время запаздывания в связи систем. Кроме того, более точно определить задержку в связи позволяет увеличение уровня шума.
3.2.2. Влияние измерительного шума на оценку времени задержки в связи. На рис. 4 показаны распределения циклической фазы ф временных рядов осцилляторов (5) при различных уровнях нелинейности и дисперсии измерительного шума.
350 300
ш
I 250 ^
о. ф
* 200 I
о. ф
I 150
го
X
100
50
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 о©
Рис. 3. Ошибка оценки времени задержки в зависимости от длины ряда и уровня шума при к = 0.7, А = 40
1200 -| 800 -400 -
0 -
а С 1200 -|
800 -
400 -
0 -
0
в
Рис. 4. Распределения циклической фазы при к = 0.9 и длине ряда 50 характерных периодов. Дисперсии измерительного шума а: 0.2 (а, б), 0.9 (в, г). Уровень нелинейности Х1>2: 0.1 (а, в), 1.0 (б,г)
Рис. 4, а и в, соответствующие фазам сигналов систем в случае слабой нелинейности (X = 0.1), демонстрируют более равномерный характер распределения ф и, как следствие, дают возможность более точно рассчитать фазу [15], чем в случае сильной нелинейности X = 1.0, рис. 4, б и г. Распределения ф, представленные на рис. 4, в и г (для а = 0.9) носят менее плавный характер, чем представленные на рис. 4, а и б (для а = 0.3), соответственно, что объясняется более высоким уровнем измерительного шума в системах. Однако существенной зависимости распределения ф от измерительного шума в широких пределах не наблюдается. Эти результаты указывают на то, что точность определения оценки запаздывания в связи в большей степени зависит от нелинейности взаимодействующих систем, чем от уровня шума.
Заключение
На практике, когда приходится иметь дело с сигналами различной природы, в случае сильной связи между системами, велика вероятность получения неточной оценки времени задержки в связи и её ошибочного толкования. В этом случае следует принимать во внимание зависимости оценки времени задержки от значения управляющих параметров исследуемых систем.
В работе на простых эталонных моделях связанных автоколебательных систем исследованы свойства оценки времени задержки в связи с помощью метода моделирования фазовой динамики для случая сильной связи. Показано, что при больших значениях коэффициента связи может быть получена точная оценка времени задержки в связи между системами. Однако при больших значениях расстройки частот, слабом уровне шума и малой длине ряда точность оценки снижается. Уменьшение уровня шума также увеличивает количество ложных выводов о преимущественном
направлении взаимодействия между системами. Соответственно, увеличение длины ряда, уровня шума и коэффициента связи способствует более точному определению направления и времени задержки в связи между осцилляторами.
Автор выражает благодарность В.И. Пономаренко и М.Д. Прохорову за плодотворные обсуждения.
Работа выполнена при поддержке грантов РФФИ № 12-02-00377 и 13-02-00227, гранта президента РФ для молодых ученых МК-4435.2012.8
Библиографический список
1. Arnhold J., Grassberger P., Lehnertz K., Elger C.E. A robust method for detecting interdependences: Application to intracranially recorded EEG // Physica D: Nonlinear Phenomena. 1999. Vol. 134. P. 419.
2. Quian Quiroga R., Kraskov A., Kreuz T., Grassberger P. Performance of different synchronization measures in real data: A case study on electroencephalographs signals // Phys. Rev. E. 2002. Vol. 65. 041903.
3. Smirnov D.A., Bodrov M.B., Perez Velazquez J.L., Wenneberg R.A., Bezruchko B.P. Estimation of coupling between oscillators from short time series via phase dynamics modeling: Limitations and application to EEG data // Chaos. 2005. Vol. 15. 024102.
4. Smirnov D.A., Andrzejak R.G. Detection of weak directional coupling: Phase dynamics approach versus state space approach // Phys. Rev. E. 2005. Vol. 71. 036207.
5. Пиковский А.С., Розенблюм М.Г., Куртс Ю. Синхронизация: фундаментальное нелинейное явление. М.: Техносфера, 2002.
6. Rosenblum M.G., Pikovsky A.S. Detecting direction of coupling in interacting oscillators // Phys. Rev. E. 2001. Vol. 64. 045202.
7. Smirnov D.A., Bezruchko B.P. Estimation of interaction strength and direction from short and noisy time series // Phys. Rev. E. 2003. Vol. 68. 046209.
8. Smirnov D.A., Barnikol U.B., Barnikol T.T., Bezruchko B.P., Hauptmann C., Buehr-le C., Maarouf M., Sturm V., Freund H.-J., Tass P.A. The generation of Parkinsonian tremor as revealed by directional coupling analysis // Europhys. Lett. 2008. Vol. 83. 20003.
9. Mokhov I.I., Smirnov D.A. El Nino Southern Oscillation drives North Atlantic Oscillation as revealed with nonlinear techniques from climatic indices // Geophysical Research Letters. 2006. Vol. 33. 024557.
10. Smirnov D.A., Bezruchko B.P. Detection of couplings in ensembles of stochastic oscillators // Phys. Rev. E. 2009. Vol. 79. 046204.
11. Смирнов Д.А., Карпеев И.А., Безручко Б.П. Выявление связи между осцилляторами по коротким временным рядам: условие применимости метода моделирования фазовой динамики // Письма в ЖТФ. 2007. Т. 33, вып. 4. С. 19.
12. Pikovsky A.S., Rosenblum M.G., Kurths J. Phase synchronization in regular and chaotic systems // Int. J. Bifurc. Chaos. 2000. Vol. 10. P. 2291.
13. Schäfer C., Rosenblum M.G., Abel H.-H., Kurths J. Synchronization in the human cardiorespiratory system //Phys. Rev. E. 1999. Vol. 60. P. 857.
14. Janson N.B., Bälanov A.G., Anishchenko VS., McClintock P.V.E. Coherence resonance versus synchronization in a periodically forced self-sustained system // Phys. Rev. E. 2002. Vol. 65. 036212.
15. Kralemann B., Cimponeriu L., Rosenblum M.G., Pikovsky A.S., and Mrowkä R. Phase dynamics of coupled oscillators reconstructed from data // Physical Rev E. 2008. Vol. 77. 066205.
Саратовский государственный Поступила в редакцию 12.11.2012
университет им. Н.Г. Чернышевского
ESTIMATION OF INTERACTION DIRECTION BETWEEN OSCILLATORY MODEL SYSTEMS IN CASE OF CLOSE COUPLING
V.S. Khorev
The task of detection statistically significant interaction, its direction and delay between time data series of two oscillatory systems in case of close coupling is investigated with nonlinear modeling approach. Numerical experiments on oscillatory model systems with different coupling function variants are used to study main dependences.
Keywords: Close coupling, interaction, coupling direction, time delay, oscillatory systems.
Хорев Владимир Срегеевич - родился в Гомеле (1987), окончил магистратуру в Саратовском государственном университете им. Н.Г. Чернышевского по специальности «биомедицинская инженерия» (2011). Аспирант кафедры динамического моделирования и биомедицинской инженерии СГУ. Принимал участие в 10 международных и всероссийских конференциях. Область научных интересов - нелинейная динамика, анализ временных рядов, математическое моделирование, приложение методов нелинейного анализа к задачам биофизики и медицины.
410012 Саратов, ул. Астраханская, 83
Саратовский государственный университет им. Н.Г. Чернышевского E-mail: [email protected]
