Модель Фіцх'ю-Нагумо

1. Походження та спрощення від Ходжкіна-Хакслі

Модель Ходжкіна-Хакслі (ХХ) (1952) точно відтворює потенціал дії гігантського аксона кальмара за допомогою чотирьох пов'язаних ОДР: мембранна напруга V та змінні воріт m, h, n для динаміки натрієвих і калієвих каналів. Однак складність моделі ускладнює геометричну інтуїцію — 4D фазовий простір не має простої картини.

Річард Фіцх'ю (1961) помітив, що в динаміці ХХ змінна m (швидка активація натрію) відстежує V майже миттєво, тоді як h і n повільніші та корельовані. Він запропонував 2D карикатуру з використанням швидкої змінної мембранного потенціалу v з кубічною нульовою кривою (замінюючи швидку динаміку V і m) та повільної відновлювальної змінної w (замінюючи комбіновану повільну динаміку h і n).

Незалежно, Дж. Нагумо, С. Арімото та С. Йошізава (1962) побудували електронну аналогову схему, що реалізовує ті ж рівняння, демонструючи генерацію потенціалів дії за допомогою тунельного діода.

2. Рівняння та параметри

Стандартні рівняння Фіцх'ю-Нагумо:

Кубічний доданок v − v³/3 створює N-подібну нульову криву за v (dv/dt = 0), що є геометричним серцем моделі. Повільне лінійне відновлювальне рівняння утримує траєкторії обмеженими та генерує зворотній шлях потенціалу дії.

Формулювання через τ

Альтернативна параметризація явно використовує сталі часу:

3. Аналіз фазової площини

Фазова площина — це площина (v, w). Кожна точка площини представляє стан системи; часова еволюція прокладає траєкторію через цю площину. Аналіз фазової площини запитує: які якісні траєкторії для всіх початкових умов?

Ключові інструменти аналізу фазової площини:

• Нульові криві: криві, де dv/dt = 0 (нульова крива за v) і dw/dt = 0 (нульова крива за w); траєкторії перетинають нульову криву за v вертикально, а нульову криву за w горизонтально.
• Нерухомі точки: перетини двох нульових кривих; де система знаходиться у спокої.
• Поле векторів: напрямок (dv/dt, dw/dt) у кожній точці; стрілки показують, куди рухаються траєкторії.
• Граничні цикли: замкнені траєкторії, що є стійкими (притягувальними); відповідають periodical-коливанням (повторним спайкам).
• Сепаратриса: межа у фазовому просторі між областями, що збігаються до нерухомої точки, і тими, що спочатку виконують великий вихід (потенціал дії).

Сила фазової площини ФН полягає в тому, що один потенціал дії з'являється як траєкторія, яка: стартує поблизу нерухомої точки спокою, виштовхується за сепаратрису струмом I_ext, мчить вгору по швидкій кубічній нульовій кривій за v, повільно проходить верхнє плато, спадає по правій гілці, повільно проходить нижнє плато і, нарешті, спірально повертається до спокою.

4. Нульові криві та нерухомі точки

Прирівнення кожного рівняння до нуля визначає нульові криві:

По мірі зростання I_ext нульова крива за v зсувається вгору. Нерухома точка ковзає по лівій гілці кубічної кривої. Коли I_ext перевищує критичне значення I_c ≈ 0,34, нерухома точка переходить на нестійку середню гілку — запускаючи надкритичну біфуркацію Хопфа та початок повторних спайків.

Стійкість нерухомої точки

Якобіан у нерухомій точці (v*, w*):

5. Граничний цикл та біфуркація Хопфа

Коли I_ext перевищує поріг Хопфа I_c, стійка нерухома точка стає нестійкою і з'являється стійкий граничний цикл — нейрон генерує потенціали дії periodically без зовнішнього введення. Біфуркація Хопфа у ФН є надкритичною: амплітуда циклу зростає безперервно від нуля при перетині I_ext порогу I_c.

Граничний цикл зростає з I_ext: більший струм → коротший період → вища частота спайків. Крива f-I (частота спайків від введеного струму) стартує з кінцевої частоти на порозі Хопфа (збудливість типу II) або зростає безперервно від 0 (збудливість типу I), залежно від типу біфуркації.

6. Класи I і II збудливості

Нейрони класифікуються за тим, як частота їх спайків залежить від прикладеного струму. Цю класифікацію запропонував Ходжкін (1948), і вона безпосередньо читається з фазової площини ФН:

Клас I — біфуркація SNIC

Нейрон Класу I переходить від спокою до повторних спайків через біфуркацію сідлового вузла на інваріантному колі (SNIC). Ключова ознака: частота спайків стартує з нуля і безперервно зростає зі струмом — f(I_c) = 0. Нейрони Класу I можуть інтегрувати входи з часом (поведінка фільтру низьких частот).

Клас II — біфуркація Хопфа

Нейрон Класу II переходить через біфуркацію Хопфа. Частота спайків на порозі є ненульовою — нейрон стрибкоподібно переходить від мовчання до спайків з мінімальною частотою. Стандартна ФН є Класом II. Нейрони Класу II є кращими резонаторами (поведінка смугового фільтра) і найсильніше реагують на входи поблизу природної частоти Хопфа.

7. Пов'язані осцилятори ФН

Два осцилятори ФН, пов'язані лінійним дифузійним (щілинний контакт) доданком, демонструють різноманітну поведінку синхронізації залежно від сили зв'язку та різниці власних частот:

Просторово розширена ФН — реакційно-дифузійні хвилі

Додавання просторової дифузії v уздовж 1D аксона дає рівняння ФН у часткових похідних:

Просторово розширена ФН є однією з основних моделей серцевої аритмії. Поломка спіральних хвиль у 2D ФН представляє перехід від реентерантної тахікардії до фібриляції. Параметричні регіони, що підтримують спіральні хвилі проти планарних хвиль проти нерегулярних візерунків, прямо відповідають клінічній класифікації аритмій.

8. JavaScript-реалізація

// Розв'язувач ОДР Фіцх'ю-Нагумо за допомогою RK4
// Будує v(t) та w(t); малює траєкторію на фазовій площині

const params = {
  a:    0.7,
  b:    0.8,
  eps:  0.08,
  Iext: 0.5,    // спробуйте 0 (мовчання), 0.34 (поріг), 0.5–1.5 (повторні спайки)
};

const dt    = 0.1;    // мс крок часу
const T_end = 300;    // мс загальний час симуляції
const tSteps = Math.floor(T_end / dt);

// Початкові умови (поблизу спокою)
let v = -1.2, w = -0.6;

// Права частина ФН
function fhnRHS(v, w) {
  const dv = v - (v * v * v) / 3 - w + params.Iext;
  const dw = params.eps * (v + params.a - params.b * w);
  return [dv, dw];
}

function rk4Step(v, w, dt) {
  const [k1v, k1w] = fhnRHS(v, w);
  const [k2v, k2w] = fhnRHS(v + dt/2*k1v, w + dt/2*k1w);
  const [k3v, k3w] = fhnRHS(v + dt/2*k2v, w + dt/2*k2w);
  const [k4v, k4w] = fhnRHS(v + dt*k3v, w + dt*k3w);
  return [
    v + dt/6 * (k1v + 2*k2v + 2*k3v + k4v),
    w + dt/6 * (k1w + 2*k2w + 2*k3w + k4w),
  ];
}

// Генерація траєкторії
const traj = [];
for (let i = 0; i < tSteps; i++) {
  traj.push([v, w]);
  [v, w] = rk4Step(v, w, dt);
}

// Малювання нульових кривих на canvas фазової площини
function drawNullclines(ctx, W, H, scale) {
  const { a, b, Iext } = params;
  const vRange = [-2.5, 2.5];

  // Нульова крива за v: w = v - v³/3 + Iext
  ctx.strokeStyle = '#6366f1';
  ctx.lineWidth   = 2;
  ctx.beginPath();
  for (let px = 0; px <= W; px++) {
    const vv = vRange[0] + (px / W) * (vRange[1] - vRange[0]);
    const ww = vv - (vv*vv*vv)/3 + Iext;
    const y = H/2 - ww * scale;
    px === 0 ? ctx.moveTo(px, y) : ctx.lineTo(px, y);
  }
  ctx.stroke();

  // Нульова крива за w: w = (v + a) / b
  ctx.strokeStyle = '#f59e0b';
  ctx.beginPath();
  for (let px = 0; px <= W; px++) {
    const vv = vRange[0] + (px / W) * (vRange[1] - vRange[0]);
    const ww = (vv + a) / b;
    const y = H/2 - ww * scale;
    px === 0 ? ctx.moveTo(px, y) : ctx.lineTo(px, y);
  }
  ctx.stroke();
}

// Малювання траєкторії на фазовій площині
function drawTrajectory(ctx, traj, W, H, scale) {
  const vRange = [-2.5, 2.5];
  ctx.strokeStyle = '#34d399';
  ctx.lineWidth   = 1.5;
  ctx.beginPath();
  traj.forEach(([vv, ww], i) => {
    const x = (vv - vRange[0]) / (vRange[1] - vRange[0]) * W;
    const y = H/2 - ww * scale;
    i === 0 ? ctx.moveTo(x, y) : ctx.lineTo(x, y);
  });
  ctx.stroke();
}

// Побудова кривої f-I: частота спайків від струму
function buildFICurve(IValues, T = 500) {
  return IValues.map(Iext => {
    params.Iext = Iext;
    let v = -1.2, w = -0.6;
    const steps = Math.floor(T / dt);
    const traj = [];
    for (let i = 0; i < steps; i++) {
      traj.push([v]); [v, w] = rk4Step(v, w, dt);
    }
    // Рахуємо спайки лише у другій половині (ігноруємо перехідний процес)
    let spikes = 0, above = false;
    traj.slice(steps / 2).forEach(([vv]) => {
      if (!above && vv > 0)  { spikes++; above = true; }
      if ( above && vv < 0)  { above = false; }
    });
    return (spikes / (T / 2)) * 1000; // Гц
  });
}