CC BY
4
Прикладные задачи
нелинейной теории колебаний и волн
Известия высших учебных заведений. Прикладная нелинейная динамика. 2022. Т. 30, № 4 Izvestiya Vysshikh Uchebnykh Zavedeniy. Applied Nonlinear Dynamics. 2022;30(4)
Научная статья УДК 530.182
DOI: 10.18500/0869-6632-2022-30-4-391-410
EDN: BPJAOD
Реконструкция интегрированных уравнений системы фазовой автоподстройки частоты под периодическим внешним воздействием по скалярному временному ряду
М.В. Сысоева1,2 М.В. Корнилов1'3'4, Л. В. Такаишвили3, В. В. Матросов4, И. В. Сысоев1'3'4
1 Саратовский филиал Института радиотехники и электроники им. В. А. Котельникова РАН, Россия 2 Саратовский государственный технический университет имени Гагарина Ю. А., Россия 3 Саратовский национальный исследовательский государственный университет имени Н.Г. Чернышевского, Россия 4Национальный исследовательский Нижегородский государственный университет им. Н. И. Лобачевского, Россия E-mail: Elbobrichek@mail.ru, kornilovmv@gmail.com, nar7187@yandex.ru, matrosov@rf.unn.ru, ivssci@gmail.com
Аннотация. Цель настоящего исследования — разработка методики реконструкции уравнений системы фазовой автоподстройки частоты, находящейся под периодическим внешним воздействием, по скалярному временному ряду одной переменной. Методы. Вместо исходной модели реконструируется модель, интегрированная по времени, что позволяет существенно снизить чувствительность метода к шумам наблюдения, поскольку не требуется оценивать вторую производную наблюдаемой численно. Внешнее периодическое воздействие аппроксимируется тригонометрическим полиномом от времени, интеграл от которого также представляет собою тригонометрический полином. Допущение о непрерывности неизвестной нелинейной функции используется для построения целевой функции и оптимизации. Результаты. Показано, что предложенный подход даёт существенное преимущество над ранее разработанным подходом к реконструкции неинтегрированных уравнений, позволяя добиться приемлемых оценок параметров при измерительном шуме порядка 10% от среднеквадратичного отклонения сигнала даже при наличии внешнего воздействия. Заключение. Описанный подход существенно расширяет возможности реконструкции систем фазовой автоподстройки частоты, позволяя реконструировать системы под произвольным периодическим воздействием и при этом существенно увеличивая устойчивость к шуму.
Ключевые слова: реконструкция, система фазовой автоподстройки частоты, периодическое внешнее воздействие.
Благодарности. Работа выполнена при поддержке Совета по грантам Президента Российской Федерации, грант МД-3006.2021.1.2 и Российского фонда фундаментальных исследований, грант № 19-02-00071.
Для цитирования: Сысоева М.В., Корнилов М.В., Такаишвили Л. В., Матросов В. В., Сысоев И. В. Реконструкция интегрированных уравнений системы фазовой автоподстройки частоты под периодическим внешним воздействием по скалярному временному ряду//Известия вузов. ПНД. 2022. T. 30, № 4. С. 391-410. DOI: 10.18500/0869-6632-2022-30-4391-410. EDN: BPJAOD
Статья опубликована на условиях Creative Commons Attribution License (CC-BY 4.0).
Поступила в редакцию 22.10.2021, принята к публикации 14.12.2021, опубликована 1.08.2022
Article
DOI: 10.18500/0869-6632-2022-30-4-391-410
Reconstruction of integrated equations of periodically driven phase-locked loop system from scalar time series
M.V. Sysoeval'2m, M. V Kornilov1'3'4, L. V Takaishvili3, V V Matrosov4, I. V Sysoev1'3'4
1 Saratov Brunch of Kotelnikov Institute of Radioengineering and Electronics of RAS, Russia 2Yuri Gagarin State Technical University of Saratov, Russia
3 Saratov State University, Russia 4Lobachevsky State University of Nizhny Novgorod, Russia E-mail: Elbobrichek@mail.ru, kornilovmv@gmail.com, nar7187@yandex.ru, matrosov@rf.unn.ru, ivssci@gmail.com Received 22.10.2021, accepted 14.12.2021, published 1.08.2022
Abstract. Purpose of this work is to develop a reconstruction technique for the equations of a phase-locked loop system under periodic external driving from a scalar time series of one variable. Methods. Instead of the original model, we reconstructed a time-integrated model. So, since it is not necessary to evaluate the second derivative of the observable numerically, the method sensitivity to observation noise has significantly decreased. The external periodic driving is approximated with a trigonometric polynomial of time, the antiderivative of which is also a trigonometric polynomial. The assumption about continuity of an unknown nonlinear function is used to construct the target function for optimization. Results. It is shown that the proposed approach gives a significant advantage over the previously developed approach to the reconstruction of non-integrated equations, allowing to achieve acceptable parameter estimates with measurement noise being about 10% of the RMS deviation of the signal even in the presence of external driving. Conclusion. The described approach significantly extends the possibilities of reconstruction of phase-locked loop systems, allowing systems to be reconstructed under arbitrary periodic driving and at the same time significantly increasing noise resistance.
Keywords: reconstruction, phase-locked loop system, periodic driving.
Acknowledgements. This study was supported in part by the President of the Russian Federation under Grant MD-3006.2021.1.2, and in part by the Russian Foundation for Basic Research, grant No. 19-02-00071.
For citation: Sysoeva MV, Kornilov MV, Takaishvili LV, Matrosov VV, Sysoev IV. Reconstruction of integrated equations of periodically driven phase-locked loop system from scalar time series. Izvestiya VUZ. Applied Nonlinear Dynamics. 2022;30(4):391-410. DOI: 10.18500/0869-6632-2022-30-4-391-410
This is an open access article distributed under the terms of Creative Commons Attribution License (CC-BY 4.0).
Введение
Теоретические возможности методов, основанных на реконструкции уравнений по временным рядам, очень велики [1]: непрямое измерение величин и даже целых функций, которые нельзя измерить напрямую; кластеризация объектов по параметрам модели; верификация различных моделей одного и того же объекта или явления; прогноз будущего поведения во времени или прогноз поведения при небольшом изменении параметров (прогноз бифуркаций). В приложении к биологическим нейронам все эти возможности, несомненно, востребованы, так как прямые измерения большинства модельных параметров невозможны, эксперименты сами по себе носят травмирующий для клетки характер, а вот опознавать клетки и по их сигналам уметь определять, к какому типу принадлежит та или иная клетка, очень важно. Однако прогресс в данной области за последние годы не очень велик, что обусловлено целым рядом факторов, главные из которых несовершенство моделей и шумы и помехи во время измерения. К настоящему времени ни для одной из известных математических моделей нейронов: Ходжкина-Хаксли [2], ФитцХью-Нагумо [3,4], Хиндмарш-Роуз [5], Моррис-Лекара [6] не удалось реконструировать уравнения по экспериментальным временным рядам клеточной активности реальных нейронов. Одна из главных причин этого — необходимость либо измерять все переменные модели в эксперименте,
что невозможно, поскольку в моделях часть из них носит условный характер. Фактически, единственная измеримая переменная — трансмембранный потенциал. Для получения оставшихся переменных (скрытых) есть несколько подходов. Либо они получаются численным дифференцированием или интегрированием, как предлагалось в [7], либо методом временных задержек, как в [8], либо необходимо задавать стартовые догадки для них и затем подгонять сами ряды вместе с неизвестными параметрами и нелинейными функциями к наблюдаемым (такой подход получил название методов работы со скрытыми переменными [9]).
Метод временных задержек подходит в двух случаях: когда уравнения для объекта совсем неизвестны и годится любой подход к реконструкции вектора состояния, либо при реконструкции систем с запаздыванием в случае, если запаздывающая переменная — единственная (если уравнение второго порядка, приходится прибегать к комбинации задержек и численного дифференцирования [10]). Для реконструкции уравнений нейронов, написанных из биофизических соображений, этот метод не годится, поскольку получить недостающие переменные сдвигом нельзя — они имеют другой смысл и другую размерность. Численное дифференцирование и интегрирование имеют ровно те же проблемы: либо используется, когда об уравнениях для объекта ничего не известно [7], либо в случае, когда нам повезло и уравнения составлены так, что одна из переменных имеет, например, физический смысл координаты, а другая — скорости (заряда — тока и т.п.). Для моделей нейронов [2,3,5,6] это неприменимо в силу их устройства, но некоторые другие модели биологических систем позволяют реконструировать часть переменных по наблюдаемой таким способом, например, модели голосовых связок [11] или нефрона [12]. Наиболее значимые успехи в данной области связаны с использованием статистических оценок [13].
Кроме моделей нейронов, построенных из биофизических соображений, есть ещё модели генераторов, способные демонстрировать широкое множество характерных для активности нейронов режимов поведения. Один из таких генераторов был изначально предложен для целей связи и представляет собою систему фазовой автоподстройки частоты (ФАПЧ) с полосовым фильтром [14]. Режимы нейроподобной генерации, включая периодические и хаотические пачки импульсов (бёрсты) были обнаружены в модели в работах [15,16], а затем продемонстрированы экспериментально [17]. Система ФАПЧ не является моделью нейрона, поэтому верифицировать её по данным активности биологического нейрона или оценить с её помощью физиологически значимые параметры нейрона не получится. Но она может быть полезна для кластеризации нейронов и для детектирования связей между ними, что уже достаточно для того, чтобы решение задачи о реконструкции этой системы по экспериментальным данным было оправдано. Её преимущество состоит в том, что для этой системы одна из двух недостающих переменных может быть получена численным дифференцированием, а вторая — численным интегрированием наблюдаемой.
Ранее на основе подхода, предложенного в [18] для осцилляторов первого порядка с запаздыванием в собственной динамике, в работе [19] была построена методика реконструкции одной автономной системы ФАПЧ с полосовым фильтром, также имеющей запаздывание по переменной у. Далее аналогичный подход был распространён на ансамбль связанных ФАПЧ [20] уже без запаздывания. Однако для практического применения подхода представленных в [19,20] наработок недостаточно, поскольку они имеют ряд существенных недостатков. Во-первых, предложенная методика довольно чувствительна к шуму: в численном эксперименте удалось реконструировать динамику в ансамбле только при наличии аддитивного шума порядка 1%. В реальном эксперименте, например, при внеклеточных записях нейронов уровень шума гораздо выше [21]. К сожалению, несмотря на применение сглаживающих полиномов для численного дифференцирования необходимость получения второй производной наблюдаемой не даёт возможности добиться работоспособности метода при заметно больших шумах. Во-вторых, предложенная в [19,20] методика оперирует при оценке нелинейной функции слагаемыми типа у/у или ь/у(Ь — т), где V — это некоторое известное выражение, зависящее от наблюдаемой у и её первой и второй
производной, а аргумент t — т означает, что в данном случае используется задержанное во времени значение наблюдаемой. Такого рода слагаемые дают сингулярность при у ^ 0, а поскольку в эксперименте у измерено с ошибками, вызванными шумом, все результаты для |у| < р,, где величину ^ приходится подбирать эмпирически, оказываются ненадёжными. К сожалению, малые значения у характерны для достаточно длительных отрезков времени, особенно в возбудимом режиме. Как следствие этого, большая часть временного ряда исключается из рассмотрения, а для части режимов нелинейную функцию вовсе удаётся реконструировать только частично. В-третьих, предложенная в [19] методика работает только для уединённой системы, что при реконструкции по биологическим данным почти бессмысленно: отдельные клетки, как правило, молчат и ничего не генерируют, а при реконструкции связанных, используя методику из [20], необходимо иметь ряды всех.
Данная работа имеет целью расширить область применения методики реконструкции ФАПЧ с полосовым фильтром для целей её дальнейшего применения к экспериментальным временным рядам реальных нейронов. Чтобы устранить приведённые выше недостатки, осуществлён переход от исходных уравнений, описанных в работах [15,16], к уравнениям, интегрированным по времени, а внешнее воздействие будем аппроксимировать тригонометрическим полиномом степени Р, как это было предложено для систем другого типа в работе [22].
1. Математическая модель системы ФАПЧ
Динамика рассматриваемого здесь фазоуправляемого генератора при наличии внешнего воздействия Iext (t) описывается неавтономной системой обыкновенных дифференциальных уравнений третьего порядка, определенной в цилиндрическом фазовом пространстве (ф mod 2п, у, z):
у,,
z, (1)
Y + Iext(t) — (£\ + £2)Z — (1 + £\ COS ф)у.
В терминах системы фазовой автоподстройки частоты ф — текущая разность фаз подстраиваемого и опорного генератора, у — начальная частотная расстройка, £1,£2 — параметры инерционности фильтров. Применительно к динамике нейрона переменную у можно интерпретировать как описывающую изменение мембранного потенциала, параметры е\ и позволяют задавать необходимый динамический режим, у — внешний ионный ток в межклеточном веществе постоянный на времени измерения (медленно меняющаяся компонента [2]), Iext — периодическое воздействие другого нейрона или внешняя электрическая стимуляция. Если Iext имеет ненулевую постоянную компоненту (сдвиг), она фактически учитывается как часть параметра у.
В [23] было показано, что данный генератор способен демонстрировать различные динамические режимы, характерные для нейронов: регулярную импульсную динамику (спайки), пачечные колебания (бёрсты) с различным числом импульсов в пачке. В работе [15] было произведено разбиение плоскости параметров модели на области существования соответствующих динамических режимов: квазирегулярные колебания, колебания с различным количеством импульсов в пачке, хаотические колебания. В работе [24] подробно исследована модель рассматриваемого генератора, находящегося в возбудимом состоянии, в которой отклики появляются только при подаче внешнего воздействия. Наличие возбудимого режима в рассматриваемой модели является очень важным, поскольку нейроны мозга большую часть времени находятся именно в возбудимом состоянии. Показано, что при у = 0 и Iext(t) = 0 в системе (1) существует континуум состояний
£1 £2
dф ~dt dy dt dz dt
равновесия (ф*, 0, 0), расположенных на отрезке ф* е [— агееов^/^); агееов^/^)], в одном из которых в конечном счёте оказывается система (1) при ненулевых начальных значениях у и г. Заставить систему генерировать импульсы при у = 0 можно с помощью подачи внешнего импульсного воздействия 1ех±(£) — гармоническое воздействие в такой системе при таких параметрах приводит к вынужденным гармоническим колебаниям. В [25] рассмотрено влияние параметров периодической стимуляции на ответ исследуемого генератора. Различные способы оценки откликов генератора на внешнее возбуждение показывали, что отклик существенно зависит от амплитуды стимулирующих импульсов и слабее зависит от периода их следования; это подтверждается результатами действия пуассоновской случайной последовательности.
2. Методика реконструкции интегрированных уравнений системы ФАПЧ
В аппаратной реализации системы ФАПЧ [17] возможно измерение единственной переменной у. Она же соответствует трансмембранному потенциалу, регистрируемому в биологическом эксперименте от отдельных нейронов [26,27]. Поэтому задача реконструкции уравнений (1) сразу была сформулирована таким образом, что измеряется только переменная у. Следуя работе [19], переменная ф восстанавливалась численным интегрированием методом Симпсона, переменная z — численным дифференцированием с использованием аппроксимирующего по т точкам полинома (фильтр Савицки-Голэя [28]). Для реконструкции по незашумлённым или слабозашумлённым реализациям подойдёт самое малое т = 3, при наличии шума т может быть подобрано, в том числе используя в качестве критерия качество реконструкции нелинейной функции, как это было сделано в [19]. Метод, предложенный в [19,20], требовал также численного получения второй производной наблюдаемой — временного ряда величины dz/dt, что являлось основною причиною низкой его устойчивости к шуму. Алгоритм, предлагаемый далее, не требует численного расчёта второй производной. Переменная ф может быть восстановлена численным интегрированием методом Симпсона. Поскольку метод использует построение парабол на интервалах, равных двум временам выборки, далее всегда будет рассматриваться нечётное количество точек в измеренном ряде (N): при необходимости всегда можно отбросить одно значение без существенной потери информации.
Чтобы повысить общность и применимость метода к экспериментальным данным, нужно сделать реконструируемую модель менее «ломкой», более общей. Для этого последнее уравнение системы (1) перепишем в виде (2), где функция f не обязана соответствовать формуле (3) (в реальной системе она почти наверняка не соответствует), то есть метод реконструкции будем писать сразу для произвольной непрерывной функции f.
dz Y 1 , , £1 + £2 t( , ^
~T7 = - + -Iext{t)--2 - f (ф)y, (2)
at £i£2 £i^2 £i^2
fi \ 1+ £1 COS Ф (3) f (Ф) = -—-. (3)
£1£2
Далее последнее уравнение отдельного осциллятора (2) проинтегрируем по времени и введём дополнительные обозначения ao и a1 аналогично [20]:
2 = aot + a1y -J f (ф) ^ф-dt + J I'ext(t)dt, (4)
Y £1 + £2 T, hxt
a0 = -, a1 =--, Iext = -, (5)
£1£2 £1£2 £1£2
где мы воспользовались тем, что J zdt = у, J (dz/dt) = z. Далее введём в рассмотрение новую
нелинейную функцию Ф(ф), которая, очевидно, будет гладкою, поскольку по определению её первая производная — /(ф), для которой мы условились, что она непрерывна.
I /(ф) ^ = 1 /(ф)йф = Ф(ф). (6)
Также отметим, что константа интегрирования в уравнении (4) может быть без потери общности внесена в Ф(ф).
Следуя идеям работы [18], будем строить целевую функцию для вычисления коэффициентов ао и а1, исходя из минимизации длины описания нелинейной функции Ф. Это позволит отказаться от явного разложения функции / и её первообразной Ф в ряд и, таким образом, одновременно повысит общность метода (он может использоваться для произвольной /) и уменьшит параметризацию (число подлежащих оценке параметров а). Это в свою очередь улучшит статистические свойства оценок оставшихся в модели коэффициентов: чем меньше параметров, тем лучше их оценки. Внешнее воздействие, следуя [22], представим в виде тригонометрического полинома, где частота □ = 2л/ТехЛ известна. Вопрос о том, как её подобрать, разберём далее. Таким образом, базовым для нас становится следующее выражение для нелинейной функции Ф, полученное из уравнения (4):
к=К
Ф(ф) = аоЬ + а1у + {а2к соя(кО,1) + а2к+1 8т(кО,Ь)) — г, (7)
к=1
где К — степень тригонометрического полинома. Постоянное слагаемое в формуле (7) не нужно: оно, как и константа интегрирования, включено в Ф; то есть фактически Ф может быть определена с точностью до константы, но это и естественно, поскольку Ф есть первообразная /.
Далее все выкладки будем вести, следуя работе [20], но имея в виду новое выражение для нелинейной функции (7). Введём сортирующее отображение Q(n), где п — номер точки в ряде, ставящее в соответствие п-й точке в исходном ряде её номер Q(n) в отсортированном по возрастанию ф. Рассмотрим также обратное отображение Q-1, вычисляющее по номеру в отсортированном ряде номер в исходном, так что (^-1((^(п)) = п. Рассмотрим точку, находящуюся в отсортированном ряде непосредственно перед (^(п)-й, тогда в исходном ряде она имеет номер Я-1(Жп) — 1), который обозначим как рп для краткости. Тогда приращение Ьп функции Ф на отрезке [ф(рга); ф(^)] будет выражаться формулою (8):
= Ф(ф(п)) — Ф(ф(рп)) = аоАг(п) + а1Ау(п) — Аг(п) +
к=К
+ (°-2кА&к + а2к+1АЕк), к=1
Аг(п) = г(п) — г('Рп),
Ау(п) = у(п) — у('Рп), (8)
Аг(п) = г(п) — г(рп),
А&к = сов(к01(п)) — соё(к01(рп)),
АЕ& = ё'т(к01(п)) — ё'т(к01(рп)).
При этом значение п такое, что Q(n) = 0 (будем здесь и далее считать, что нумерация начинается с нуля, то есть п = 0,1,... — 1) не позволено, поскольку для него нет соответствующего рп (нет предыдущей точки в отсортированном ряде — наше значение самое маленькое). Тогда в качестве целевой функции рассмотрим величину Ь:
Ь(ао, аь..., а2к+1) = ^2 . (9)
п
Из (8) очевидно, что Ь будет зависеть от 1к=0 2к+1 квадратично. Значит, формулу (9) можно рассматривать как формулировку задачи на наименьшие квадраты для аппроксимации величин А 2 (п). Выражения при коэффициентах суть базисные функции (вернее, их значения, доступные в моменты измерений), а Ьп соответствуют невязкам. При верном выборе а^ целевая функция Ь будет много меньше, чем при неверном, когда коэффициенты определены неверно и функция / имеет разрывы почти в каждой точке.
Поскольку даже при верном выборе в общем случае Ь > 0, то при конечном N полученные таким образом оценки должны быть смещёнными — их математическое ожидание не будет в точности совпадать с истинными значениями. Так как функция Ф не только непрерывна, но и дифференцируема (по определению она является первообразною непрерывной функции /) на полуинтервале 0 ^ ф < 2п, то 8га ^ (Ф(ф>('рп)), то 8га является аппроксимацией дифференциала функции Ф в точке ф('га) справа, или в точке ф(п) слева, а для дифференцируемой функции в пределе это одно и то же. Тогда должно выполняться следующее соотношение:
2п 2п
^=/«2=/ш <«2 = (£)'
0 0
2п
0. (10)
Таким образом, оценки предложенным методом являются асимптотически несмещёнными.
Предложенный подход позволяет оценить значения ао и ai, являющиеся комбинациями исходных параметров у, ei и £2, но не сами эти параметры. Следует отметить, что если рассмотреть третье уравнение системы (1), то можно видеть, что при замене 1 + £1 cos ф = f (ф) оценить отдельно у, £1 и е2 будет нельзя, потому что все члены уравнения, включая правую часть, стоят с неизвестными коэффициентами (или функциями). Таким образом, система уравнений для определения коэффициентов, которую можно было бы составить, оказалась бы вырождена (нет свободного члена). Это значит, что и никаким другим способом оценить у, £1 и £2 отдельно по данным нельзя, то есть данное ограничение является свойством системы, а не недостатком метода.
Изложенный алгоритм предполагает знание периода внешнего воздействия. В действительности, даже небольшая ошибка в его значении сильно повлияет на результаты реконструкции, как это было показано в [29]. Самое простое и эффективное решение данной проблемы при современных вычислительных средствах — перебор значений пробного периода воздействия Т в некотором, довольно большом диапазоне с шагом h порядка шага выборки, хотя целочисленное соотношение между Т и h не обязательно. По минимуму построенной зависимости L(T) можно достаточно точно оценить истинный период. Полагаясь на результаты [22], можно рассчитывать, что высокая чувствительность целевой функции к периоду воздействия в типичном случае обеспечит верное определение периода даже при невысоких значениях степени тригонометрического полинома Р и недостаточном сглаживании m в борьбе с шумами измерений.
Предложенный подход позволяет также реконструировать внешнее воздействие на основе оценок коэффициентов a^|fc=2 2^+1 с точностью до множителя 1/(е 1е2), который в общем случае не может быть определён. Для этого нужно продифференцировать входящий в (7) тригонометрический полином:
к=К
I'ext(t) = ^ ( -kQa2k sin(kQt) + kQa2k+1 cos(kQt)). (11)
k=1
Поскольку значения коэффициентов |fc=2 2K+1 и значение частоты Q известны либо оценены, внешнее воздействие по формуле (11) можно рассчитать для любого момента времени t.
После получения значений можно таблично рассчитать по формуле (7) функцию Ф. Таблично, потому что в формулу (7) входят значения у и z, а эти переменные доступны только
0
в моменты измерения. Но сама по себе функция Ф, скорее всего, не очень интересна, для реконструкции модели интересна её первообразная /(ф). Как аналитически проинтегрировать (7) по ф — не ясно, поэтому проще выразить / из формулы (2), подставив туда выражения для ао и а1 из (4):
/(ф) = (ао + а^ + 1'ех^) — У. (12)
Помимо того, что расчёт по формуле (12) возможен только для ф, соответствующих наблюдаемым у (иначе не будет значений у, г и йг/А£), имеются два ограничения, которых лишён весь остальной алгоритм. Эти ограничения были присущи исходному методу, предложенному в работе [19]: во-первых, необходимо численно рассчитать вторую производную наблюдаемой, поскольку (Ы/(И = = й'2у/(И'2"; во-вторых, при у ~ 0 рассчитанные значения будут сильно неточными, особенно при наличии шумов. При этом следует отметить, что даже несмотря на эти ограничения, полученные оценки будут заметно точнее, чем при использовании оригинального подхода из [19], поскольку коэффициенты ак были рассчитаны, не опираясь на оценку (Ы/(И и без ограничений для у « 0.
3. Результаты
3.1. Генерация рядов и восстановление вектора состояния. Уравнения (1) решались методом Эйлера. Интервал дискретизации А£ = 1 /32, равный шагу интегрирования, подбирался эмпирически таким образом, чтобы получить устойчивое решение. Подробное исследование этого вопроса было проведено в [30].
Для тестирования метода система ФАПЧ рассматривалась в двух режимах: возбудимом (у = 0,£1 = 4,е2 = 10) и колебательном (у = 0.075, е1 = 4.5, е2 = 10 — спайковый режим согласно классификации, данной в [20]). Подавались три вида воздействия: прямоугольные импульсы, гауссовские импульсы и гармоническое воздействие: импульсное воздействие обоих типов в обоих режимах (см. рис. 1, c, e), а гармоническое - только для колебательного режима, см. рис. 1, a. Выбор типов воздействия обусловлен тем, что гармоническое воздействие легче всего для описания и на нём проще всего тестировать метод, гауссовские импульсы требуют для описания тригонометрических полиномов высоких порядков, но, будучи гладкими, они всё же могут быть удовлетворительно описаны, а прямоугольные импульсы типичны в радиотехнике и рассматривались для данной системы ФАПЧ в работе [25]. Если система находится в подпо-роговом режиме, под гармоническим воздействием она демонстрирует вынужденные линейные колебания — режим как бессмысленный с точки зрения нейродинамики, так и неудобный для реконструкции, поскольку содержит слишком мало информации.
В соответствии с [25] введём обозначения: А8р — амплитуда стимулирующего импульса, Тер — длительность входного импульса, Т8р — период стимуляции. Чтобы сохранить одинаковую мощность внешнего воздействия для разных его форм, были использованы следующие параметры:
1. для прямоугольных импульсов А8р = 0.26, Т8р = 100, х8р = 10, скважность = Т8р/х8р = 10;
2. для гауссовских импульсов А^р = Л^ = 0.26, Тёр = Т8р = 100, хёр = т8р/л/гс;
3. для гармонического воздействия А^ = Asp\/2/Qsp, Т^ = Твр = 100.
Ко всем временным рядам добавлялся измерительный шум со среднеквадратичным отклонением в 10% от среднеквадратичного отклонения сигнала, на рис. 1 временные ряды и фазовые портреты приведены без шума наблюдений. Реконструкция проводилась как по рядам с шумом, так и по рядам без шума, но далее приводятся результаты только только по рядам с шумом, если это не оговорено дополнительно, поскольку метод ориентирован на работу с зашумлёнными данными.
0.05 0.00 -0.05 -0.10
1000
У z
Pulses
0.04 0.02 0.00 N -0.02 -0.04 -0.06 -0.08
1000
Рис. 1. Временные ряды переменных у (наблюдаемая) и г модели (1) и внешнего воздействия (а, с, е). Фазовые портреты в проекции (у, z) системы (1) (b, d, f). Части а, Ъ соответствуют гармоническому воздействию на систему, находящуюся в колебательном режиме, части с, d — воздействию прямоугольными импульсами на систему в колебательном режиме и части е, f — воздействию прямоугольными импульсами на систему в возбудимом режиме (цвет online)
Fig. 1. Time series of variables у ( observable) and г from the model (1) and external driving (a, c, e). Phase portraits in the projection (y, z) of the system (1) (b, d,f). Parts a, b correspond to the system under harmonic driving in the oscillatory regime. Parts c, d correspond to the system under square pulse driving in the oscillatory regime. Parts e, f correspond to the system under square pulse driving in the excitable regime (color online)
1000
0.050 0.025 0.000
N
-0.025 -0.050 -0.075 -0.100
Ь
3.2. Реконструкция параметров ФАПЧ. Оценить качество реконструкции проще всего, сопоставив восстановленные значения эффективных параметров ао и щ со значениями, рассчитанными аналитически по формуле (5). Для ао при импульсном воздействии введена дополнительная поправка, вызванная ненулевым средним внешнего воздействия, которое компенсируется эффективным изменением параметра у. Результаты реконструкции для всех пяти режимов представлены на рис. 2. Видно, что в колебательном режиме для всех трёх типов воздействия оценки а! оказываются довольно точными: их относительная погрешность не превышает 2%, причём для импульсного воздействия такой точности удаётся добиться уже при К = 3 — дальнейшее увеличение числа гармоник не приводит к уточнению оценок. Для возбудимого
« 3 о
|Н о в 2
1-
-0.28 --0.30 --0.32 --0.34 -
0_G
0_Н
0_S
S_G
S_S
Situe
Отгие
10
15 P
20
25
30
10
15 P
20
25
30
Рис. 2. Зависимости восстановленных значений коэффициентов ао и ai от степени использованного тригонометрического полинома К. Три кривые, помеченные как О (OjG для гауссовых импульсов, 0_Н для гармонического воздействия и 0_S для прямоугольных импульсов), относятся к колебательному режиму, истинное значение для них обозначено горизонтальною линией Отше- Две кривые, помеченные как S (S_G для гауссовых импульсов, S_S для прямоугольных), относятся к возбудимому режиму, для них истинное значение обозначено как Зтгие- Величина ао для удобства дана умноженною на 103 (цвет online)
Fig. 2. The dependences of the reconstructed values of the coefficients ao and ai on the order of the used trigonometric polynomial K. The three curves marked as О (0_G for Gauss pulses, 0_H for harmonic driving and 0_S for square pulses) belong to the oscillatory regime, the true value for them is indicated by the horizontal line Отше- Two curves marked as S (S_G for Gauss pulses, S_S for square ones) belong to the excitable regime, for them the true value is denoted as Зтгие-The value of ao is given multiplied by 103 for the convenience (color online)
режима относительная ошибка оценок а! несколько больше и составляет порядка 4%, при этом лучшие результаты достигаются при использовании К = 5.
Ошибки оценок ао в целом существенно выше, чем аь что может быть вызвано различными причинами. Одна из них заключается в том, что в формуле (7) этот коэффициент стоит при времени, а оно не является ограниченною переменною в отличие от у, г и ф. Неограниченный линейный рост I с увеличением длины ряда, вообще говоря, является проблемою при численном счёте, ухудшая статистические характеристики в первую очередь оценки ао. В работе [19] основным источником ошибок была производная г!г/сИ, в модифицированном алгоритме, изложенном в данной работе, нет необходимости расчёта второй производной, но в связи с интегрированием исходного уравнения (1) появился член, содержащий линейную зависимость от времени, так что ошибки определения ао в некотором смысле можно считать платою или компенсацией за отсутствие необходимости двойного численного дифференцирования. Ещё одна возможная причина невысокой точности восстановления параметра ао — малость величины у в рассмотренных режимах, что автоматически означает и малость её вклада в общую динамику, и высокие относительные погрешности оценок при малых абсолютных.
3.3. Реконструкция внешнего воздействия. Вторым важным критерием успешности процедуры реконструкции является восстановление внешнего воздействия при использовании тригонометрических полиномов различной степени. Результаты для колебательного режима при некоторых использованных значениях степени тригонометрического полинома К приведены на рис. 3. В численных экспериментах использовались и другие: большие и меньшие значения К, результаты для которых не отображены на рис. 3, чтобы не загромождать его.
Гармоническое воздействие хорошо восстанавливается при любом К, что естественно. Существенное увеличение К относительно оптимального значения К = 1 ведёт лишь к небольшим искажениям синусоиды, что хорошо видно на рис. 3, а. Для гауссовых импульсов К = 15 позволяет добиться приличной аппроксимации — такой, что визуально осцилляции внутри
периода, неизбежно присутствующие при аппроксимации конечным тригонометрическим полиномом (см. рис. 3, Ъ — оранжевая кривая, соответствующая К = 5), становятся незаметны. Точная аппроксимация прямоугольных импульсов теоретически достижима только в пределе, учитывая,
2 -
o -
100
2oo
300
100
b
t
K=5
2oo
300 K=25
100
2oo
300
Рис. 3. Реконструированное нормированное внешнее воздействие I'xt для колебательного режима: a — импульсное гауссово, b — импульсное прямоугольное, c — гармоническое. Различные кривые на каждом рисунке соответствуют разным степеням тригонометрического полинома К (цвет online)
Fig. 3. Reconstructed normalized external driving I'xt for the oscillatory regime: a — Gauss pulses, b — square pulses, c — harmonic driving. The different curves in each subplot correspond to different orders of the trigonometric polynomial К (color online)
o
o
o
a
c
что они имеют вертикальный фронт, а значит, содержат все частоты, включая бесконечно большие (или по крайней мере столь большие, что они неразрешимы при текущей частоте выборки). Рис. 3, c содержит несколько различных вариантов: вариант при К = 5 очень визуально похож на реконструированное при том же К = 5 воздействие гауссовыми импульсами. Реконструкция при К = 15 имеет больше локальных минимумов на периоде, что логично, а сама амплитуда колебаний меньше; на самом импульсе уже два экстремума, а фронт заметно круче, что отличает его от аппроксимации гауссовых импульсов. Реконструкция при К = 25 имеет уже три максимума на самом импульсе и ещё более крутой фронт.
Результаты реконструкции воздействия для возбудимого режима (рис. 4) качественно не отличаются от того, что представлено на рис. 3, Ь, c для колебательного режима. В целом можно констатировать, что предложенный подход позволяет довольно точно восстановить форму внешнего периодического воздействия, хотя для некоторых типов требует очень высокой степени тригонометрического полинома. Интересно, что точность аппроксимации воздействия слабо влияет на точность реконструкции параметров ФАПЧ ао и а^ Например, по рядам ФАПЧ под воздействием прямоугольных импульсов при К = 5 значения параметров восстанавливаются не хуже, чем при больших значениях, но само воздействие как прямоугольное идентифицировать нельзя.
б -4 -
Gauss Pulses
5c
icc
i5c t
K=5
200 K=15
250
300
b
Square Pulses
50
100 K=1
150 t
- K=5
200
250 300 K=25
Рис. 4. Реконструированное нормированное внешнее воздействие I'xt для возбудимого режима: a — импульсное гауссово, b — импульсное прямоугольное. Различные кривые на каждом рисунке соответствуют разным степеням тригонометрического полинома К (цвет online)
Fig. 4. Reconstructed normalized external driving I'xt for the excitable regime under different driving types: a — Gauss pulses, b — square pulses. The different curves in each subplot correspond to different orders of the trigonometric polynomial К (color online)
0
0
a
3.4. Реконструкция нелинейной функции. Третий критерий успешности реконструкции — идентификация нелинейной функции / (ф), которую по рассчитанным коэффициентам ак можно рассчитать по формуле (12). При этом надо понимать, что реконструкция нелинейной функции / — относительно слабое место алгоритма, поскольку формула (12) предполагает деление на у, а значения у ^ 0 часто встречаются, особенно для возбудимого режима.
Результаты реконструкции нелинейной функции для колебательного режима представлены на рис. 5, результаты для всех трёх типов воздействия приведены для К = 15. Разное число периодов реконструированного синуса объясняется различиями динамических режимов из-за наличия воздействия — фаза успевает вырасти на 2п разное число раз за одно и то же время наблюдения. Вертикальные «шумные» линии из точек соответствуют векторам состояния, когда у ~ 0, их достаточно много в ряде из-за особенности динамики системы. Чтобы получить эти результаты, пришлось использовать довольно значительное усреднение: т = 151. В целом, можно сказать, что реконструкция нелинейной функции успешна: удаётся восстановить её периодический характер и она похожа на гармоническую.
Результаты реконструкции в возбудимом режиме показаны на рис. 6 для воздействия гауссовыми импульсами — a и прямоугольными импульсами — Ь. Вся наблюдаемая динамика
-0.2
Gauss Pulses
Square Pulses
Harmonic
10 20 Ф
30
b
10 20 Ф
30
10
20
Ф
Рис. 5. Реконструированная нелинейная функция f (ф) для колебательного режима: a — при импульсном гауссовом воздействии, b — при импульсном прямоугольном, c — при гармоническом. Для общности все результаты приведены при К = 15, гп = 151
Fig. 5. Reconstructed nonlinear function f (ф) for the oscillatory regime under different driving types: a — Gauss pulses, b — square pulses, c — harmonic driving. For generality, all results are given for К = 15, rn = 151
Gauss Pulses
0.2
0.1 -
0.0 -
-0.1 -
-0.2
•
1
• * 1
* & • 1 1
K=15
3
Ф
0.2
0.1 -
0.0 -
-0.1 -
-0.2
Square Pulses
i
"ЛЧ - j
K=15
12
b
3 4 Ф
Рис. 6. Реконструированная нелинейная функция f (ф) для возбудимого режима: a — при импульсном гауссовом воздействии, b — при импульсном прямоугольном. Для общности все результаты приведены при К =15, гп = 151
Fig. 6. Reconstructed nonlinear function f (ф) for the excitable regime under different driving types: a — Gauss pulses, b — square pulses. For generality, all results are given for К = 15, rn = 151
0
0
0
a
c
0
1
2
4
5
6
0
5
6
7
a
была сосредоточена на одном периоде нелинейной функции. При этом наблюдаются существенные искажения и функция в целом довольно сильно зашумлена, что в первую очередь объясняется длинными интервалами у ~ 0 в этом режиме. Это значит, что для значительного числа ф во временном ряде нет или очень мало существенно отличных от нуля значений у, что приводит к ошибкам при реконструкции.
3.5. Подбор периода воздействия. До сих пор все результаты были приведены для случая, когда период воздействия Т известен. В действительности, конечно, это не так. В лучшем случае, период может быть известен приблизительно. Однако, как показывает опыт реконструкции неавтономных систем [22,29], даже небольшая ошибка в определении Т может быть весьма критична, так как она быстро приводит к накоплению фазового сдвига во времени на основной гармонике, а на высших — кратно быстрее. Самый простой способ решить эту проблему в эксперименте — перебор пробного Т (обозначим его как Т). Перебор можно осуществлять с любым разумным значением шага ДТ, поскольку никакая перевыборка значений с этим шагом не производится, а значение □ = 2п/Т просто используется как множитель в формуле (7). Тем не менее, для простоты в приведённом ниже примере мы ограничились шагом, равным шагу выборки ДТ = Д . При выборе диапазона перебора можно ориентироваться на расположение локальных экстремумов в сглаженном временном ряде (в зашумлённом ряде они будут стоять почти подряд) — маловероятно, что период воздействия будет много меньше наименьшего из них или много больше наибольшего. В представленном случае перебор производился на отрезке Т е [2; 320].
Для самого сложного в аппроксимации воздействия прямоугольными импульсами результаты подбора — зависимость Ь(Т) приведена на рис. 7. Видно, что даже при использовании относительно небольшой степени тригонометрического полинома К = 5, не позволяющей, как показано на рис. 3, Ь и рис. 4, Ь, удовлетворительно описать форму импульса, зависимость Ь(Т) демонстрирует чёткий минимум на значении Т = Т = 100 (см. Таблицу), соответствующем истинному значению, по крайней мере, с точностью до Д£. Минимум достаточно глубокий, несмотря на наличие шума, — для его визуализации пришлось прибегнуть к логарифмическому масштабу.
a
100
200
300
т
7.0 5.0
3.0 2.0
1.0
0.2
K
15 5
100
200
300
b
т
Рис. 7. Зависимость целевой функции (9) от пробного периода воздействия при реконструкции по временным рядам с воздействием прямоугольными импульсами с использованием тригонометрического полинома степени К = 5, К = 15 и К = 25 при сглаживании гп = 151. Рис. a соответствует колебательному режиму, рис. b — возбудимому (цвет online)
Fig. 7. The dependence of the target function (9) on the trial driving period at reconstruction from time series generated under square pulses using a trigonometric polynomial with order К = 5, К = 15 and К = 25 with smoothing rn = 151. Fig. a corresponds to the oscillatory regime, fig. b corresponds to the excitable regime (color online)
0
0
Таблица. Локальные минимумы целевой функции (9) Table. Local minima of the target function (9)
Колебательный режим Возбудимый режим
К Т = 100 Т = 200 Т = 300 Т = 100 Т = 200 Т = 300
25 0.204 0.208 0.210 0.393 0.409 0.423
15 0.210 0.218 0.243 0.400 0.434 0.545
5 0.251 1.195 2.257 0.552 3.148 5.158
При этом увеличение К приводит к углублению не только и не столько главного минимума, сколько его побочных кратных вариантов — Т = 2Т иТ = 3Т. Эти побочные глубокие минимумы очевидно проистекают из возможности описать частоту воздействия с помощью второй и третьей гармоник, соответственно (при этом период сокращается как раз в 2 или 3 раза); такая возможность улучшается при увеличении степени полинома. Тем не менее, минимум Т = 100 всегда является глобальным. Для гауссовых импульсов зависимости оказались очень сходными, только побочные минимумы были более выражены для относительно малых К. Для гармонического воздействия результаты существенно не менялись при изменении К ^ 3, при использовании К = 1 и К = 2 отсутствовали локальные экстремумы на двойном и тройном, или только тройном, значении периода.
Полученные результаты свидетельствуют, что при использовании аппроксимации воздействия тригонометрическим полиномом даже относительно невысокой степени даже для самых сложных импульсов с вертикальными фронтами легко можно детектировать истинное значение периода воздействия по зависимости целевой функции от пробного периода Ь(Т) как в колебательном, так и в возбудимом режиме. Наличие 10% измерительного шума не может помешать процедуре. Таким образом, незнание периода воздействия никак не ограничивает применимость описанного в работе метода.
Заключение
В работах [19,20] уже были предложены подходы к реконструкции системы ФАПЧ с полосовым фильтром, предложенной в [14]. Они опирались, как и представленный в данной работе алгоритм, на неявную аппроксимацию одной из нелинейных функций модели; такой подход был предложен изначально для систем с запаздыванием [18], а затем развит и для других типов уравнений, например, для обобщённых осцилляторов ван дер Поля [31]. Основная проблема предложенных в [18,19,31] решений заключалась в том, что они требовали численного расчёта второй производной, что сильно ограничивает применимость метода при наличии измерительных шумов: приемлемый уровень шума, как правило, составляет не более 1-2% от среднеквадратичного отклонения сигнала. Как решить эту проблемы в общем случае, пока не ясно, но для ряда систем результат может быть достигнут переходом к интегрированным по времени уравнениям, как это было показано в [32] для осцилляторов ван дер Поля и Рэлея. Рассматриваемая в данной работе модель системы фазовой автоподстройки частоты также допускает подобную модификацию, чем мы и воспользовались.
При этом была решена и вторая, специфичная именно для реконструкции рассматриваемой системы ФАПЧ проблема: чувствительность метода к близким к нулю значениям наблюдаемой. В работах [19,20] векторы состояния, соответствовавшие у ~ 0 вовсе исключались из рассмотрения, поскольку в таком случае сама целевая функция оказывалась плохо определена, что приводило к необходимости исключать до половины всего ряда в колебательных режимах,
а в возбудимом режиме при у = 0 реконструкция вовсе не проводилась, так как в таком режиме значения у ~ 0 соответствуют большей части временного ряда, как это показано в [15,30]. Представленный здесь алгоритм лишён этого недостатка при реконструкции параметров ФАПЧ и внешнего воздействия. Решение естественным образом следовало из перехода к интегрированной системе — соответствующая наблюдаемой переменная у ушла из знаменателя при определении нелинейной функции интегрированной системы и переопределении целевой функции алгоритма. При этом следует понимать, что такая модификация алгоритма возможна не для всякой системы и её применимость зависит от вида уравнений.
Аппроксимация периодического внешнего воздействия была произведена проверенным ранее для осцилляторов ван дер Поля-Тоды [22] способом путём разложения в тригонометрический полином. При этом частота воздействия подбиралась просто перебором, а не с помощью алгоритма нелинейной оптимизации, как в [22], где использовался по сути метод Ньютона, поскольку мощности современных компьютеров допускают такой подход, при этом отсутствует вероятность разбегания алгоритма. Такая аппроксимация показала высокую эффективность и устойчивость: с одной стороны, можно описать достаточно сложное воздействие с близкими к вертикальным фронтами коротких импульсов большой скважности, с другой стороны, даже при существенном увеличении числа гармоник сверх необходимого, метод остаётся устойчивым.
Представляется, что предложенный подход может быть непосредственно использован для реконструкции по экспериментальным данным неавтономного генератора, описанного в [17]. При этом для случая нерегулярного импульсного воздействия представленный подход не годится. Можно предположить, что частично проблему может решить использование радиальных базисных функций для аппроксимации импульсов, где каждому импульсу будет сопоставлена отдельная функция.
Список литературы
1. Безручко Б. П., Смирнов Д. А. Математическое моделирование и хаотические временные ряды. Саратов: ГосУНЦ «Колледж», 2005. 320 с.
2. Hodgkin A. L., Huxley A.F. A quantitative description of membrane current and its application to conduction and excitation in nerve // J. Physiol. 1952. Vol. 117, no. 4. P. 500-544.
DOI: 10.1113/jphysiol.1952.sp004764.
3. FitzHugh R. Impulses and physiological states in theoretical models of nerve membranes // Biophysical Journal. 1961. Vol. 1, no. 6. P. 445-466. DOI: 10.1016/S0006-3495(61)86902-6.
4. Nagumo J., Arimoto S., Yoshizawa S. An active pulse transmission line simulating nerve axon // Proc. IRE. 1962. Vol. 50, no. 10. P. 2061-2070. DOI: 10.1109/JRPR0C.1962.288235.
5. Hindmarsh J.L., Rose R.M. A model of neuronal bursting using three coupled first order differential equations // Proc. R. Soc. Lond. B. 1984. Vol. 221, no. 1222. P. 87-102.
DOI: 10.1098/rspb.1984.0024.
6. Morris C., Lecar H. Voltage oscillations in the barnacle giant muscle fiber // Biophysical Journal. 1981. Vol. 35, no. 1. P. 193-213. DOI: 10.1016/S0006-3495(81)84782-0.
7. Gouesbet G., Letellier C. Global vector-field reconstruction by using a multivariate polynomial L2 approximation on nets // Phys. Rev. E. 1994. Vol. 49, no. 6. P. 4955-4972.
DOI: 10.1103/PhysRevE.49.4955.
8. Packard N. H., Crutchfield J. P., Farmer J.D., Shaw R. S. Geometry from a time series // Phys. Rev. Lett. 1980. Vol. 45, no. 9. P. 712-716. DOI: 10.1103/PhysRevLett.45.712.
9. Baake E., Baake M., Bock H. G., Briggs K. M. Fitting ordinary differential equations to chaotic data // Phys. Rev. A. 1992. Vol. 45, no. 8. P. 5524-5529. DOI: /10.1103/PhysRevA.45.5524.
10.
11.
12.
13.
14.
15.
16.
17.
18.
19.
20
21
22
23
24
25
26
406
Пономаренко В. И., Прохоров М.Д. Оценка порядка и реконструкция модельного уравнения системы с запаздыванием // Письма в ЖТФ. 2006. Т. 32, № 17. С. 73-80. Ishizaka K., Flanagan J. L. Synthesis of voiced sounds from a two-mass model of the vocal cords // Bell. Syst. Tech. J. 1972. Vol. 51, no. 6. P. 1233-1268. DOI: 10.1002/j.1538-7305.1972.tb02651.x. Barfred M., Mosekilde E., Holstein-Rathlou N.-H. Bifurcation analysis of nephron pressure and flow regulation // Chaos. 1996. Vol. 6, no. 3. P. 280-287. DOI: 10.1063/1.166175. Huys Q. J.M., Ahrens M. B., Paninski L. Efficient estimation of detailed single-neuron models // Journal of Neurophysiology. 2006. Vol. 96, no. 2. P. 872-890. DOI: 10.1152/jn.00079.2006. Шалфеев В. Д. Исследование динамики системы фазовой автоподстройки частоты с разделительным конденсатором в цепи управления // Известия высших учебных заведений. Радиофизика. 1968. Т. 11, № 3. С. 397-406.
Мищенко М. А., Шалфеев В. Д., Матросов В. В. Нейроноподобная динамика в системе фазовой синхронизации // Известия вузов. ПНД. 2012. Т. 20, № 4. С. 122-130. DOI: 10.18500/08696632-2012-20-4-122-130.
Matrosov V.V., Mishchenko M.A., Shalfeev V.D. Neuron-like dynamics of a phase-locked loop // The European Physical Journal Special Topics. 2013. Vol. 222, no. 10. P. 2399-2405. DOI: 10.1140/epjst/e2013-02024-9.
Мищенко М. А., Большаков Д. И., Матросов В. В. Аппаратная реализация нейроподобного генератора с импульсной и пачечной динамикой на основе системы фазовой синхронизации // Письма в ЖТФ. 2017. Т. 43, № 13. С. 10-18. DOI: 10.21883/PJTF.2017.13.44806.16737. Sysoev I. V., Prokhorov M.D., Ponomarenko V.I., Bezruchko B.P. Reconstruction of ensembles of coupled time-delay systems from time series // Phys. Rev. E. 2014. Vol. 89, no. 6. P. 062911. DOI: 10.1103/PhysRevE.89.062911.
Сысоева М. В., Сысоев И. В., Пономаренко В. И., Прохоров М.Д. Реконструкция уравнений нейроподобного осциллятора, моделируемого системой фазовой автоподстройки частоты с запаздыванием, по скалярному временному ряду // Известия вузов. ПНД. 2020. Т. 28, № 4. С. 397-413. DOI: 10.18500/0869-6632-2020-28-4-397-413.
Sysoeva M. V, Sysoev I. V., Prokhorov M.D., Ponomarenko V.I., Bezruchko B.P. Reconstruction of coupling structure in network of neuron-like oscillators based on a phase-locked loop // Chaos, Solitons & Fractals. 2021. Vol. 142. P. 110513. DOI: 10.1016/j.chaos.2020.110513. Luttjohann A., Pape H.-C. Regional specificity of cortico-thalamic coupling strength and directionality during waxing and waning of spike and wave discharges // Scientific Reports. 2019. Vol. 9, no. 1. P. 2100. DOI: 10.1038/s41598-018-37985-7.
Смирнов Д. А., Сысоев И. В., Селезнев Е. П., Безручко Б. П. Реконструкция моделей неавтономных систем с дискретным спектром воздействия // Письма в ЖТФ. 2003. Т. 29, № 19. С. 69-76.
Мищенко М.А. Нейроноподобная модель на основе системы фазовой автоподстройки частоты // Вестник Нижегородского университета им. Н. И. Лобачевского. 2011. Т. 5, № 3. С. 279-282.
Мищенко М. А., Жукова Н. С., Матросов В. В. Возбуждение фазоуправляемого генератора импульсным воздействием // Известия вузов. ПНД. 2018. Т. 26, № 5. С. 6-19. DOI: 10.18500/08696632-2018-26-5-6-19.
Мищенко М.А., Ковалева Н.С., Половинкин А. В., Матросов В. В. Возбуждение фазоуправляемого генератора импульсной последовательностью // Известия вузов. ПНД. 2021. Т. 29, № 2. С. 240-253. DOI: 10.18500/0869-6632-2021-29-2-240-253.
Zhang J., Davidson R. M., Wei M., Loew L. M. Membrane electric properties by combined patch clamp and fluorescence ratio imaging in single neurons // Biophysical Journal. 1998. Vol. 74, no. 1. P. 48-53. DOI: 10.1016/S0006-3495(98)77765-3.
27. Jin L., Han Z., Platisa J., Wooltorton J. R. A., Cohen L. B., Pieribone V.A. Single action potentials and subthreshold electrical events imaged in neurons with a fluorescent protein voltage probe // Neuron. 2012. Vol. 75, no. 5. P. 779-785. DOI: 10.1016/j.neuron.2012.06.040.
28. Savitzky A., Golay M.J.E. Smoothing and differentiation of data by simplified least squares procedures // Anal. Chem. 1964. Vol. 36, no. 8. P. 1627-1639. DOI: 10.1021/ac60214a047.
29. Bezruchko B. P., Smirnov D. A. Constructing nonautonomous differential equations from experimental time series // Phys. Rev. E. 2001. Vol. 63, no. 1. P. 016207.
DOI: 10.1103/PhysRevE.63.016207.
30. Сысоев И. В., Сысоева М. В., Пономаренко В. И., Прохоров М.Д. Нейроподобная динамика в системе фазовой автоподстройки частоты с запаздывающей обратной связью // Письма в ЖТФ. 2020. Т. 46, № 14. C. 36-38. DOI: 10.21883/PJTF.2020.14.49665.18267.
31. Sysoev I. V. Reconstruction of ensembles of generalized Van der Pol oscillators from vector time series // Physica D. 2018. Vol. 384-385. P. 1-11. DOI: 10.1016/j.physd.2018.07.004.
32. Sysoev I. V., Bezruchko B. P. Noise robust approach to reconstruction of van der Pol-like oscillators and its application to Granger causality // Chaos. 2021. Vol. 31, no. 8. P. 083118. DOI: 10.1063/5.0056901.
References
1. Bezruchko BP, Smirnov DA. Extracting Knowledge From Time Series: An Introduction to Nonlinear Empirical Modeling. Berlin: Springer; 2010. 410 p. DOI: 10.1007/978-3-642-12601-7.
2. Hodgkin AL, Huxley AF. A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 1952;117(4):500-544.
DOI: 10.1113/jphysiol.1952.sp004764.
3. FitzHugh R. Impulses and physiological states in theoretical models of nerve membranes. Biophysical Journal. 1961;1(6):445-466. DOI: 10.1016/S0006-3495(61)86902-6.
4. Nagumo J, Arimoto S, Yoshizawa S. An active pulse transmission line simulating nerve axon. Proc. IRE. 1962;50(10):2061-2070. DOI: 10.1109/JRPROC.1962.288235.
5. Hindmarsh JL, Rose RM. A model of neuronal bursting using three coupled first order differential equations. Proc. R. Soc. Lond. B. 1984;221(1222):87-102. DOI: 10.1098/rspb.1984.0024.
6. Morris C, Lecar H. Voltage oscillations in the barnacle giant muscle fiber. Biophysical Journal. 1981;35(1):193-213. DOI: 10.1016/S0006-3495(81)84782-0.
7. Gouesbet G, Letellier C. Global vector-field reconstruction by using a multivariate polynomial L2 approximation on nets. Phys. Rev. E. 1994;49(6):4955-4972. DOI: 10.1103/PhysRevE.49.4955.
8. Packard NH, Crutchfield JP, Farmer JD, Shaw RS. Geometry from a time series. Phys. Rev. Lett. 1980;45(9):712-716. DOI: 10.1103/PhysRevLett.45.712.
9. Baake E, Baake M, Bock HG, Briggs KM. Fitting ordinary differential equations to chaotic data. Phys. Rev. A. 1992;45(8):5524-5529. DOI: /10.1103/PhysRevA.45.5524.
10. Ponomarenko VI, Prokhorov MD. Evaluating the order and reconstructing the model equations of a time-delay system. Tech. Phys. Lett. 2006;32(9):768-771.
DOI: 10.1134/S1063785006090100.
11. Ishizaka K, Flanagan JL. Synthesis of voiced sounds from a two-mass model of the vocal cords. Bell. Syst. Tech. J. 1972;51(6):1233-1268. DOI: 10.1002/j.1538-7305.1972.tb02651.x.
12. Barfred M, Mosekilde E, Holstein-Rathlou NH. Bifurcation analysis of nephron pressure and flow regulation. Chaos. 1996;6(3):280-287. DOI: 10.1063/1.166175.
13. Huys QJM, Ahrens MB, Paninski L. Efficient estimation of detailed single-neuron models. Journal of Neurophysiology. 2006;96(2):872-890. DOI: 10.1152/jn.00079.2006.
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
408
Shalfeev VD. Investigation of the dynamics of a system of automatic phase control of frequency with a coupling capacitor in the control loop. Radiophysics and Quantum Electronics. 1968;11(3): 221-226. DOI: 10.1007/BF01033800.
Mishchenko MA, Shalfeev VD, Matrosov VV. Neuron-like dynamics in phase-locked loop. Izvestiya VUZ. Applied Nonlinear Dynamics. 2012;20(4):122-130 (in Russian). DOI: 10.18500/0869-6632-2012-20-4-122-130.
Matrosov VV, Mishchenko MA, Shalfeev VD. Neuron-like dynamics of a phase-locked loop. The European Physical Journal Special Topics. 2013;222(10):2399-2405. DOI: 10.1140/epjst/e2013-02024-9.
Mishchenko MA, Bolshakov DI, Matrosov VV. Instrumental implementation of a neuronlike generator with spiking and bursting dynamics based on a phase-locked loop. Tech. Phys. Lett. 2017;43(7):596-599. DOI: 10.1134/S1063785017070100.
Sysoev IV, Prokhorov MD, Ponomarenko VI, Bezruchko BP. Reconstruction of ensembles of coupled time-delay systems from time series. Phys. Rev. E. 2014;89(6):062911. DOI: 10.1103/PhysRevE.89.062911.
Sysoeva MV, Sysoev IV, Ponomarenko VI, Prokhorov MD. Reconstructing the neuron-like oscillator equations modeled by a phase-locked system with delay from scalar time series. Izvestiya VUZ. Applied Nonlinear Dynamics. 2020;28(4):397-413 (in Russian). DOI: 10.18500/0869-66322020-28-4-397-413.
Sysoeva MV, Sysoev IV, Prokhorov MD, Ponomarenko VI, Bezruchko BP. Reconstruction of coupling structure in network of neuron-like oscillators based on a phase-locked loop. Chaos, Solitons & Fractals. 2021;142:110513. DOI: 10.1016/j.chaos.2020.110513. Luttjohann A, Pape HC. Regional specificity of cortico-thalamic coupling strength and directionality during waxing and waning of spike and wave discharges. Scientific Reports. 2019;9(1):2100. DOI: 10.1038/s41598-018-37985-7.
Smirnov DA, Sysoev IV, Seleznev EP, Bezruchko BP. Reconstructing nonautonomous system models with discrete spectrum of external action. Tech. Phys. Lett. 2003;29(10):824-827. DOI: 10.1134/1.1623857.
Mishchenko MA. Neuron-like model on the basis of the phase-locked loop. Vestnik of Lobachevsky University of Nizhni Novgorod. 2011;5(3):279-282 (in Russian).
Mishchenko MA, Zhukova NS, Matrosov VV. Excitability of neuron-like generator under pulse stimulation. Izvestiya VUZ. Applied Nonlinear Dynamics. 2018;26(5):6-19 (in Russian). DOI: 10.18500/0869-6632-2018-26-5-6-19.
Mishchenko MA, Kovaleva NS, Polovinkin AV, Matrosov VV. Excitation of phase-controlled oscillator by pulse sequence. Izvestiya VUZ. Applied Nonlinear Dynamics. 2021;29(2):240-253 (in Russian). DOI: 10.18500/0869-6632-2021-29-2-240-253.
Zhang J, Davidson RM, Wei M, Loew LM. Membrane electric properties by combined patch clamp and fluorescence ratio imaging in single neurons. Biophysical Journal. 1998;74(1):48-53. DOI: 10.1016/S0006-3495(98)77765-3.
Jin L, Han Z, Platisa J, Wooltorton JRA, Cohen LB, Pieribone VA. Single action potentials and subthreshold electrical events imaged in neurons with a fluorescent protein voltage probe. Neuron. 2012;75(5):779-785. DOI: 10.1016/j.neuron.2012.06.040.
Savitzky A, Golay MJE. Smoothing and differentiation of data by simplified least squares procedures. Anal. Chem. 1964;36(8):1627-1639. DOI: 10.1021/ac60214a047. Bezruchko BP, Smirnov DA. Constructing nonautonomous differential equations from experimental time series. Phys. Rev. E. 2001;63(1):016207. DOI: 10.1103/PhysRevE.63.016207.
30. Sysoev IV, Sysoeva MV, Ponomarenko VI, Prokhorov MD. Neural-like dynamics in a phase-locked loop system with delayed feedback. Tech. Phys. Lett. 2020;46(7):710-712.
DOI: 10.1134/S1063785020070287.
31. Sysoev IV. Reconstruction of ensembles of generalized Van der Pol oscillators from vector time series. PhysicaD. 2018;384-385:1-11. DOI: 10.1016/j.physd.2018.07.004.
32. Sysoev IV, Bezruchko BP. Noise robust approach to reconstruction of van der Pol-like oscillators and its application to Granger causality. Chaos. 2021;31(8):083118.
DOI: 10.1063/5.0056901.
Сысоева Марина Вячеславовна — родилась в Саратове (1987). Окончила с отличием факультет нано- и биомедицинских технологий Саратовского государственного университета по направлению «Биомедицинская инженерия» (магистр, 2011). Защитила диссертацию на соискание учёной степени кандидата физико-математических наук по специальностям «Биофизика» и «Радиофизика» (2015, СГУ). Доцент кафедры «Радиоэлектроника и телекоммуникации» и старший научный сотрудник лаборатории радиотехнических устройств и видеотехники СГТУ имени Гагарина Ю. А., старший научный сотрудник СФ ИРЭ им. В. А. Котельникова РАН. Научные интересы — нейронаука, математическое моделирование, радиотехническое моделирование. Автор 34 статей, индексируемых в базах данных Scopus и Web of Science.
Россия, 410019 Саратов, ул. Зелёная, 38
Саратовский филиал Института радиотехники и электроники
им. В. А. Котельникова РАН
Россия, 410054 Саратов, Политехническая, 77
Саратовский государственный технический университет имени Гагарина Ю. А. E-mail: bobrichek@mail.ru ORCID: 0000-0003-4228-9124 AuthorID (eLibrary.Ru): 651344
Корнилов Максим Вячеславович — родился в Саратове (1988). Окончил Лицей математики и информатики (2005), механико-математический факультет (2009, бакалавр математики) и факультет нано- и биомедицинских технологий (2011, магистр техники и технологии) Саратовского государственного университета. Защитил диссертацию на соискание учёной степени кандидата физико-математических наук по специальности «Радиофизика» (2015, СГУ). Доцент кафедры системного анализа и автоматического управления СГУ имени Н.Г. Чернышевского, старший научный сотрудник СФ ИРЭ им. В. А. Котельникова РАН. Научные интересы — анализ временных рядов, математическое моделирование биологических процессов, автоматическое распознавание изображений, компьютерное зрение. Автор 6 статей, индексируемых в базах данных Scopus и Web of Science.
Россия, 410019 Саратов, ул. Зелёная, 38
Саратовский филиал Института радиотехники и электроники
им. В. А. Котельникова РАН
Россия, 410012 Саратов, Астраханская, 83
Саратовский национальный исследовательский государственный
университет имени Н. Г. Чернышевского
Россия, 603950 Нижний Новгород, пр. Гагарина, 23
Нижегородский государственный университет им. Н. И. Лобачевского
E-mail: kornilovmv@gmail.com
ORCID: 0000-0003-1634-9007
AuthorID (eLibrary.Ru): 681980
Такаишвили Лев Вячеславович — родился в Саратове (1998). Окончил МОУ СОШ № 77 (2016). С отличием окончил Институт физики Саратовского государственного университета имени Н.Г. Чернышевского (магистр, 2022). Научные интересы — радиофизические модели биологических нейронов.
Россия, 410012 Саратов, Астраханская, 83
Саратовский национальный исследовательский государственный университет имени Н. Г. Чернышевского E-mail: nar7187@yandex.ru
Матросов Валерий Владимирович — родился в 1960 году. Окончил Горьковский (Нижегородский) государственный университет им. Н.И. Лобачевского по специальности «Прикладная математика» (1982). Доктор физико-математических наук (2007), профессор (2009). Профессор кафедры теории колебаний и автоматического регулирования (с 2007), заведующий этой кафедрой (с 2013), декан радиофизического факультета ННГУ им. Н.И. Лобачевского (с 2014). Под его руководством защищено 4 кандидатских диссертации. Член диссертационного совета 24.2.340.03 при ННГУ (радиофак); научно-методического совета исследовательской школы «Колебательноволновые процессы в природных и искусственных средах». Почетный работник сферы образования РФ, награжден знаком НТОРЭС им. А. С. Попова «За заслуги в развитии радиоэлектроники и связи». Автор 52 статей, индексируемых в базах данных Scopus и Web of Science.
Россия, 603950 Нижний Новгород, пр. Гагарина, 23
Нижегородский государственный университет им. Н. И. Лобачевского
E-mail: matrosov@rf.unn.ru
ORCID: 0000-0003-3146-111X
AuthorID (eLibrary.Ru): 24397
Сысоев Илья Вячеславович — родился в Саратове (1983). Окончил факультет нелинейных процессов Саратовского государственного университета имени Н. Г. Чернышевского по специальности «Физика» (СГУ, 2004). Защитил диссертацию на соискание учёной степени кандидата физико-математических наук (2007) и доктора физико-математических наук (2019) по специальности «Радиофизика». Профессор кафедры системного анализа и автоматического управления СГУ, старший научный сотрудник СФ ИРЭ им. В. А. Котельникова РАН, ведущий научный сотрудник кафедры теории колебаний и автоматического регулирования ННГУ им. Н. И. Лобачевского. Зам. главного редактора журнала «Известия вузов. ПНД». Научные интересы — исследование сигналов биологической природы методами нелинейной динамики, исследование эффективности и модернизация подходов к анализу сигналов. Автор 69 статей, индексируемых в базах данных Scopus и Web of Science.
Россия, 410019 Саратов, ул. Зелёная, 38
Саратовский филиал Института радиотехники и электроники
им. В. А. Котельникова РАН
Россия, 410012 Саратов, Астраханская, 83
Саратовский национальный исследовательский государственный
университет имени Н. Г. Чернышевского
Россия, 603950 Нижний Новгород, пр. Гагарина, 23
Нижегородский государственный университет им. Н. И. Лобачевского
E-mail: ivssci@gmail.com
ORCID: 0000-0002-9203-5504
AuthorID (eLibrary.Ru): 128457
