CC BY
30
Радиофизика
Вестник Нижегородского университета им. Н.И. Лобачевского, 2010, № 4 (1), с. 67-72
УДК 537.86+621.371.399+519.21
СТАТИСТИКА ВЫХОДНОЙ СПАЙКОВОЙ ПОСЛЕДОВАТЕЛЬНОСТИ НЕЙРОННОЙ МОДЕЛИ СЛУХОВОГО АНАЛИЗАТОРА
© 2010 г. Ю.В. Ушаков
Нижегородский госуниверситет им. Н.И. Лобачевского
nomah@list.ru
Поступила в редакцию 22.12.2009
Изучается система из трёх нейроподобных элементов под воздействием двух гармонических сигналов и шума. Представлен вероятностный анализ механизма генерации импульсов (спайков), позволяющий аналитически описать распределение межспайковых интервалов на выходе системы. Теоретические кривые подтверждаются численным экспериментом.
Ключевые слова: нейрон, сигнал, шум, статистика, радиофизика.
Введение
Удивительные функциональные возможности нервной системы человека, животных и даже простых организмов уже многие годы привлекают внимание специалистов различных ес-тественно-научных областей. Нервные клетки (нейроны) очень детально исследованы физиологами и описаны разнообразными математическими моделями, которые в настоящее время активно изучаются в том числе и методами радиофизики [1].
Представленная в данной работе система уравнений моделирует входную часть слухового анализатора млекопитающих. Согласно исследованиям, базилярная мембрана среднего уха выполняет Фурье-анализ звукового сигнала [2, 3] и далее каждая частотная компонента действует на отдельный сенсорный нейрон (сенсор), преобразующий синусоиду в последовательность коротких электрических импульсов (спайков), которые отправляются на вход промежуточных нейронов (интернейронов - ИН). Для предупреждения недоразумений скажем здесь, что в дальнейшем изложении понятие нейрон применяется, вообще говоря, к математической модели нейрона, т.е. нейроподобному элементу.
Изучаемый в данной статье ансамбль составлен из двух сенсоров и одного ИН. Поскольку в реальных нервных системах шумы играют весьма существенную роль [1], введено действие аддитивного белого гауссова шума на каждый элемент модели. Целью работы было выяснение статистических закономерностей преобразования сигнала с двумя частотами в спектре зашумлённой нелинейной системой.
Современные исследования подобных моделей (см., например, [4-6]) обычно направлены на выявление различных аспектов конструктивной роли шума [7] и, в основном, сосредоточены на поведении таких частных характеристик, как, например, отношение сигнал-шум для одной из спектральных компонент входного сигнала, в зависимости от шумовой интенсивности. В представленной работе получены выражения для более общей характеристики: распределения межспайковых интервалов (РМСИ) выходного сигнала ИН. Таким образом, результаты могут использоваться во многих исследованиях нейроподобных и других пороговых систем с шумом, которыми традиционно интересуется радиофизика.
Далее в статье приводится детальное описание исследуемой модели и её параметров, затем следует теоретический анализ модели, и в заключении обсуждаются полученные результаты.
Модель
Для исследования заявленных вопросов использовалась математическая модель, описанная в общих чертах выше. В качестве её базового элемента выбран нейрон типа «пороговый интегратор с утечкой», который при всей простоте достаточно полно передаёт необходимые для данного исследования принципы функционирования живой нервной клетки [8]. Все описанные далее динамические переменные и параметры следует считать безразмерными, хотя они имеют прямую связь с наблюдаемыми в экспериментах физическими величинами.
Итак, математически изучаемую систему можно записать следующим образом:
Vi = -ЦМ + A COS £ty + -Щ1,\ (О,
• v2 = -H2v2 + Л COS Q2r + ^2 (0, (1)
v = -|j,v + + k2s2{t) +
где V;(/) моделирует мембранный потенциал /-го сенсора; ц, - параметр релаксации; At и Q, -амплитуды и частоты гармонических входных сигналов, соответственно; Dt — интенсивность шума сенсора; <f,(/) - 5-коррелированный белый гауссов шум /-го сенсора с нулевым средним. В третьем уравнении системы (1) (уравнение для ИН) используются аналогичные символы без индексов. ИН принимает на вход спайковые последовательности сенсоров S\ 2(t), взвешенные С ПОМОЩЬЮ коэффициентов k\ И &2- Описание этих последовательностей следует ниже. Все источники шума в системе - £1,2(0 и Щ) — считаются статистически независимыми.
Необходимо пояснить, что пороговый интегратор с утечкой, вообще говоря, не имеет никакого механизма генерации спайков. При расчетах вводится порог «зажигания» vth, по достижении которого мембранный потенциал
нейрона сбрасывается в значение v®2 для сенсоров или v° для ИН, и мы говорим, что нейрон сгенерировал спайк. Таким образом, последовательность спайков, приходящих с выхода сенсора на вход ИН, определяется моментами времени генерации этих спайков, и каждый из них представляется в виде 8-функции Дирака. В данном исследовании порог зажигания всегда равен единице для всех трёх нейронов. Поскольку генерация спайка требует достаточной для достижения порога величины внешнего воздействия, система является нелинейной по определению. В модели (1) это скрыто за сокращенной записью спайковых последовательностей сенсоров
N,(0 , ч
Sj(t) = где tij — это моменты време-
)= О
ни зажигания сенсоров (/ = 1, 2). Именно они, а также их количество iV)(f), отсчитываемое при t > 0, нелинейным образом связаны с vi>2(0- В силу наличия шумов в уравнениях сенсоров значения Ц являются ещё и случайными. Поэтому получить аналитическую зависимость tjj(yi(t),vth, А„ П„ Df) не представляется возможным. В указанной форме записи уравнений автоматически определяется возбудимость связей между элементами, то есть ИН испытывает влияние не всей временной динамики мембранных потенциалов сенсоров, а только сгенерированных ими спайков.
Представленные ниже численные и аналитические результаты получены при следующих постоянных значениях некоторых параметров: Ц1 = ц2 = 1, И- = 0.3665, = /)2 = И = 1.6-10 3,
к\ = к2 = 0.97, V® =у2 =0, V0 = -1. Частоты входных сигналов находились в соотношении несократимой рациональной дроби т/п.
Начальные условия при компьютерном моделировании и теоретическом анализе системы (1) полагаются следующими: у12 (0) = у°2 , у(0) — у0.
Потенциал сброса ИН равен минус единице в целях моделирования такого реально существующего и важного эффекта, как период рефрактерности: сразу после генерации спайка нейрон в течение определенного промежутка времени тте/ не способен повторить импульс и вообще не реагирует на внешние воздействия. Учет рефрактерности позволяет исключить из рассмотрения слишком частое генерирование спайков ИН. В численном исследовании реф-рактерность отслеживается по уровню у(7), что определяет значение гге/ « 6.28, которое используется в теории. Поскольку временная дистанция между импульсами сенсоров и без того задается периодами внешних синусоид, введение для них интервала рефрактерности не является необходимым.
Подытоживая описание модели, необходимо сказать об ограничениях применяемых параметров.
1) Входные сигналы являются подпороговы-
ми, то есть в отсутствие шума они не способны зажигать сенсоры. Из решения уравнений сенсоров в системе (1) [9] для подпорогового синусоидального сигнала имеем: Д / +1 < %.
2) Параметры релаксации сенсоров и частоты входных сигналов связаны таким образом, что на каждом максимуме входной синусоиды может быть сгенерировано не более одного спайка. Однако в то же время «спайк на каждом периоде» - это наиболее вероятная ситуация, по сравнению, например, с генерацией спайка через период.
3) Как и входные сигналы сенсоров, отдельные входные спайки для ИН считаются подпо-роговыми, а два спайка одновременно оказывают сверхпороговое воздействие. Формально это записывается следующим образом:
1\г <1=
+к2 >1.
Итак, мембранный потенциал сенсора совершает подпороговые гармонические колеба-
ния вслед за входным сигналом с тем же периодом. Согласно решению уравнения для сенсора, к этому движению аддитивно добавляется «окрашенный» гауссов шум [9] той же интенсивности Д. Вблизи максимума синусоиды шумовая флуктуация может обеспечить достижение порога, то есть генерацию спайка. Как уже было сказано, для выбранных параметров спайк на каждом периоде - это наиболее вероятный вариант. С меньшей вероятностью спайк будет сгенерирован через два периода от предыдущего, с ещё меньшей - через три периода, и так далее. Типичная реализация у,{^) и плотность РМСИ сенсора представлены на рис. 1.
V. V к/ \А/~
0 10 20 30 40 1
Рі 1 л |Л ..1 Л1 ...
10
20
30
40 І
Рис. 1. Типичная реализация мембранного потенциала сенсора и характерная гистограмма РМСИ, совпадающая с распределением времени первого достижения порога генерации. А{= 1.165, П, - 0.6
ИН испытывает воздействие только 5-им-пульсов (стайков сенсоров) и белого шума, поэтому реализация у(/) - это аналогичный окрашенный шум с вертикальными всплесками (толчками от входных спайков), затухающими экспоненциально во времени, и сбросами при зажигании.
Теоретический анализ
Возвращаясь к статистике времён генерации спайков сенсорами, следует указать, что для них задача о нахождении РМСИ решена различными способами с удовлетворительной степенью точности в работах предшественников [8]. К тому же такие распределения быстро оцениваются с помощью численного моделирования. В результате будем полагать РМСИ ргдОО сенсоров известными. То есть начиная с момента времени 1 = 0 нам известна статистика входных спайков ИН. Так как эта статистика непу-ассоновская, процесс у(/) не является марковским [10]. Поэтому для его анализа выполняется детальное вероятностное рассмотрение.
Вероятностные гипотезы
Найдём вероятность генерации спайка интернейроном в дифференциально малом промежутке времени (ґ, 1+ск). Для указанных выше ограничений на параметры системы имеем четыре вероятностные гипотезы такого события.
1. Интернейрон генерирует импульс по приходу отдельного спайка первого сенсора. Отдельным спайком будем называть импульс, который поступает на вход ИН по истечении периода рефрактерности после зажигания или переходного процесса, следующего за толчком от предыдущего входного спайка. Поскольку спайк в нашей системе моделируется 5-функцией, отпадает необходимость учитывать производную процесса у(?) до момента его прихода. Вероятность рассматриваемой ситуации, таким образом, есть вероятность того, что в данный момент времени сумма полученного толчка и собственного шума достигнет порога генерации. А вероятность того, что всё это произойдет в указанном интервале времени, совпадает с вероятностью генерации спайка первым сенсором: рі (*)<#. В результате получаем полный вклад данной вероятностной гипотезы, который запишем следующим образом: рі(0Фоі(£і)*#, где
Фоі(*і) = ^ег£с|^^(у,й -Лі)| - это и есть вероятность того, что спайк, пришедший от первого сенсора, в сумме с собственным шумом ИН окажется не меньше порога генерации. Выражение получается из стационарной плотности вероятности упомянутого окрашенного шума. Стационарный её вариант можно использовать потому, что за каждым спайком ИН следует период реф-рактерности хгф который в данном случае длиннее характерного масштаба установления постоянного значения дисперсии шума. Функция егіс(х) - это дополнительный интеграл ошибки.
2. Интернейрон генерирует импульс по приходу отдельного спайка второго сенсора. Ситуация полностью аналогична предыдущей.
3. Интернейрон генерирует импульс по приходу спайка первого сенсора на фоне воздействия спайка второго сенсора, который пришёл чуть раньше. Иначе говоря, в какой-то момент Ґ в прошлом спайк от второго сенсора не зажёг ИН, однако потенциал у(І') совершил при этом вертикальный скачок до уровня к2+0]П^г)> после чего экспоненциально затухает по закону е~^. В текущий момент времени { > Ґ приходит спайк первого сенсора. Если за время і-Ґ эффект от предыдущего толчка не затух до уровня шума, то мы не можем считать
данный спайк первого сенсора отдельным, учесть это феноменологически, с помощью та-поскольку система «помнит», что было до это- *
го. Таким образом, данная гипотеза представ- кого множителя: 1- I . Во-вторых, в си-
ляет собой два независимых события: спайк второго сенсора не зажёг ИН в момент /'ив данный момент і пришел спайк первого сенсора. Вероятность гипотезы есть произведение вероятностей этих событий: (1 - ФогС^г)) Рі(0 Ш. Домножая её на условную вероятность генерации спайка ИН, получаем следующее слагаемое:
Рі(0^(і -Ф02(^2)) \ Р2(ОФі(*Ь*2»*'- *')<*'»
гаях.{і-Т2,хг#}
где Тг = (1/ц.)1п(^2/^1/2) - это период затухания импульса второго нейрона до уровня В п (амплитуды шума),
■'чГ
лу наличия периода рефрактерности некоторые области РМСИ сенсоров могут просто выпасть из формулы (2), что приведет здесь к ненормированному распределению. Поэтому использовать указанный множитель можно только после нормировки. В итоге получаем следующую запись распределения вероятности первого спайка ИН во времени:
Р(0 = Р(0
'■гкГ
(3)
Фі(^2^-0 = |егй|^(уА-Аі-^ ** °)| -
вероятность того, что спайк первого сенсора зажжёт ИН, если подействует на него спустя время І Ґ после прихода спайка второго сенсора. Интеграл учитывает все возможные Ґ, попадающие в период Т2 «памяти» ИН. Очевидно, ЧТО Фі(*ь &2, *“0 ~> ФоіС*і) при (*“0 > °с, то есть данная ситуация превращается в первую, когда система «забывает» о том, что в прошлом приходил недостаточный толчок от другого сенсора.
4. Интернейрон генерирует спайк по приходу импульса второго сенсора на фоне импульса от первого сенсора. Рассуждения совпадают с предыдущим случаем.
Убедившись, что все четыре ситуации действительно исключают друг друга, то есть не могут быть реализованы одновременно, суммируем их вклады, получая плотность распределения вероятности генерации первого спайка интернейроном во времени:
РСО=Рі (ОФ01&)+Р2 (0^02 (^2)+
I
+Рі(/Хі-ф02(^)) (2)
і
+Р2(0(і-Фоі&)) /рі(ОФг(*2 А,*
где р(0 = 0 При І < їгф
Необходимо, однако, учесть еще два фактора, влияющих на это распределение. Во-первых, нужно помнить, что ненулевая вероятность генерации спайка ИН в некотором интервале (/, і-<М) снижает вероятность найга его спайк в следующем промежутке времени. Постараемся
где р(0 - это нормированный результат выражения (2).
Множество состояний интернейрона
Для каждого сенсорного нейрона РМСИ совпадает с плотностью распределения времени первого достижения порога мембранным потенциалом. Действительно, на каком бы периоде ни произошло зажигание сенсора, сразу после сброса можно считать момент генерации
Состояние 0
0.5
1
0.5
0.5
л С0 (9)
" VI/ ■ , ^ і : (-]: и ; А: ©<
10
20
30
40 ♦
р,( 0
і
10
20
30
10
;©
20
Состояние 1
30
40 і
Р2(1)
0
I
1
0.5
0
I
1
0.5
0
..Л®...;.
■ і ® Л ; Л. ;л.
10
20
30
1.3
©
_Л
10
20
30
р (0
©
®
А.
10
20
30
40 *
Рис. 2. Плотности распределения вероятностей генерации спайков трёх нейронов в разных состояниях для отношения частот входных сигналов 01/02 = 4/3. Пример возможного перехода между двумя состояниями. См. пояснения в тексте
спайка нулевым и наблюдать идентичное распределение времени достижения порога в следующий раз.
Иначе обстоит дело с ИН. Известные плотности распределения спайков сенсоров во времени (р1(0 и рг(0) определяют плотность распределения во времени входных спайков ИН. Каждый из этих спайков с определённой вероятностью может зажечь ИН, после чего время последнего отсчитывается от момента зажигания. В этом новом времени плотность распределения входных спайков смещена влево, по сравнению с предыдущей ситуацией. Будем говорить, что каждое зажигание ИН переключает его в новое состояние, которое определяется формой плотности распределения входных спайков во времени. Если частоты синусоид, действующих на сенсоры, относятся как целые числа т!п, то каждый /и-й пик Р1<7) совпадает с каждым п-м пиком рг(0 во времени, т.е. вся схема расстановки пиков Ри(0 периодически повторяется, следовательно, интернейрон имеет конечный набор состояний.
Посмотрим на рис. 2. В исходном состоянии (состояние 0) входящие в ИН импульсы распределены во времени в соответствии с пиками р1,2(0- Если, например, ближайший возможный спайк первого сенсора (состояние 0, Р1(0, пик 1) будет сгенерирован (вероятность такого события равна площади под первым ПИКОМ Р](/) в нулевом состоянии) и вызовет зажигание ИН, то после его сброса распределение входных спайков во времени будет смещено (состояние 1), что необходимо учитывать в выражении (2). Теперь ближайший спайк с определенной вероятностью придет от второго сенсора (состояние 1, рг(0> пик
2). Однако он, вероятнее всего, при данных параметрах системы попадёт в период реф-рактерности, равный »6.28. Поэтому следующим импульсом, способным зажечь ИН, может быть спайк опять первого сенсора (состояние 1, ра (0, пик 3), хотя вероятность его генерации тоже не равна, конечно, единице. И так далее. Вполне очевидно, что после каждого зажигания ИН оказывается в разных состояниях по отношению к распределению вероятности входящих спайков во времени.
Согласно наиболее простому предположению, ИН «обходит» весь набор состояний равномерно. В таком случае результирующее РМСИ на выходе системы получается усреднением распределений времени первого достижения порога мембранным потенциалом ИН по всему набору:
Р-(0 = ^=Ь------------, (4)
1# £ р0)(0
о у'=о
где М — число возможных состояний интернейрона.
На каждом периоде биений двух синусоидальных сигналов число несовпадающих пиков Р:(0 и р2(0 есть (т - 1) + (и - 1). Существует также один общий пик. Отсюда заключаем, что ИН имеет М = т + п — 1 различных состояний. Так, например, отношение частот в ситуации на рис. 2 есть = 4/3, поэтому М = 4+3-1 = 6. Цифры в окружностях, от 0 до 5, расставлены около всех видимых пиков распределений для установления однозначного соответствия между каждым пиком и состоянием ИН.
Алгоритм вычисления РМСИ интернейрона
Резюмируя изложенные рассуждения, представим их в виде пошаговой процедуры. Для нахождения плотности РМСИ интернейрона при заданных параметрах системы (1) необходимо выполнить следующее.
1) Получить РМСИ р!(/) и р2(0 сенсоров путём численного моделирования системы (1) с исключенным ИН или с помощью теоретических процедур, описанных в [8].
2) Определить все возможные состояния ИН по отношению к пикам распределений Р1,г(0» то есть найти смещения этих распределений (см. рис. 2) для подстановки в выражение (2).
3) С помощью выражений (2) и (3) посчитать
р ^ (?) для каждого возможного состояния.
4) Наконец, получить роц((?) по формуле (4).
В качестве подтверждения правильности
сделанных предположений и выводов на рис. 3 изображено сравнение результатов работы приведенного алгоритма с результатами, полученными путем прямого численного моделирования системы (1). Налицо хорошее совпадение кривых даже для довольно сложных форм плотности распределения.
Заключение
На примере математической модели, описывающей нейронный ансамбль на входе слухового анализатора млекопитающего, в работе проанализирован механизм генерации импульсов зашумлённым пороговым элементом под действием независимых случайных спайковых последовательностей. Развит аналитический подход к нахождению плотности распределения межепайко-
Л А л
§ 10 15 20 T Кварта. Отношение частот: 4/3.
out - A jt J
5 10 15 20 25 T Малая терция. Отношение частот: 6/5.
м out J
5 10 15 20 Г
Малая секунда. Отношение частот: 9/8,
Рис. 3. Гистограммы РМСИ интернейрона для различных отношений частот с названиями этих отношений в музыкальной терминологии. T - длина межспайкового интервала в выходном сигнале интернейрона. Сплошная линия - теоретический результат. Прерывистая получена с помощью компьютерного симулирования модели (1). Параметры: (кварта) A1 = 1.165, П1 = 0.6, А2 = 1.085, П2 = 0.45; (малая терция) A1 = 1.125, П1 = 0.54, А2 = 1.085, П2 =
= 0.45; (малая секунда) А1 = 1.2, ^ = 0.675, А2= 1.165, fi2 = 0.6
вых интервалов в выходном сигнале. Таким образом, поведение этого распределения может быть детально изучено при варьировании параметров источников шума, входных сигналов и самой системы, что требуется, например, при выявлении разного рода конструктивного шумового влияния. Использованный метод получения аналитических выражений позволяет без труда расширять их для более сложных входных сигналов и более широких областей параметров системы, что обеспечивает возможность приложения результатов данной работы в других нейронных и радиофизических системах. Исследование актуально также в контексте многих нерешенных проблем, которые связаны с восприятием гармоничных и диссонансных музыкальных созвучий млекопитающими [2, 11], поскольку экспериментальные данные о функционировании мозга
обнаруживают фундаментальность механизмов, лежащих в основе эффектов сенсорного восприятия [3, 12].
Благодарности
Выражаю благодарность доц. А.А. Дубкову за консультации при получении и интерпретации результатов работы.
Работа выполнена при финансовой поддержке РФФИ, грант № 08-02-01259-а.
Список литературы
1. Gerstner W., Kistler W.M. Spiking neuron models. Single neurons, populations, plasticity. Cambridge University Press, 2002. 400 p.
2. Helmholtz H.L.F. On the sensations of tone. New York: Dover, 1954. 576 p.
3. Скучик Е. Основы акустики. Т. II / Под ред. Ю.М. Сухаревского. М.: Ин. лит., 1959. 566 с.
4. Balenzuela P. Neural mechanism for binaural pitch perception via ghost stochastic resonance // Chaos. 2005. V. 15. P. 023903-1-023903-8.
5. Lopera A., Buldu J.M., Torrent M.C. et al. Ghost stochastic resonance with distributed inputs in pulse-coupled electronic neurons // Phys. Rev. E. 2006. V. 73. P. 021101-1-021101-6.
6. Анищенко В.С., Ануфриева М.В., Вадивасова Т.Е. Стохастический резонанс в бистабильной системе под воздействием хаотического сигнала // Письма в ЖЭТФ. 2006. Т. 32. № 20. С. 12-17.
7. Анищенко В.С., Вадивасова Т.Е., Астахов В.В. Нелинейная динамика хаотических и стохастических систем. Саратов: Изд-во Саратовского ун-та, 1999. 368 с.
8. Burkitt A.N. A review of the integrate-and-fire neuron model: II. Inhomogeneous synaptic input and network properties // Biol. Cybern. 2006. V. 95. P. 97-112.
9. Gardiner C.W. Handbook of stochastic methods: for physics, chemistry and natural sciences. 2 ed. New York: Springer, 1996. 442 p.
10. Тихонов В.И., Миронов М.А. Марковские процессы. М.: Сов. радио, 1977. 488 с.
11. Plomp R., Levelt W.J.M. Tonal consonance and critical bandwidth // J. Acoust. Soc. of America. 1965. V. 38. P. 548-560.
12. Fishman Y.I., Volkov I.O., Noh M.D., et al. Consonance and dissonance of musical chords: neural correlates in auditory cortex of monkeys and humans // J. Neurophysiol. 2001. V. 86. P. 2761-2788.
THE OUTPUT SPIKE TRAIN STATISTICS OF THE AUDITORY ANALYZER NEURAL MODEL
Yu.V. Ushakov
The system of three neural-like elements stimulated by the mixture of two sinusoidal signals and a white Gaussian noise is considered. The probabilistic analysis of the spike generation mechanism is presented which makes it possible to get analytical expressions for the interspike interval probability distribution at the system output. Theoretical results are proved by numerical simulations.
Keywords: neuron, signal, noise, statistics, radiophysics.
