РАДИОФИЗИКА, ЭЛЕКТРОНИКА, АКУСТИКА
Спектральная плотность мощности условного марковского
импульсного процесса
Ю. В. Ушаков0, A.A. Дубков6
Нижегородский государственный университет имени Н. И. Лобачевского, радиофизический факультет, кафедра математики. Россия, 603950, г. Нижний Новгород, просп. Гагарина, д. 23.
E-mail: "nomah@list.ru, ь dubkov@rf.unn.ru
Статья поступила 16.02.2010, подписана в печать 28.04.2010
На примере модели нейронной системы, генерирующей условную марковскую последовательность дельта-импульсов, показана процедура вывода выражения для оценки спектральной плотности мощности такого сигнала.
Ключевые слова: немарковский процесс, скрытый марковский процесс, условный марковский процесс, нейрон, спайк, спектральная плотность мощности.
УДК: 537.86, 621.371.399, 519.21. PACS: 05.40.-a, 02.50.Еу, 02.50.Ga.
Введение
Многие физические и биологические процессы могут быть описаны как точечные, т. е. в виде последовательности случайных точек на оси времени. Типичными примерами являются импульсные сигналы нервных клеток (нейронов), функционирующих в существенно зашумленной среде нервной системы [1], а также скачкообразное поведение интенсивности излучения лазера в некоторых практически интересных режимах работы [2]. Актуальность спектрального анализа таких процессов не вызывает сомнений. Однако аналитический подход в этой области вызывает серьезные затруднения в том случае, когда точечный процесс не является процессом восстановления. В настоящей работе рассматривается немарковская последовательность импульсов действия (спайков) модели нейрона, которая, согласно результатам статьи [3], может быть описана скрытой марковской цепью, т. е. является условным марковским процессом [4]. Предложено аналитическое выражение для оценки спектральной плотности мощности (СИМ) такого рода процессов.
В первом разделе статьи дано описание исследуемой модели и необходимых для дальнейшего изложения результатов работы [3]. Во втором разделе описывается хорошо известная процедура вывода общего выражения для СИМ точечных процессов, дальнейшее развитие которой в третьем разделе приводит к выводу формулы СИМ процесса со скрытой марковской цепью и результатам ее использования в изучаемой системе.
1. Модель
Рассмотрим кратко результаты работы [3]. Изучаемая система является моделью нейронного ансамбля слухового анализатора млекопитающих. Она состоит из трех возбудимых элементов, два из которых моделируют периферийные сенсорные нейроны (будем называть их сенсорами), находящиеся под воздействием гармонических сигналов Л^ собО^?, а третий элемент моделирует интернейрон (ИН), принимающий спайки
от сенсоров и генерирующий аналогичные спайки, которые передаются другим нейронам. Существенную роль в функционировании нервных клеток играет шумовое воздействие огромного количества «соседей» (~104 на нейрон), поэтому в рассматриваемой модели присутствуют аддитивные источники шума. Основным изучаемым в работе объектом является спайковая последовательность ИН. В силу наличия шума в системе и одинаковой формы спайков этот сигнал рассматривается как последовательность случайных межспай-ковых интервалов (МСИ). Она является немарковской, поскольку на вход ИН поступают спайки сенсоров с неэкспоненциальной плотностью распределения МСИ (РМСИ). Система схематично изображена на рис. 1. В качестве ее базового элемента использовалась безразмерная математическая модель нейрона «пороговый интегратор с утечкой», которая записывается в виде стохастического уравнения Ланжевена: и = -¡ли + /ext(0 + £(0< гДе у(0 мембранный потенциал нейрона, ¡л — релаксационный параметр, /ext(0 — внешний ток, Ç(f) — белый гауссов шум. Необходимо добавить, что для этой модели нейрона существует
m
m
Рис. 1. Исследуемая модель из трех нейронов
граничное условие: при достижении мембранным потенциалом заданного порога yth подразумевается генерация спайка, и v(t) сбрасывается в некоторое фиксированное значение. Спайки моделируются дельта-функциями Дирака.
В статье [3] показано, что в том случае, когда частоты гармонических сигналов, действующих на сенсоры, находятся в отношении несократимой дроби: =т/п, ИН можно считать системой с M состояниями: M = т + п — I. В момент генерации спайка ИН переключается из одного состояния в другое, и время для него начинает отсчитываться заново («память» о прежнем движении «стирается»). Зная плотности РМ-СИ сенсоров, удается получить распределение времени первого достижения (РВПД) порога генерации мембранным потенциалом ИН в каждом возможном состоянии. Любое из этих РВПД устроено таким образом, что каждый его пик, а точнее временной промежуток, в центре которого расположен максимум пика, соответствует переходу в другое состояние и однозначно его определяет. В результате можно изобразить РВПД состояний ИН, как показано на рис. 2, где с помощью обозначения (/ -¥ k) отмечены временные промежутки, при попадании в которые межспайкового интервала ИН переключается из состояния i в состояние k. Таким образом, для вероятности такого переключения ж ¡и
записываем
Я"/* :
p(i){t)dt
(j—ik)
и замечаем, что она зависит только от РВПД /э(,)(0 в текущем состоянии, т. е. данный, вообще говоря, немарковский процесс может быть описан с помощью скрытой марковской цепи с матрицей переходов {7Гд,}. Наблюдаемая последовательность МСИ называется в таком случае условным марковским процессом (см. [4, 5]).
Вполне понятно, что вероятностное РМСИ интернейрона /9ои1(т) (т будем использовать для значения случайного МСИ) получается усреднением РВПД по состояниям, для чего необходимо найти стационарные вероятности состояний р,-, решив систему уравнений [6]
(м-1
= 6=1,2.....М- 1,
;=о
М-1
;=о
М-1
В итоге находим рт[(т)= ^ Р1Р{1>(т), которое удобно
1=0
сравнить с РМСИ ИН, полученным численно, и найти вполне приемлемое сходство (рис. 2,д), подтвержда-
ем
JL
о
1
10
-А
20
—г 30
40 t
Рис. 2. РВПД состояний {а, б, в, г) и результирующее РМСИ (д) интернейрона при П] =0.4, Пг = (3/2)0]. Тонкой пунктирной линией показано соответствие пиков РМСИ и РВПД. Цифры в окружностях показывают, в какое состояние переходит ИН, если в текущем состоянии генерирует МСИ, попадающий в область помеченного пика. Цифры со стрелкой в скобках (/ -»• к) также обозначают некоторые промежутки, в которые должен попасть МСИ для переключения ИН из состояния I в состояние к. На графике д тонкая сплошная кривая изображает теоретический результат, а жирная пунктирная — результат численного моделирования
20 ВМУ. Физика. Астрономия. „Y> 5
ющее определенную степень справедливости предыдущих рассуждений.
2. Общее выражение для СПМ
Рассмотрим несложную процедуру вывода выражения для СПМ импульсного процесса на конкретном примере исследуемого сигнала, являющегося последовательностью дельта-функций:
Г— 1
где 1Г — момент генерации г-го спайка ИН.
Для этого вводится комплексная амплитуда сигнала г
СЦш) = / сИ, и тогда СПМ может быть найдена
о
как предел [7]
= Нш
Т—>оо
№) 1 т
(1)
ЯМ = т-т
1 + Нш — Ие
Л/—>оо N
(2)
Рассмотрим одну реализацию достаточно большой длительности Т, содержащую большое число спай-ков N. Для изучаемого процесса справедливо полагать длину среднего МСИ (т) = Т/Ы [8]. Тогда выражение (1) без труда сводится к следующему:
■ N г— 1
ЕЕ
,г=2 ?=1
В случае независимых МСИ (процесс восстановления), что имеет место, например, для МСИ сенсоров, среднее от комплексной экспоненты в выражении (2) распадается на произведение средних для отдельных интервалов, и усреднение выполняется с помощью известного РМСИ р(т) сенсора. Двойное суммирование последовательно вычисляется по формуле для суммы геометрической прогрессии, и после несложных преобразований [7] получается формула
2 1 - |0М|
(3)
э/ц){1Г-1Я)\ —
0\ш{тг +ТГ+1
£>/"■4'г ' ' г+1 1 ^
' X р(тг, 7>+1, . . . , Тц_1) йтг... йтч_ 1, (4)
где г, = — и — случайный МСИ, а р(тг, тг+\, ■ ■ ■ ... — совместная плотность вероятности выпа-
дения подряд межспайковых интервалов т>, тг+\ и т.д.
Рассмотрим для наглядности частный случай
е^р{п) йтх х
е>ип'2р{т2/т1)с1т2
е^р(т3/тит2)с1т3. (5)
В силу стационарности процесса усреднение не должно зависеть от того, в каком месте на оси времени берется первый МСИ в степени комплексной экспоненты, т. е. среднее зависит только от количества МСИ
усреднения: (е-
/з(ш), или в общем случае
(6)
Таким образом, будем считать, что генерация случайного значения т\ происходит в некотором неизвестном состоянии /. А так как все состояния и сами имеют определенные вероятности выпадения р; в скрытой цепи Маркова, необходимо провести соответствующее усреднение. Поэтому
М-1
р(т1)=^р1р{Чт1) = рш(т1).
(7)
(=0
Если случайное значение т\ выпадает, когда ИН находится в /-м состоянии, то длина интервала т\ определяет состояние к, в которое переключается ИН, согласно р(1){п). После переключения в к-е состояние случайная величина второго МСИ т2 определяется плотностью р(к){т2). Таким образом, для первой условной плотности вероятности в (5) имеем
р(т2/п) = <
р<0>(т2), при пб(г^-О), р<г>(т2), при п е(/ —> 1),
>р<м-1)(т2), при пб(г -¥М — V)
(8)
где 0(ш) = / р(т)е1шт йт — характеристическая функция последовательности МСИ сенсора, а 5т_^т)(ш) — СПМ флуктуаций сигнала. Знак интеграла без пределов здесь и далее означает интегрирование от 0 до оо.
3. Условный марковский процесс
В первом разделе была показана возможность рассмотрения процесса генерации МСИ интернейрона как условного марковского процесса. В каждом из скрытых состояний известно время первого достижения порога генерации. При этом генерация того или иного интервала т означает переключение в соответствующее состояние. В общем случае усреднение в выражении (2) проводится следующим образом:
(для примера см. рис. 2).
Вообще говоря, для некоторого р(1){т) промежутки значений т, переключающие ИН в некоторое состояние циклически повторяются. Однако в интересующих частных случаях число ненулевых пиков всякого р(1){т) меньше или равно числу состояний М. Поэтому предполагается однозначное соответствие пиков и состояний. В результате выражение (5) с учетом (6)-(8) переписывается в следующем виде:
М-1
1,к=0
/зМ=Х> е,ипхр{1)(т\) йт\ х
((—>А)
е^р(К\т2)йт2
е^р(т3/тит2)с1т3. (9)
Условная плотность вероятности р{т3/т\,т2) рассматривается аналогичным образом, из чего следует независимость переменных интегрирования в выражении (9) и становится возможным записать следующее:
М-1
/з е/ьУ°( О^х
1X1=о
((—>А)
е>ы'рт(0 сИ
е^'р(!){1)сИ. (10)
(А'—>/)
Очевидно, что усреднение комплексной экспоненты с большим числом МСИ в показателе выполняется просто добавлением необходимого числа операций суммирования и соответствующих интегралов.
По аналогии с определением характеристической функции введем обозначения
в;{ш) =
М-1
Ыи) =
р(1){0е>ы> сИ.
¡и!
(/—>А")
Очевидно, 9,{ш) = ^ В итоге выражение (10)
А'—0
для произвольного числа МСИ п коротко записывается в векторно-матричной форме
Ш=рв^1{ш)в{ш),
(11)
где р — вектор-строка с элементами po.pi.....рм-ь
в(ш) — матрица М хМ с элементами 0;а.(ш), которая возводится в степень п—1, а в(ш) — вектор-столбец с элементами 0о(^),#1(ш).....6м-1(ш)- И выражение (2) с учетом (6) приводится к следующему:
5(ш) = — П + 2 Игл Ие
(т) 1 А/—>оо
'N-1
п=1
N
>. (12)
При ш ф 0 по крайней мере одна из норм матрицы
0{ш), а именно ||0(а>)|| = шах^ \9ц{ш)\, очевидно, не
<' / '
превышает единицы. Этого условия достаточно [9] для сходимости матричного ряда в выражении (12) при N -¥ оо. Применяя формулу суммирования убывающей геометрической прогрессии, получаем конечное выражение СПМ спайковой последовательности ИН
ЭД = -1-{1 + 2/>Ее[£~
(г)
емг'ем}, (13)
где Е — единичная матрица и [Е — @{ш)]~
матри-
ца, обратная к (Е — 0{ш)).
На рис. 3-5 сплошными кривыми изображены СПМ спайковой последовательности ИН для различных отношений частот П^Пг- На тех же графиках пунктирными кривыми показаны результаты, полученные по формуле (3), т.е. в предположении, что РМСИ ИН рои1(т) имеет место для независимых одинаково распределенных МСИ и в(ш) = § р01й(1)е1и1 с1.1. Хорошо видно, что полученная таким образом кривая отражает усредненную зависимость 5(ш) и в некоторых случаях (рис. 3 и 4) «упускает» многие детали СПМ в прак-
Рис. 3. СПМ спайковой последовательности ИН. П] =0.4, Пг = (3/2)0]
Рис. 4. СПМ спайковой последовательности ИН. П] =0.4, Пг = (5/4)0]
21 ВМУ. Физика. Астрономия. „М' 5
S( со)
0.30 0.25 0.20 0.15 0.10 0.05 0
Рис. 5. СИМ спайковой последовательности ИН. fij =0.4, О2 = (16/15)0]
тически важной низкочастотной области. Однако для случая, изображенного на рис. 5, сплошная и точечная кривые очень схожи: в этом построении частоты fii и fi2 мало отличаются, и пики СПМ «расставлены» на частотах, кратных средней частоте входных гармонических сигналов (fii +Q,-2)/2.
Заключение
В работе представлена процедура вывода выражения для оценки СПМ условного марковского процесса для вполне конкретной математической модели. Однако такая модель может быть использована в достаточно широком круге приложений, поскольку без труда и принципиальных сложностей дополняется большим числом входных и промежуточных элементов. Подобные модели находят применение при изучении сенсорных нейронных систем [10], а скрытая марковская модель вообще активно используется для распознавания речи, графических символов и других задач цифровой обработки сигналов [11]. Таким образом, результаты проведенного в настоящей работе анализа могут быть весьма полезны, по крайней мере в указанных областях исследований.
Power spectral density of the conditional Markov pulse process Yu.V. Ushakov11, A. A. Dubkov''
Department of Mathematics, Faculty of Radiophysics, N. I. Lobachevsky Nizhnii Novgorod State University, 23 Gagarin ave., Nizhnii Novgorod 603950, Russia. E-mail: "nomah@list.ru, b dubkov@rf.unn.ru.
The procedure of power spectral density expression obtaining is provided for the example of the neural system model generating conditional Markov series of delta-pulses.
Keywords: non-Markov process, hidden Markov process, conditional Markov process, neuron, spike, power spectral density.
PACS: 05.40.-a, 02.50.Ey, 02.50.Ga. Received 16 February 2010.
English version: Moscow University Physics Bulletin 5(2010).
Сведения об авторах
1. Ушаков Юрий Владимирович — ассистент (ННГУ); тел.: (831) 462-32-81, e-mail: nomah@list.ru.
2. Дубков Александр Александрович — канд. физ.-мат. наук, доцент, зав. кафедрой (ННГУ); тел.: (831) 462-32-81, e-mail: dubkov@rI.unn.ru.
Работа выполнена при финансовой поддержке РФФИ (грант 08-02-01259-а).
Список литературы
1. Gerstner IF., Kistler \V.M. Spiking Neuron Models. Single Neurons, Populations, Plasticity. Cambridge, 2002.
2. Arecchi F.T., Meucci R. 11 Eur. Phys. J. B. 2009. 69. P. 93.
3. Ушаков 10.В. 11 Вестн. ННГУ. Радиофизика. 2010. В печати.
4. Стратонович Р.Л. Условные марковские процессы и их применение к теории оптимального управления. М., 1965.
5. Rabiner L.R. // Proc. IEEE. 1989. 77, N 2. P. 257.
6. Тихонов В.И., Миронов М.А. Марковские процессы. М., 1977.
7. Стратонович Р.Л. Избранные вопросы теории флюктуа-ций в радиотехнике. М., 1961.
8. Тихонов В.И. Статистическая радиотехника. М., 1966.
9. Приклонский В.И. Численные методы. М., 1999.
10. Lopera A., Buldu J.M., Torrent М.С. et al. 11 Phys. Rev. E. 2006. 73. P. 021101.
11. Королёв A.B., Силаев A.M. 11 Изв. вузов. Радиофизика. 2005. XLVIII, № 4. С. 358.