ДИСКРЕТНАЯ МОДЕЛЬ НЕЙРОННОЙ АКТИВНОСТИ
В.И. Некоркин, Л.В. Вдовин
В работе представлена новая модель, описывающая хаотические спайк-берстовые колебания нейронов, заданная в виде двумерного разрывного отображения. Модель получена на основе дискретной модификации модели нейрона ФитцХью - Нагумо и разрывного отображения типа Лоренца. Исследована динамика модели, найдены значения параметров, при которых в системе возникает хаотический аттрактор, соответствующий спайк-берстовым колебаниям, изучены его свойства и характеристики. Показано, что модель позволяет описывать и другие режимы нейронной активности: подпороговые колебания, одиночные периодические и хаотические спайки.
Введение
Нервная клетка (или нейрон) является [1] базовым элементом мозга. Посредством электрических и химических синапсов разной полярности нейроны формируют разнообразные пространственно распределенные крупномасштабные сети, состоящие из большого числа элементов. Поэтому изучение ключевых функциональных свойств мозга тесно связано с исследованием коллективной активности разнообразных нейронных сетей. Одним из эффективных подходов исследования такой активности является моделирование нейронных сетей с помощью нелинейных динамических систем [2]. Очевидно, что при таком подходе прежде всего необходимо построить подходящие динамические системы, описывающие отдельный нейрон. Предъявляемые к таким моделям требования отчасти являются противоречивыми. С одной стороны, они должны описывать сложные многокомпонентные физико-химические процессы в клетке, вызывающие разнообразные формы активности (например, под-пороговые, спайковые и берстовые колебания), наблюдаемые в «живых» нейронах. Благодаря этому, модели, достаточно подробно описывающие процессы в нервной клетке, представляют собой сложные нелинейные системы. Классическим представителем моделей такого типа является система Ходжкина - Хаксли [3], основанная на детальном анализе ионного транспорта через мембрану клетки. С другой стороны, в реальных нейронных сетях число элементов, охваченных многочисленными синаптическими связями, достигает сотен тысяч. Поэтому при моделировании таких сложных систем приходится упрощать модели индивидуальных нейронов, насколько это возможно в рамках рассматриваемой задачи. Например, в искусственных нейронных сетях используются так называемые формальные нейроны - элементы с простейшей бинарной логикой [4]. Для разрешения отмеченного противоречия при
моделировании нейронных сетей очень часто используют феноменологические модели, описывающие основные особенности динамики нейронов, однако, детально не учитывающие действие ионных каналов мембраны. К таким моделям относятся, например, системы ФитцХью - Нагумо [5], Хиндмарш - Розе [6], Морриса - Ле-кара [7] и др. Эти модели представляют собой системы нелинейных обыкновенных дифференциальных уравнений. В последнее время получил популярность и другой класс феноменологических моделей - системы с дискретным временем в виде нелинейных точечных отображений [8-11]. Высокую эффективность такие модели показали при исследовании коллективных процессов в крупномасштабных нейронных сетях. Например, в [12,13] проведено исследование одномерных и двумерных многоэлементных динамических сетей точечных отображений, моделирующих процессы активации и коллективных колебаний кортикальных нейронных сетей. Установлено хорошее качественное соответствие между динамическими структурами в этих сетях и режимами активности, полученными ранее с помощью различных вариантов модели Ходжкина - Хаксли.
В широком ареале режимов нейронной активности важное место занимают так называемые спайк-берстовые колебания - последовательности из трех и более потенциалов действия (спайков), обычно возникающих на волне деполяризации. Во-первых, такие колебания достаточно широко представлены в природе - они характерны для пирамидальных нейронов гиппокампа [14], кортикальных СН [15] и 1В [16] нейронов, таламокортикальных и других нейронов. Во-вторых, считается [17,18], что спайк-берстовые колебания играют ключевую роль в процессах обработки информации нейронными сетями. Наиболее известной моделью, описывающей такую колебательную активность, является модель Хиндмарш - Розе, которая задается системой дифференциальных уравнений третьего порядка, включающей две нелинейные функции. Спайк-берстовые колебания также моделировались и с помощью нелинейных точечных отображений. В работе [19] для одномерных и двумерных отображений достаточно общего вида проведена классификация возможных бифуркационных механизмов возникновения берстовых колебаний. В [20] рассматривается кусочно-линейное одномерное разрывное точечное отображение. Показано, что с помощью такой модели можно описывать спайк-берстовую хаотическую активность нейронов. Кроме того, для взаимодействующей пары отображений проведено исследование различных режимов синхронизации, включая режим хаотической синхронизации. В [11] двумерное точечное отображение было использовано для исследования динамических механизмов рождения и регуляризации хаотических берстов в ансамбле элементов, находящихся в режиме синхронных берстовых колебаний. Еще одна модель в виде двумерного точечного отображения была введена в [21]. Отображение содержит быструю и медленную переменную и две линии разрыва. Изучены [21,22] динамические механизмы генерации спайков, берстов, и режимы их синхронизации в случае системы двух взаимосвязанных отображений. Дальнейшее развитие эта модель получила в [23,24]. В работе [25] изучена синхронизация и распространение берстов в сети отображений [21], моделирующей нейронную сеть, состоящую из нейронов в режиме берстовой генерации, связанных посредством химических и электрических синапсов. В [24] представлена некоторая модификация модели [21] в виде дополнительной квадратичной нелинейности, введенной в отображение. С помощью модифицированного отображения стало возможным моделировать так называемые подпороговые колебания - колебания ниже порога возбуждения потенциала
действия. Подобные колебания наблюдаются, например, в оливомозжечковой системе, отвечающей за формирование универсального ритма (около 10 Гц) в центральной нервной системе [26,27].
В данной работе предлагается новая, достаточно универсальная, дискретная модель нейронной активности. Модель представляет собой двумерное точечное отображение с одной линией разрыва, реализующее динамические механизмы, характерные для канонического в теории динамических систем отображения Лоренца и для широко известной в нейродинамике модели ФитцХью - Нагумо. Будет показано, что, несмотря на сравнительную простоту, модель позволяет описывать многие режимы нейронной активности: спайк-берстовые колебания, в том числе хаотические, под-пороговые колебания, режим одиночной, периодической и хаотической генерации спайков. Однако главное внимание уделяется моделированию хаотических спайк-берстовых колебаний.
Статья организована следующим образом. В разделе 1 вводится исследуемая модель. В разделе 2 рассматривается релаксационная динамика модели. Исследуются быстрые и медленные движения и доказывается существование хаотического релаксационного аттрактора, соответствующего спайк-берстовым колебаниям модели. Раздел 3 посвящен исследованию модели в случае, когда динамика модели не является релаксационной. Получены условия существования инвариантной области, содержащей хаотический аттрактор, определяющий спайк-берстовые колебания модели. В разделе 4 обсуждаются другие возможные аттракторы модели и отвечающие им режимы нейронной активности. В Заключении дается краткое обсуждение полученных результатов.
1. Модель
Обозначим через f : К2 — К2 отображение (х,у) — (х,у), задаваемое следу ющей системой:
'х = х + Г(х) - у - ря(х - ,1) (1)
у = у + е(х - 3),
где х-переменная качественно описывает динамику мембранного потенциала нервной клетки, у - совокупное действие всех ионных токов, проходящих через мембрану нейрона и отвечающих за восстановление покоя мембраны, а функции ¥ (х) и Н(х — ,) имеют вид
¥ (х) = х(х - а)(1 - х), 0 <а< 1, (2)
х> 0,
Н (х) = { (3)
[0, х < 0.
Параметр е (е > 0) определяет характерный временной масштаб изменения восстанавливающей переменной у, параметр .] контролирует уровень деполяризации мембраны (3 < ,), а параметры в (в > 0) и ^ (, > 0), как будет показано далее, характеризуют порог возбуждения берстовых колебаний.
В случае в = 0 система (1) представляет собой дискретную версию модели ФитцХью - Нагумо. Поэтому некоторые черты динамики этой модели характерны
и для системы (1). С другой стороны, очевидно, что при в > 0 за счет наличия пороговой функции Н(х — й) при подходящем выборе параметров й и в, в динамике системы (1) могут появляться принципиально новые свойства. Потребуем, чтобы параметры й и в удовлетворяли следующим условиям:
3шт < й < 3шах> (4)
в (3шах) — F(3шт), (5)
где
3ш
3ш
1 + а — V1 — а + а
2
1 + а + л/1 — а + а
2
3
При таком выборе й линия разрыва х = й расположена на возрастающем участке функции F(х), что, как будет показано ниже, имеет принципиальное значение для хаотизации отображения /. Условие (5) устанавливает взаимное расположение ветвей изоклины вертикальных наклонов слева и справа от линии разрыва х = й. Поскольку система (1) содержит функцию Н, отображение / является разрывным.
Множество точек разрыва имеет вид
^
А =и /-г8, 8 = {(х,у): х = !}. (6)
г=0
Поскольку для всех е > 0
^ = 1 + Р (х)+ е > 0, х е М+, (7)
0(х,у)
где
М+ = {х : г- < х < г+}
г
±
1 + а ± V 4 — а + а
2
2
отображение / является [28] взаимно-однозначным, если х е К+. Будем рассматривать траектории отображения /, не имеющие общих точек с А и целиком расположенные в области
Р+ = {(х,у): х е М+, х = й, у>Р (3шах) — в} .
Заметим, что ограничение по у, присутствующее в Р+, введено с целью устранения несущественного описания поведения траекторий на фазовой плоскости вне этой области.
Для удобства изложения представим отображение / в виде
{/г, х <й, /2, х > й,
где
/г : (х,у) ^ (х + Р(х) — у, у + е(х — 3)), /2 : (х,у) ^ (х + Р(х) — у — в, у + е(х — 3)).
3
2. Релаксационная динамика системы и спайк-берстовые колебания
Рассмотрим динамику отображения f при е << 1. В этом случае переменная y с течением времени изменяется медленно по сравнению с переменной x. Поэтому в системе (1) присутствует два разных масштаба времени и скорости. Быстрые движения определяются, главным образом, «замороженной» системой, в которой параметр е = 0. Медленные же движения описываются уравнением для переменной y при е = 0, зависящим, однако, от состояния быстрой переменной, которое является квазиравновесным. В результате под действием медленной и быстрой систем в (1) может возникнуть режим, состоящий из чередующихся временных интервалов, в течение которых система находится в квазиравновесном состоянии, и «мгновенных» скачков, переводящих систему из одного квазиравновесного состояния в другое. Такие режимы называются релаксационными колебаниями. Изучим такие колебания в системе (1).
2.1. Быстрые движения. В первом приближении быстрые движения в отображении f (1) описываются следующей системой
x = g(x), y = y0 = const, (8)
где
{gi(x) = x + F(x) — y0, x < d, g2(x) = x + F(x) — y0 — в, x > d.
0.2 0.4 0.6 0.8 1.0 1.2 Р Рис. 1. Разбиение плоскости параметров (Р,уо) на области, соответствующие различной динамике отображения (8) для й = 0.45, а = 0.01. Серым тоном выделена область Уь
Заметим, что в отображении (8) уо играет роль нового параметра. Зафиксируем в (8) параметры а и й, а параметры у0 и в будем рассматривать как контрольные. Найдем условия, при которых для значений уо, удовлетворяющих неравенству
Уо > F(,7тах) - в, (9)
из условия х(0) € М+, (х(0) = й) следует, что х(п) € М+ для п € Z+. Анализируя вид функции последования, устанавливаем, что это будет иметь место, если
lim g\(x) < r+, lim g2(x) > r . (10)
x yd x\d
Неравенства (5), (9) и (10) выделяют на плоскости (ß,yo) некоторую область Y (рис. 1), имеющую следующие границы:
Во = {(в,Уо) : в = Р(7тах) - Р(7т1п),
Р(7т1п) < Уо < й + Р(й) - Г- - Р(7тах) + Р(7т1п)}, В1 = {(в,уо): Уо = й + Р(й) - г- - в,
Р(7тах) - Р(7т1п) < в < Г+ - Г-}, В2 = {(в,Уо): У0 = й + Р(й) - Г+,
Р(7тах) - й - Р(й) + Г+ < в < Г+ - Г-},
Т1 = {(в,Уо): Уо = Р(7тах) - в,
Р(7тах) - Р(7т1п) < Уо < Р(7тах) - й - Р(й) + Г+}.
Рассмотрим динамику отображения (8) для (в, Уо) € К. Прежде всего, заметим, что функция последования отображения д является монотонно возрастающей при всех
х € М+, х = й.
2.1.1. Неподвижные точки. Для значений параметров из области У отображение (1) в области х > й неподвижных точек не имеет, а в области х < й может иметь одну или две неподвижные точки. Прямые
Ао = {(в,Уо) : Уо = Р(й),
Р(7тах) - Р(7т1п) < в < й - Г-}, Т2 = {(в, Уо) : Уо = Р(7т1п),
Р(7тах) - Р(7т1п) < в < й - Г- + Р(й) - Р(7т1п)}
делят У на три подобласти - У1 ,У2,У3 (см. рис. 1), где
У1 = Ур|{Уо >Р(й)},
У2 = Ур|{Р(7т1п) <Уо (й)}, У3 = Ур|{Уо (7т1п)}-Для (в, Уо) € У2 отображение (8) имеет две неподвижные точки х = х1 и х = х2, где о 1 + а 2(1 - а + а2)1/2 / а п
хт =---сов---I ,
1 3 3 43 3/'
о 1 + а 2(1 - а + а2)1/2 /а . п
х° =---сов —|— ,
2 3 3 \3 3/'
27
сов а = —-
2(1 - а + а2)3/2
(1 + а)(2а2 - 5а + 2) Уо--27-
Неподвижная точка х = х^ является устойчивой, а х = х° - неустойчивой. При переходе из У2 в подобласть У1 неподвижная точка х = исчезает, попадая на линию
разрыва х = й при (в,Уо) € О0. Так что в У1 существует единственная устойчивая неподвижная точка х = х^. Переход из У2 в подобласть У3 сопровождается касательной бифуркацией, при которой неподвижные точки х = х0 и х = х0 сливаются при (в, Уо) € Т и исчезают. В подобласти У3 отображение (8) неподвижных точек не имеет.
2.1.2. Динамика отображения (8) для подобласти У\. В случае
(в, Уо) € У\ отображение (8) имеет единственную устойчивую неподвижную точ-
о. 1.
ку х = хо. При этом справедливы неравенства
92 (й) <91 (й) <й. (11)
Из условия (11) и монотонного роста функции последования д(х) следует, что для (в, Уо) € У1 все траектории отображения стремятся к неподвижной точке х = х® (рис. 2, а).
2.1.3. Динамика отображения (8) для подобласти У2. Прежде всего, заметим, что для (в, Уо) € У2 отображение (8) имеет две неподвижные точки, а функция 9(х) удовлетворяет условиям
91 (й) > й, 92(й) < й. (12)
Положим для краткости Ь = д2(й), с = д1(й). Рассмотрим образы точек х = Ь и х = с под действием (8). Поскольку при х > й отображение (8) не имеет непо-
Рис. 2. Динамика отображения (8) при а = 0.125, й = 0.45 для различных подобластей У^: а - (в, у0) € Ух: в = 0.4, уо = 0.125; б - (в, уо) € У +: в = 0.24, уо = -0.2; в - (в, уо) € У2\{У2+}: в = 0.6, уо = 0.015; г - (в, уо) € Уз: в = 0.5, уо = 0.2
движных точек, а функция g2 (x) является монотонно возрастающей, то g2 (x) < x и, следовательно,
g2(c) < c. (13)
При x < d имеем g\(x) > x, если x £ (x2, d] и поэтому
g^b) >b, b>x°2. (14)
Неравенство (14) выделяет в Y2 некоторую подобласть Y+, граница которой задается прямыми T2, Во (см. рис. 1) и линией Г
Г = {(ß, уо) : d + F(d) - yo - ß - x2(yo) = 0}-
Для (ß, yo) £ Г отображение (8) имеет гомоклиническую к x = x°j траекторию [29]. Для точек области Y2+ одновременно выполняются условия (13), (14) и, следовательно, интервал I = {x : b < x < c} является инвариантным для отображения (8). Таким образом, для (ß,y0) £ Y2+ отображение (8) обладает бистабильными свойствами - все траектории с начальными условиями x < x2 стремятся к неподвижной точке x = x2, ас начальными условиями x2 < x < b приходят с течением времени в интервал I(y0) и никогда его не покидают (рис. 2, б).
2.1.4. Динамика отображения (8) для подобласти Y3. При переходе из Y2 в Y3 через прямую T2 (см. рис.1) неподвижные точки x = x1 и x = x°j исчезают через касательную бифуркацию, и, следовательно, для (ß, y0) £ Y3 неравенства (14) автоматически выполняются. Отсюда, поскольку (13) для Y3 также выполняется, вытекает существование инвариантного интервала I(y0) для (ß, yo) £ Y3. Для этой области параметров отображение (8) теряет свойство бистабильности и все его фазовые траектории с начальными условиями вне I с течением времени приходят в этот интервал и остаются в нем (рис. 2, г).
Таким образом, отображение (8) имеет инвариантный интервал I(y0) для значений параметров из области Y2+ |J Y3.
2.1.5. Инвариантное множество интервала I. Рассмотрим структуру инвариантного множества, содержащегося в I(yo). Покажем, что в Y2+ |J Y3 существует подобласть, для точек которой отображение (8) на интервале I(yo) действует подобно известному отображению Лоренца. Согласно [28], отображение (8) на I(yo) будет отображением типа Лоренца, если выполняются следующие условия:
1. g'(x) > 0 для всех x £ I(yo)\{d};
2. множество прообразов точки разрыва x = d, D = (J g-n(d) принадлежит
n>o
I (yo);
3. lim g(x) = b, lim g(x) = c.
x\d x yd
Из отмеченных выше свойств функции g(x), очевидно, следует, что условия (1), (3) выполняются. В [28] показано, что условие (2) будет выполненным, если
g'(x) > q> 1, x £ I(yo). (15)
Нетрудно видеть, что условие (15) для отображения (8) выполняется при в,Уо, удовлетворяющих следующим неравенствам:
Уо <й + Р(й) - в - 3 Уо > й + Р(й) - ,1тах-
(16)
На плоскости (в,Уо) неравенства (16) выделяют некоторую область, границы которой задаются двумя прямыми
¿1 = {(в, Уо) : Уо = й + Р(й) - в - 3т1п,
й + Р(й) - Г+ < Уо Утш)}, ¿2 = {(в, Уо) : Уо = й + Р(й) - 7тах),
Р (^тах) й Р (й) + 3тах < в < -Г + 3тах} ■
Для рассматриваемых значений параметра й (см. (4)) прямая Ь2 на плоскости (в, Уо) расположена между прямыми Т2 и В2, а прямая ¿1 параллельна прямой В1 и расположена левее этой прямой (см. рис. 1). Потребуем, чтобы прямая Ь1 на плоскости (в, Уо) была расположена правее прямой Т1. Это требование будет выполнено для значений параметров, удовлетворяющих неравенству
й + Р(й) - Р(3тах) - 3т1п > 0. (17)
При выполнении (17) в области У существует подобласть У включающая в себя У+ и часть подобласти У3 (на рис. 1 подобласть У ь отмечена серым фоном). Для (в, Уо) € Уь отображение (8) на I(уо ) действует подобно отображению Лоренца. Неравенства (4), (17) выделяют на плоскости (а, й) область А+ (рис. 3, а), соответствующую существованию Уь. При этом параметр в должен удовлетворять (см. рис. 1) условиям
Рис. 3. Области параметров (а) и В+ (б), соответствующие существованию Уь
Р (3тах) Р (3т1п)
<в< 3тах 3п
(18)
которые определяют на плоскости (а, в) некоторую область В + (рис. 3, б).
Обозначим через (уо) инвариантное множество отображения (8) на 1(уо). Структура множества существенно зависит от величины параметра ц [28]. Нетрудно показать, что
' д'(с), Уо < й + Р(й) - в/2 - (1 + а)/3, д'(Ь), уо > й + Р(й) - в/2 - (1 + а)/3
Рис. 4. а - Бифуркационная диаграмма отображения (8); цветом показано распределение изображающих точек внутри инвариантного интервала для й = 0.45, а = 0.1, в = 0.3. б, в - Распределение точек внутри инвариантного интервала для различных значений уо: 0.022 (б); -0.043 (в)
и, в свою очередь, зависит от параметров отображения (8) и, в частности, от параметра Уо. Следовательно, при изменении уо структура (уо) также меняется. На рис. 4, а представлена бифуркационная диаграмма отображения (8) где параметр уо играет роль контрольного. Плотность вероятности распределения траекторий множества показана градациями серого цвета. Интенсивность цвета пропорциональна значению плотности вероятности. Диаграмма показывает, что структура (уо) существенно зависит от уо. При некоторых у о траектории множества (у о) полностью покрывают интервал I(уо) (рис. 4, б), а при других уо - локализованы в нескольких непересекающихся интервалах (рис. 4, в). В первом случае (см. рис. 4, б) отображение (8) на I(уо) является строго транзитивным [28], поскольку для любого интервала
к
1о С I(уо) \ {й} существует к > 0 такое, что У дг1о Э 1^1 (уо). Следовательно [28],
г=о
в этом случае периодические точки плотны в I(уо). Множество является хаотическим аттрактором, поскольку все его траектории являются неустойчивыми (д > 1) и не покидают I(уо).
2.2. Медленные движения. Рассмотрим систему (1) в области
3 > 3т1п, (19)
в которой неподвижная точка 0(х = 3,у = Р(3)) отображения f является неустойчивой. Уравнение медленного движения системы (1) представляет собой уравнение эволюции медленной переменной У при условии, что быстрая переменная х поддерживается в квазиравновесном состоянии, то есть
У = У + е(х - 3), (20)
У = Р (х). (21)
Уравнение (21) задает на фазовой плоскости (х,у) так называемую кривую медленных движений [30], являющуюся множеством неподвижных точек для различных Уо отображения (8), определяющего быстрые движения. Кривая медленных движений системы (1) состоит из двух компонент. Одна из них,
Ж?(0) = {(X, у) : у = Р(х), X < 3т1п} ,
образована семейством устойчивых неподвижных точек х = хо(уо) отображения (8) (рис. 5, а) и, следовательно, является устойчивой относительно быстрых движений. Вторая компонента
^(0) = {(х, у) : у = Р (х), 3т1п <х<й} .
образована семейством неустойчивых неподвижных точек отображения (8) х = х2(уо) и является неустойчивой относительно быстрых движений. Поскольку компонента ^^(О) целиком расположена в полуплоскости х < 3, переменная у на ^®(0) под действием отображения (20) монотонно убывает, что приводит к монотонному возрастанию переменной х при движении по кривой ^^(0).
2.3. Релаксационные хаотические спайк-берстовые колебания. Пусть параметры системы (1) принадлежат областям В+ и В + (см. рис. 3). Кроме того, выберем параметр 3 таким, чтобы для всех уо интервал I(уо) на плоскости (х, уо) был целиком расположен правее прямой х = 3. Нетрудно показать, что это требование выполняется в случае выполнения условия
3 + Р(3) - й - Р(й) + в < 0. (22)
При 0 < е << 1 характерные черты разбиения фазовой плоскости (х,у) на траектории в значительной степени сохраняются по отношению к случаю е = 0, однако в случае е = 0 медленные движения происходят не только на кривых медленных движений, но и в тонких слоях порядка еа, 0 < а < 1 , в окрестностях этих кривых. Быстрые движения имеют место вне слоев медленных движений. Обозначим через МБ(е) слой, образованный траекториями медленных движений в окрестности кривой
^(е) = {(х,у): у = Р (х) + + ■■■, х < 3^,
У
0.03 0.02 0.01 0.0
а -о.1 о.1 о.з о.5 х
Рис. 5. а - Хаотический аттрактор А; б - хаотические берстовые колебания. Параметры: J = 0.1, е = 0.001, в = 0.3, й = 0.45, а = 0.1
вид которой был получен стандартными методами теории возмущений. Заметим, что кривая Ш1(е) является устойчивой по отношению к быстрым движениям и, следовательно, быстрые траектории притягиваются слоем Мв(е).
Рассмотрим поведение траекторий системы (1) с начальными условиями из слоя Мв(е). Пусть Ь - одна из таких траекторий. Сначала эволюция Ь определяется, в основном, медленными движениями, которые, главным образом, описываются уравнениями (20), (21). Поэтому траектория Ь движется в слое Мв(е) до тех пор, пока она не достигнет окрестности точки С(х = ,1т\п, у = ¥(1т;п)). В этой точке многообразие Ш1(е) теряет устойчивость, и далее движение траекторий Ь задается системой быстрых движений. При этом эволюция переменной х в значительной степени определяется отображением (8), а переменная у медленно увеличивается в соответствии со вторым уравнением в системе (8). При выходе из окрестности точки С в отображении (8) параметр уо
^^ ¥ (1тт)
и, следовательно, (в,Уо) е Уь (см. рис. 1). Поэтому динамика переменной х соответствующей траектории Ь подчиняется отображению, представленному либо на рис. 2, в, либо на рис. 2, г. Следовательно, после прохождения окрестности точки С переменная х, соответствующая Ь, «быстро» (скачкообразно) попадает в интервал I(уо), в котором эволюционирует подобно траекториям множества , а переменная у изменяется незначительно. Далее переменная у начинает медленно и монотонно расти, поскольку для всех уо интервал I(у0) расположен в области х > ,1. При этом для каждого у переменная х
остается локализованной в интервале, близком к I(уо). Это сохраняется до тех пор, пока переменная у не превысит значение, соответствующее линии Г. После чего в отображении (8) параметр у0 переходит в область У2 \ У2+ и динамика переменной х соответствует отображению, представленному на рис. 2, в. Согласно этому отображению, переменная х асимптотически приближается к устойчивой неподвижной точке. При этом переменная у незначительно меняется на фоне значения уо. Такая эволюция переменных х и у означает, что траектория Ь возвращается в слой Мв(е) и продолжает свое движение в этом слое. Многократное повторение описанных выше быстрых и медленных движений приводит к формированию аттрактора А (см. рис. 5, а) на фазовой плоскости (х,у). На рис. 5, б представлена временная реализация одной из траекторий этого аттрактора, имеющая форму спайк-берстовых колебаний.
Для понимания структуры аттрактора была вычислена его фрактальная размерность df(А) (рис. 6, а), которая для различных 3 изменяется в пределах от 1.12 до 1.6, что свидетельствует о хаотичности этого аттрактора. При увеличении 3 от значения 3т;п фрактальная размерность df(А) монотонно возрастает и при 3 = 3а, За = 0-35 достигает максимального значения, df (А) = 1.6. Это возрастание связано с тем, что увеличение 3 сопровождается уменьшением скорости изменения переменной у внутри «лоренцева» участка аттрактора. Обозначим для удобства этот уча-
„ , ,, т , „ сток через Сь. Дальнейшее увеличение
Рис. 6. Зависимость от параметра J фрактальной г "
размерности аттрактора А (а) и вероятности попа- 3 приводит к резкому уменьшению фрак-
дания точж аттракт°ра в области С- и с+ (б).Па- тальной размерности df (А) и потере мо-
раметры: е = 0.001, 6 = 0.2, й = 0.45, а = 0.1 '
г г 1 нотонности ее изменения. Это объясня-
ется тем, что при достижении некоторого значения Зс > За происходит нарушение условия (22). При этом участок Сь будет разбит на две подобласти: С- = {(х,у) : (х,у) £ Сь,х < 3}, внутри которой векторное поле будет направлено в сторону убывания переменной у, и С+ = {(х,у) : (х,у) £ Сь,х > 3}, внутри которой переменная у возрастает. Вероятность попадания точек аттрактора в эти подобласти в зависимости от параметра 3 представлена на рис. 6, б. При 3 < 3а вклад области С- в динамику системы будет незначителен по сравнению с вкладом области С+, и поэтому фрактальная размерность продолжает возрастать. Однако при достижении критического значения параметра 3а этот баланс нарушается и переменная у начинает совершать хаотические колебания вблизи некоторого значения. В этом случае большинство точек аттрактора будет локализовано в полосе, принадлежащей Сь. Поскольку распределение точек внутри Сь является неравномерным (см. рис. 4), фрактальная размерность df (А) чувствительна к изменению 3 при 3 > 3а.
а^л; 1.6
1.4
1.2
1.0
р-
X т=т ■
£ ''и
ш
3=й
N
а 0.1 0.2 0.3 0.4 3
р■
0.8 Р{С+)
0.4
0.0 б Р(в-) у/
0.1 0.2 0.3 0.4 3
3. Хаотический аттрактор и спайк-берстовые колебания
Покажем, что отображение f имеет хаотический аттрактор не только при £ << 1, но и в некотором конечном интервале значений параметра е. Для этого построим на фазовой плоскости (х, у) область 5, инвариантную относительно отображения f.
3.1. Инвариантная область. Напомним, что 5 будет f -инвариантной областью в том случае, если из условия
(х,у) е 5, (х,у) е А (23)
следует, что (х,у) е 5. Пусть Г^ - граница области 5. Будем искать Г^ в виде кусочно-линейного замкнутого контура.
Прежде всего потребуем, чтобы векторное поле (1) на Г^ было ориентировано внутрь 5 и область 5 С М+ (см. (7)). Из этих условий находим уравнения прямых, образующих Г^,
Г1 = {(х,у): у= —а1, 1 < х < 1тах})
Г2 = {(х,у): у= —к(х — -1) — а1, —1о < х < 1},
Гз = {(х,у): х= — 1о, к(1о + 1) — а1 < у < ¥(1тах) },
Г4 = {(х,у): у= ¥ (1тах 0, — 1о < х < 1},
Г5 = {(х,у): у= — л/£(х 1) + ¥(1max)) 1 < х < 1тах})
Гб = {(х,у): х= 1тах) —а1 < у <л/£(1тах — 1) + ¥ (1тах) },
где
а2 , 2^/1- а + а2
21 — а + а2 — (1 + а)
к = 2 + Ут — ' 1 = "-^-^ (24)
Заметим, что введенная таким образом область 5 определена при выполнении (4), (19) и
£ < а2/4. (25)
Найдем, какими должны быть параметры, чтобы выполнялось условие (23). Представим 5 в виде объединения трех областей, то есть 5 = 5^ 52 У 53 (рис. 7). Пусть (х, у) е 5\. Потребуем, чтобы (х,у) е 5. Это требование эквивалентно условиям
—1о < х < 1тах) (26)
у < шш{¥(1тах), V£(х — 1) + ¥(1тах)}, (27)
у> тах{—а1, —к(х — 1) — а1}. (28)
Рассмотрим условие (26). Подставляя х,у из (1) в (26), получим
у<х + ¥ (х) + 1о, (29)
у>х + ¥ (х) — 1тах- (30)
Рис. 7. Инвариантная область для значений параметров 3 = 0.327, е = 0.008, в = 0.196, d = 0.5, а = 0.25
Неравенства (29) и (30) определяют на фазовой плоскости (х, у) область, граница которой задается прямыми х = — 1о, х = 1 и двумя монотонными линиями
К = {(х,у) : у = х + ¥ (х) + 1о, —1о < х<1},
К2 = {(х,у) : у = х + ¥ (х) — 1тах, —1о < х<1}.
Линия К проходит через точку Б (см. рис. 7) и, следовательно, неравенство (29) выполняется для всех точек области 51. Потребуем теперь, чтобы линия К2 была расположена ниже прямой Г2. Очевидно, что это условие будет выполнено, если К2 либо пересечет прямую х = 1 ниже точки В, либо пройдет через эту точку. Эти требования приводят к следующим условиям на параметры:
1 + ¥(1) — 1тах < —а1.
(31)
Следовательно, в области параметров, выделенной условиями (31), неравенство (26) выполняется.
Рассмотрим неравенства (27), (28), из которых при учете (1) получаем
у > тах
у< —£(х — 1)+ ¥ (1тах),
—(£ + к)(х — 1) — к¥ (х) — а1
—а1 + £(х — 1)
(32)
(33)
1 — к
На фазовой плоскости этим условиям удовлетворяют точки, заключенные между линией
Кз = {(х, у): у = £(х — 1) + ¥ (1тах) — 1о < х < 1}
и ближайшей к ней одной из линий
К4 = {(х,у) : у = е(х — 3) — а3 — 30 < х < 3}, К5 = {(х,у):(1 — к) у = —(е + к)(х — 3) — кЕ (х) — аЗ, —Зо < х < 3}.
Линия Кз - отрезок прямой, проходящей через точку Е (см. рис. 7), расположенный над прямой Г4. Следовательно, для любой точки области неравенство (32) выполняется. Аналогично, К4 также представляет собой отрезок прямой, проходящей через точку В, расположенный ниже Г2. Рассмотрим теперь расположение К5 относительно 5*1. Потребуем, чтобы К5 на плоскости (х,у) располагался ниже 5*1, то есть ниже прямой Г2. Нетрудно убедиться, что это условие будет выполнено, если
ф(х) < 0, —Зо < х < 3, (34)
где
к2 + е
ф(х) = —-—(х — 3) + Е (х) + а3. к
При х € (—3о,3) функция ф(х) имеет единственный минимум в точке х = 0, ордината которого ф(0) = 0. Следовательно, условие (34) выполняется и, значит, неравенства (27), (28) также справедливы.
Таким образом, для значений параметров, удовлетворяющих (31), выполняется условие
¡1 (51) С 5, (х,у) € А. (35)
Совершенно по аналогичной схеме устанавливаются условия на параметры, при которых образы областей 52 и 53 принадлежат 5. Опуская соответствующие выкладки, представим здесь лишь сами условия, а именно
¡1(52) С 5,
если
если
й + Е(й) + а3 < 3тах,
/е<Е (3тах )/(2(й — 3)), (36)
^ /е < (Е(3тах) — Е(й))/(2(й — 3)),
¡25) С Б,
в < Е(3тах) + а3, в < й + Е(й) — (1 — а)3,
< (Е(3тах) — Е(й) + в)/(2(й — 3)),
/е < в/(2(3тах — 3)).
и
Наконец, рассматривая совместно условия (4), (19), (25), (36) и (37), устанавливаем область параметров, для которой область S является f -инвариантной,
Jmin <J< min i Ё—J), d_±m-l, Jm„ - d - F (fl 1 , (38)
a 1 — a a J
/£< min< a F(Jmax) F(Jmax) — F(d) + ß ß 1
2' 2(d — J) ' 2(d — J)) ' 2(Jmax — J)j '
aJmin + F (Jmax ) < ß <d + F(d) — (1 — a) Jmin, d + F (dd) < Jmax aJmin-
(39)
3.2. Хаотический аттрактор. Для каждого фиксированного а е (0,1) неравенства (39) на плоскости (!, в) выделяют область А+ (рис. 8, а), существующую при d > d0(a) (1т;п < d0(a)). Зафиксируем параметры d и в, выбрав их из
области А+. Для этих ! и в неравенства (38) определяют на плоскости (1, л/ё) область Е+, для точек которой 5 является инвариантной областью. Рассмотрим поведение траекторий отображения f с начальными условиями из области 5. Изоклина горизонтальных наклонов - прямая х =л 1 - делит 5 на две части: 51 и 52 У 53 (см. рис. 7). В силу второго уравнения в системе (1) в области 51 выполняется неравенство у < у, а в 52 и 53 - неравенство у > у. Поэтому переменная у вдоль любой траектории с начальными условиями из области 51 монотонно убывает и траектория покидает Л1. После попадания траектории в 52 У 53, наоборот, переменная у начинает монотонно возрастать до тех пор, пока траектория не покинет эту область. Затем процесс повторяется и в области 5 возникает вращательное движение вокруг неподвижной точки О. В результате на фазовой плоскости (х, у) возникает аттрактор А (рис. 9, а). На рис. 9, б представлена зависимость фрактальной размерности от параметра £. При £ < 0.036 фрактальная размерность является дробной величиной (максимальное значение равно 1.4). Заметим, что при увеличении £ от нуля значение df (А) резко изменяется. Это свидетельствует о том, что динамика переменной у существенно влияет на структуру аттрактора А. Начиная с £ = 0.036, аттрактор А перестает быть хаотическим, df (А) = 1. Это объясняется тем, что с увеличением £ увеличивается скорость прохождения «лоренцева» участка аттрактора, и «времени жизни» на этой части аттрактора недостаточно для хаотизации колебаний.
Рис. 8. Области параметров, соответствующие существованию Г: а - A+ (а = 0.25); б - E+ (а = 0.25, ß = 0.19547, d = 0.5012)
Рис. 9. а - Хаотический аттрактор А для значений параметров а = 0.25, J = 0.327, е = 0.00793, в = 0.19547, d = 0.5012, принадлежащих областям А+ и Е+; б - Зависимость фрактальной размерности df (А) от параметра е для значений параметров а = 0.2, J = 0.14, в = 0.265, d = 0.45
4. Галерея других аттракторов и режимов нейронной активности
Покажем, что с помощью отображения f кроме спайк-берстовых колебаний можно моделировать многие другие режимы нейронной активности. Для этого откажемся от ограничений, задаваемых неравенствами (5), (18) и условием У > Е(7тах) - в (см. раздел 2).
4.1. Генерация одиночных спайков и берстов. Как известно [1], одним из основных свойств нейронов, изначально находящихся в состоянии покоя, является их способность к генерации потенциала действия при превышении некоторого порога в результате действия внешнего стимула. Состоянию покоя нейрона в системе (1) отвечает неподвижная точка О, которая является устойчивой в области параметров, выделяемой неравенством
т 1 + а^/Г—а+Ж+Ъё
3 <-з-= Н(а, е). (40)
Рассмотрим отклик нейрона на внешнюю стимуляцию достаточно короткими одиночными импульсами. В системе (1) представим действие такого стимула в виде аддитивного неавтономного слагаемого е(п) в первом уравнении (рис. 10, а). Пред-положим,что е << 1 и, следовательно, динамика системы (1) является релаксационной. Кроме того, считаем, что длительность внешних импульсов настолько мала, что их действие сводится к мгновенному изменению переменной х на величину, равную амплитуде импульса. Обозначим через х = хе, у = уе (уе & Е(3)) значения переменных х, у после окончания действия стимула, а через Е - траекторию системы (1), соответствующую этим начальным условиям (хе,уе).
(1) Если амплитуда стимула является недостаточной (рис. 10, а (1)) для преодоления первого порога, определяемого слоем медленных траекторий в окрестности
неустойчивой инвариантной кривой
^1и(£) = {(х,у) : у = ¥ (х) + + •••' 1тт <х<^,
то наибольшее значение переменной х(п) вдоль траектории Е незначительно превосходит амплитуду стимула (рис. 10, в (1)). В этом случае генерации потенциала действия не происходит (рис. 10, б (1)).
(И) Увеличим амплитуду (см. рис. 10, а (и)) внешнего стимула так, чтобы был преодолен первый порог, но не достигнут второй порог, определяемый слоем медленных траекторий в окрестности инвариантной кривой
^и(£) = {(х,у): у = ¥(х) + £(х—1) — в + ! < х < 1тах} .
В этом случае быстрые движения системы (1) близки к траекториям отображения д на интервале I(уо), где уо е У ь. Поэтому траектория Е совершает некоторое количество нерегулярных колебаний около прямой х = ! (см. рис. 10, в (и)) и возвращается к неподвижной точке О. Такой траектории отвечает режим возбуждения одиночного берста с нерегулярной последовательностью спайков (см. рис. 10, б (и)).
(ш) Если амплитуда внешнего стимула достаточна для преодоления второго порога, то точка х = хе, у = уе попадает в область притяжения устойчивой инвариантной кривой
Wi(£) = {(х,у): у = ¥(х) + — в + х > 1тах} .
и траектория Е приходит в слой медленных движений в окрестности этой кривой (см. рис. 10, в (ш)). Движение траектории Е в окрестности Ш^(£) заканчивается в окрестности точки (х & 1тах, у = ¥(1тах) — в) и далее продолжается в соответствии с быстрой системой, близкой к (8), с уо е Уь Достигнув окрестности ^®(£), траектория Е продолжает свое движение в этой окрестности и асимптотически приближается к неподвижной точке О. Такому поведению траектории Е соответствует режим генерации одиночного спайка (см. рис. 10, б (ш)).
0 150 300 450 п 0.1 0.4 0.7 1.0 х
б в
Рис. 10. Режим однократной нейронной активности для трех внешних стимулов: а - внешняя стимуляция; б - отклик системы на внешний стимул; в - фазовая плоскость системы (1). Параметры 3 = 0.15, е = 0.002, в = 0.07, d = 0.7, а = 0.9
4.2. Колебательные режимы нейронной активности.
4.2.1. Замкнутая инвариантная кривая и подпороговые колебания. При
нарушении неравенства (40) неподвижная точка О теряет устойчивость. Нетрудно видеть, что при 3 = Н(а, е) ее мультипликаторы равны
(1,2 = е-, ю = у е - 3.
Следовательно, смена устойчивости точки О происходит через бифуркацию Андро-нова-Хопфа-Неймарка-Сакера [28]. Приводя систему (1) при ю = 0, п, 2п/3, п/2 в окрестности О к соответствующей этой бифуркации нормальной форме, устанавливаем, что первая ляпуновская величина задается следующим образом:
3 1 - а + а2 + Зе
¿1 \.1=Ща,е) - Щ—) К—)
При е < 4 первая ляпуновская величина является отрицательной и, следовательно, смена устойчивости точки О сопровождается рождением устойчивой замкнутой инвариантной кривой С^ь (рис. 11, а). Колебания, соответствующие С\ь, происходят ниже порога возбуждения потенциала действия и потому называются подпороговы-ми (рис. 11, б).
4.2.2. Хаотический «двухрусловый» аттрактор и хаотические спайко-вые колебания. Рассмотрим релаксационную (е << 1) динамику системы (1) в случае, когда неподвижная точка О является неустойчивой, то есть 3 > Н(а, е), и выполнены также следующие условия на параметры
Р(3тш) > Р(й) - в, (41)
Р(3тах) > Р(й). (42)
В этом случае в окрестности прямой разрыва х = й инвариантная кривая W21(е) разделяет быстрые траектории системы (8) на два потока (рис. 11, в). Первый поток образуют траектории, совершающие колебаний в окрестности х = й. Их динамика близка к динамике системы (8) на инвариантном интервале I(уо) с уо € Уь (см. раздел 2). Второй поток состоит из траекторий, преодолевших второй порог и движущихся затем в окрестности кривой WS(е) (рис. 11, д). При этом траектория системы (1), стартующая, например, из окрестности W:Ls(е), хаотически переходит из одного потока в другой. В результате такой динамики на фазовой плоскости формируется хаотический аттрактор А^. Аттрактору А^ соответствует режим нейронной активности, при котором на фоне подпороговых колебаний возникают хаотические спайковые колебания. (рис. 11, г).
4.2.3. Замкнутая инвариантная кривая и периодические спайковые колебания. Пусть параметры системы (8) удовлетворяют тем же условиям, что и в п. 4.2.2, за исключением неравенства (41). Кроме того, предположим, что выполняется неравенство, противоположное (41). В этом случае параметр в относительно мал и траектории с начальными условиями из окрестности WS2(е) не меняют своего направления движения при прохождении прямой разрыва (см. рис. 11, д). Поэтому
в системе (1) возможны движения между слоями медленных движений локализованных около кривых 2(£). В результате таких движений на фазовой плоскости формируется замкнутая инвариантная кривая С8р (см. рис. 11, д), определяющая режим периодической спайковой активности.
Рис. 11. а - Замкнутая инвариантная кривая С^ь, б - подпороговые колебания, соответствующие С^ь, параметры: 3 = 0.115, е = 0.01, в = 0.04, d = 0.5, а = 0.25. в - «Двухрусловый» хаотический аттрактор ЛьЬ, df (Л^) = 1.1544, г - спайк-подпороговые колебания, параметры: 3 = 0.15, е = 0.005, в = 0.018, d = 0.26, а = 0.25. д - Замкнутая инвариантная кривая Сзр, е - периодические спайковые колебания, параметры: 3 = 0.15, е = 0.005, в = 0.04, d = 0.5, а = 0.25
Заключение
В работе предложена новая дискретная модель нейрона в виде двумерного точечного отображения, построенная на основе дискретной версии системы ФитцХью - Нагумо. Проведено исследование динамики модели и показано существование как хаотических, так и регулярных аттракторов различных видов. Построена зависимость фрактальной размерности хаотических аттракторов от основных управляющих параметров. В основе динамического механизма хаотизации колебаний модели лежат свойства широко известного отображения типа Лоренца. С одной стороны, предложенная модель обладает относительно простой и ясной структурой, а с другой - описывает достаточно большое число сложных режимов нейронной активности и, на наш взгляд, представляется достаточно универсальной.
Несмотря на идеализации, принятые при построении модели, описываемые ею режимы находятся в хорошем качественном соответствии с их прототипами, полученными в реальных нейрофизиологических экспериментах. Например, подпоро-говые квазисинусоидальные колебания типичны для нейронов нижних олив [31], входящих в оливомозжечковую систему, которая в значительной степени определяет процесс координации движений, осуществляемый нервной системой. Для некоторых видов нейронов нижних олив морской свинки характерна [32] спайковая активность на фоне хаотических подпороговых колебаний (см. рис. 11, г). Спайк-берстовую активность демонстрируют многие нейроны, в частности, пирамидальные нейроны гиппокампа [34], таламические нейроны [35] и другие. Модель (1) позволяет не только воспроизводить широкий спектр спайк-берстовых колебаний, но и осуществлять «настройку» специфических свойств этих колебаний, характерных для различных типов нейронов.
Работа выполнена при поддержке Российского фонда фундаментальных исследований (грант 06-02-16137) и гранта Президента РФ для поддержки ведущих научных школ (НШ-7309.2006.2).
Библиографический список
1. Kandel E.R., Schwartz J.H., Jessell T.M. Principles of neural science. Prentice-Hall Int. Inc., 1991.
2. Rabinovich M.I., Varona P., Selverston A.I., Abarbanel H.D.I. Dynamical principles in neuroscience // Reviews of Modern Physics. 2006. Vol. 78, 4. P. 1213.
3. Hodgkin A.L., Huxley A.F. A quantitative description of membrane current and application to conduction and excitation in nerve // J. Physiol. Vol. 117. P. 500.
4. Потапов А.Б., Али М.К. Нелинейная динамика обработки, используемая в нейронных сетях // В кн. Новое в синергетике: взгляд в третье тысячелетие. М.: Наука, 2002.
5. FitzHugh R. Mathematical models of the threshold phenomena in the nerve membrane // Bull. Math. Biophys. 1955. Vol. 17. P. 257.
6. Hindmarsh J.L., Rose R.M. A model of neuronal bursting using three coupled first order differential equations // Philos. Trans R. Soc. London, Ser. B221. 1984. P. 87.
7. Morris C., Lecar H. Voltage oscillations in the barnacle giant muscle fiber // Biophys. J. 1981. Vol. 25. P. 87.
8. Chialvo D.R. Generic excitable dynamics on the two-dimensional map // Chaos Solutions Fract. 1995. Vol. 5. P. 461.
9. Kinouchi O., Tragtenberg M. Modeling neurones by simple maps // Int. J. Bifurcation Chaos. 1996. Vol. 6, № 12a. P. 2343.
10. Kuva S., Lima G., Kinouchi O. Tragtenberg M., Roque A. A minimal model for excitable and bursting element // Neurocomputing 2001. Vol. 38-40. P. 255.
11. de Vries G. Bursting as an emergent phenomenon in coupled chaotic maps // Phys. Rev. E. 2001. Vol. 64. P. 051914.
12. Bazhenov M., Rulkov N.F., Fellous J.-M., Timofeev I. Role of network dynamics in shaping spike timing reliability // Phys. Rev. E. 2005. Vol. 72. P. 041902.
13. Rulkov N.F., Timofeev I., Bazhenov M. Oscillations in large-scale cortical networks: map-based model // J. of Computational Neuroscince. 2004. Vol. 17. P. 203.
14. Su H., Alroy G., Kirson E.D., Yaari Y Extracellular calcium modulates persistent sodium current-dependent burst-firing in hippocampal pyramidal neurones // Journal of Neuroscience. 2001. Vol. 21. P. 4173.
15. Gray C.M., Mc Cormick D.A. Chattering cells: superficial pyromidal neurones contributing to the generation of synchronous oscillations in the visual cortex // Science. 1996. 274 (5 284). P. 109.
16. Connors B.W., Gutnick M.J. Intrinsic firing patterns of diverse neocortical neurones // Trends in Neuroscience. 1990. Vol. 13. P. 99.
17. Lisman J. Bursts as a unit of neural information: making unreliable synapses reliable. Trends in Neuroscience. 1997. Vol. 20. P. 38.
18. Izhikevich E.M., Desai N.S., Walcott E.C., Hoppensteadt F.C. Bursts as s unit of neural information: selective communication via resonance // Trends in Neuroscience. 2003. Vol. 26. P. 161.
19. Izhikevich E.M., Hoppensteadt F. Classification of bursting mappings // Int. J. Bifurcation and Chaos. 2004. Vol. 14, № 11. P. 3847.
20. Cazellis B., Courbage M., Rabinovich M. Anti-phase regularization of coupled chaotic maps modeling bursting neurons // Europhysics Letters. 2001. Vol. 56 (4). P. 504.
21. Rulkov N.F. Modeling of spiking-bursting neural behavior using two-dimensional map // Phys. Rev. E. 2002. Vol. 65. P. 0.41922.
22. Shilnikov A.L., Rulkov N.F. Origin of chaos in a two dimensional map modeling spiking-bursting neural activity // Int. J. Bifurc. Chaos. 2003. Vol. 13, № 11. P. 3325.
23. Rulkov N.F. Regularization of synchronized chaotic bursts // Phys. Rev. Lett. 2001. Vol. 86. p. 183.
24. Shilnikov A.L., Rulkov N.F. Subthreshold oscillations in a map-based neuron model // Physics Letters A. 2004. Vol. 328. P. 177.
25. Tanaka G. Synchronization and propagation of bursts in networks of coupled map neurons. Chaos, 16, 2006, 013113
26. Llinas R. I of vortex. From Neurones to Self.The MIT Press, Massachusettes, 2002.
27. Llinas R. and Yarom Y. Oscillatory properties of guinea-pig inferior olivary neurones and their pharmacological modulation: An in vitro study // J. Physol., Lond. 1986. Vol. 315. P. 569.
28. Afraimovich VS., Sze-Bi Hsu. Lectures on Chaotic Dynamical Systems. American Mathematical Society. Int. Press, 2003, 354 p.
29. Afraimovich VS., Shilnikov L.P. Strange attractors and quasiattractors // In the book «Nonlinear Dynamics and Turbulence» / eds. G.I. Barenblatt, G. Iooss, D.D. Joseph. Pitam, Boston, 1983. P. 1.
30. Arnold V.I., Afraimovich KS., Ilyashenko Yu.S., Shilnikov L.P. Bifurcation Theory, Dyn. Sys. V. Encyclopaedia Mathematics Sciences, Springer, Berlin, 1994.
31. Llinas R. Rebound excitation as the physiological basis for tremor: a biological study of the oscillatory properties of mammalian central neuron in vitro. Movements Disorders: Tremor / Eds L.J. Findley and R. Capildeo. London: Macmillan, 1984. P. 135.
32. Bernardo L.S., Foster R.P. Oscillatory behavior in inferior olive neurones: mechanism, modulation, cell aggregates // Brain Res. Bull. 1986. Vol. 17. P. 773.
33. Traub R.D., Jefferys J.G.R., Whittington M.A. Fast Oscillations in Cortical Circuits. The MIT Press, Massachusetts, 1999.
34. R.S.K. Wang and Prince D.A. Afterpotential generation in hippocampal pyramidal cells // J. Neurophysiol. 1981. Vol. 45. P. 86.
35. Deschenes M., Roy J.P. and Steriade M. Thalamic bursting mechanism: an inward slow current revealed by membrane hyperpolarization // Brain Res. 1982. Vol. 239. P. 289.
Институт прикладной физики РАН, Поступила в редакцию 05.06.2007
Н.Новгород
MAP-BASED MODEL OF THE NEURAL ACTIVITY
V.I. Nekorkin, L.V. Vdovin
A two-dimensional model exhibiting the chaotic spiking-bursting activity of real neurons is proposed. The model is given by the discontinuous two-dimensional map. It is constructed on the basis of the discrete modification of the FitzHugh-Nagumo model and one-dimensional Lorenz-type map. We have studied the dynamics of the system, found the conditions on the parameters under which chaotic attractor exists. The structure and properties of the attractor is studied. This attractor mimics spiking-bursting oscillations. We have also showed that map is applicable for simulation of other regimes of neural oscillatory activity such as subthreshold oscillations and periodic and chaotic spiking or it could be used for modeling of threshold excitability property of the neurons.
Некоркин Владимир Исаакович - родился на Украине (1948), окончил радиофизический факультет Нижегородского государственного университета (1971). Защитил диссертацию на соискание ученой степени кандидата физико-математических наук (1981, ННГУ) и доктора физико-математических наук (1992, СГУ). С 1971 года работал в Нижегородском институте прикладной математики, с 1983 по 1986 год - в Нижегородском политехническом институте, с 1986 года и по настоящее время работает на радиофизическом факультете ННГУ, профессор. Область научных интересов - теория колебаний и волн, структурообразование и нелинейные волны в непрерывных и дискретных средах, синхронизация и пространственно-временной хаос. Имеет более 100 научных публикаций в указанных направлениях, в том числе двух монографий (в соавторстве). В качестве приглашенного профессора читал лекции в Мадридском и Стэнфордском университетах. Работал в Калифорнийском университете. E-mail: vnekorkin@appl.sci-nnov.ru
Вдовин Лев Вячеславович - родился в 1983. Окончил радиофизический факультет Нижегородского государственного университета (2006). Аспирант первого года обучения ИПФ РАН. Область научных интересов - теория нелинейных колебаний и волн, нейродинамика.