Нейрон Ходжкіна-Гакслі — потенціал дії з іонних струмів

1. Нейрон як електричне коло

Плазматична мембрана розділяє заряд через свою товщину ~7 нм, діючи як конденсатор з ємністю C_m ≈ 1 мкФ/см². Іонні канали — білкові пори, що відкриваються й закриваються у відповідь на напругу — діють як потенціал-залежні резистори. Рушійна сила для кожного виду іонів — це різниця між мембранним потенціалом V та нернстівським (рівноважним) потенціалом E.

У своїх експериментах із фіксацією потенціалу на гігантських аксонах кальмара Ходжкін та Гакслі прямо виміряли провідності g_Na та g_K і отримали числові параметри, що відповідають даним. Набір параметрів нижче відтворює оригінальні результати 1952 року.

2. Чотири рівняння Ходжкіна-Гакслі

Повна система — це чотири зв'язані ЗДР: одне для напруги V і три для ворітних змінних m, h, n — кожна з яких представляє ймовірність того, що певні ворота в білку іонного каналу перебувають у відкритому стані.

Кожна ворітна змінна підкоряється кінетиці першого порядку: вона релаксує до свого стаціонарного значення m_∞(V) = α_m/(α_m + β_m) зі сталою часу τ_m(V) = 1/(α_m + β_m). Швидкості α та β — емпіричні функції V, підігнані до даних фіксації потенціалу.

Член m³ у I_Na відображає те, що натрієвий канал потребує одночасного відкриття трьох незалежних активаційних воріт. Множник h — це ворота інактивації: вони починають відкритими (~0,6), але закриваються в міру деполяризації мембрани, перекриваючи струм Na⁺ і завершуючи потенціал дії. Член n⁴ у I_K відображає чотири однакові субодиниці.

3. Ворітні змінні m, h, n — альфа/бета-кінетика

Функції швидкості нижче дійсні, коли напруга V виміряна в мВ. Вони містять усувні сингулярності за певних напруг (коли знаменник → 0), які слід обробляти чисельно за правилом Лопіталя.

// Усі напруги в мВ; швидкості в мс⁻¹
// V — це відхилення від спокою: використовуйте V_shifted = V_membrane − (−65)

function alpha_m(V) {
  const dv = 10 - V;
  return Math.abs(dv) < 1e-7
    ? 1                                    // границя за Лопіталем
    : 0.1 * dv / (Math.exp(dv / 10) - 1);
}
function beta_m(V)  { return 4    * Math.exp(-V / 18); }

function alpha_h(V) { return 0.07 * Math.exp(-V / 20); }
function beta_h(V)  { return 1    / (Math.exp((30 - V) / 10) + 1); }

function alpha_n(V) {
  const dv = 10 - V;
  return Math.abs(dv) < 1e-7
    ? 0.1
    : 0.01 * dv / (Math.exp(dv / 10) - 1);
}
function beta_n(V)  { return 0.125 * Math.exp(-V / 80); }

4. Чисельне інтегрування з RK4

Потенціал дії триває ~1 мс, а найшвидша ворітна змінна (m) має τ_m ≈ 0,1 мс на порозі. Безпечний крок часу — dt = 0,01 мс. Метод Ейлера нестійкий для цієї системи за такого кроку; RK4 і стійкий, і має четвертий порядок точності.

function derivatives(V, m, h, n, Iext) {
  const INa = 120 * m**3 * h * (V - 50);
  const IK  =  36 * n**4     * (V + 77);
  const IL  = 0.3             * (V + 54.4);
  return {
    dV: Iext - INa - IK - IL,             // Cm = 1
    dm: alpha_m(V) * (1 - m) - beta_m(V) * m,
    dh: alpha_h(V) * (1 - h) - beta_h(V) * h,
    dn: alpha_n(V) * (1 - n) - beta_n(V) * n,
  };
}

function stepRK4(V, m, h, n, Iext, dt) {
  const k1 = derivatives(V, m, h, n, Iext);
  const k2 = derivatives(V + dt/2*k1.dV, m + dt/2*k1.dm,
                          h + dt/2*k1.dh, n + dt/2*k1.dn, Iext);
  const k3 = derivatives(V + dt/2*k2.dV, m + dt/2*k2.dm,
                          h + dt/2*k2.dh, n + dt/2*k2.dn, Iext);
  const k4 = derivatives(V + dt*k3.dV, m + dt*k3.dm,
                          h + dt*k3.dh, n + dt*k3.dn, Iext);
  return {
    V: V + dt/6 * (k1.dV + 2*k2.dV + 2*k3.dV + k4.dV),
    m: m + dt/6 * (k1.dm + 2*k2.dm + 2*k3.dm + k4.dm),
    h: h + dt/6 * (k1.dh + 2*k2.dh + 2*k3.dh + k4.dh),
    n: n + dt/6 * (k1.dn + 2*k2.dn + 2*k3.dn + k4.dn),
  };
}

5. Фази потенціалу дії

Подайте короткий струмовий імпульс (I_ext = 10 мкА/см² протягом 0,5 мс), і мембранний потенціал пройде через чотири окремі фази:

Властивість «усе або нічого» виникає з порогу: слабкі стимули спричиняють підпорогову відповідь, що згасає; стимули вище ~−55 мВ запускають повний спайк. Це емерджентна властивість нелінійного зв'язку між m, h та V.

6. Повна реалізація на JavaScript

class HHNeuron {
  constructor() {
    this.V  = -10;   // відхилення від спокою: −65 мВ абсолютно → 0 тут
    this.m  = this._minf(this.V);
    this.h  = this._hinf(this.V);
    this.n  = this._ninf(this.V);
    this.dt = 0.01;   // мс
  }

  _alpha_m(V) { const d=10-V; return Math.abs(d)<1e-7 ? 1 : 0.1*d/(Math.exp(d/10)-1); }
  _beta_m(V)  { return 4    * Math.exp(-V/18); }
  _alpha_h(V) { return 0.07 * Math.exp(-V/20); }
  _beta_h(V)  { return 1    / (Math.exp((30-V)/10)+1); }
  _alpha_n(V) { const d=10-V; return Math.abs(d)<1e-7 ? 0.1 : 0.01*d/(Math.exp(d/10)-1); }
  _beta_n(V)  { return 0.125 * Math.exp(-V/80); }

  _minf(V) { const a=this._alpha_m(V); return a/(a+this._beta_m(V)); }
  _hinf(V) { const a=this._alpha_h(V); return a/(a+this._beta_h(V)); }
  _ninf(V) { const a=this._alpha_n(V); return a/(a+this._beta_n(V)); }

  _deriv(V, m, h, n, I) {
    return {
      dV: I - 120*m**3*h*(V-50) - 36*n**4*(V+77) - 0.3*(V+54.4),
      dm: this._alpha_m(V)*(1-m) - this._beta_m(V)*m,
      dh: this._alpha_h(V)*(1-h) - this._beta_h(V)*h,
      dn: this._alpha_n(V)*(1-n) - this._beta_n(V)*n,
    };
  }

  step(Iext = 0) {
    const dt = this.dt;
    let {V, m, h, n} = this;
    const k1 = this._deriv(V, m, h, n, Iext);
    const k2 = this._deriv(V+dt/2*k1.dV, m+dt/2*k1.dm, h+dt/2*k1.dh, n+dt/2*k1.dn, Iext);
    const k3 = this._deriv(V+dt/2*k2.dV, m+dt/2*k2.dm, h+dt/2*k2.dh, n+dt/2*k2.dn, Iext);
    const k4 = this._deriv(V+dt*k3.dV, m+dt*k3.dm, h+dt*k3.dh, n+dt*k3.dn, Iext);
    this.V = V + dt/6*(k1.dV+2*k2.dV+2*k3.dV+k4.dV);
    this.m = m + dt/6*(k1.dm+2*k2.dm+2*k3.dm+k4.dm);
    this.h = h + dt/6*(k1.dh+2*k2.dh+2*k3.dh+k4.dh);
    this.n = n + dt/6*(k1.dn+2*k2.dn+2*k3.dn+k4.dn);
    return this.V + (-65);  // перетворення в абсолютні мВ
  }
}

// Запускаємо 30-мс симуляцію з імпульсом 10 мкА/см² від t=1 до t=2 мс
const neuron = new HHNeuron();
const trace  = [];
for (let t = 0; t <= 30; t += 0.01) {
  const I = (t >= 1 && t <= 2) ? 10 : 0;
  trace.push({ t, V: neuron.step(I) });
}
// trace тепер містить [{ t: 0, V: -65 }, ..., { t: 5.2, V: +38 }, ...]

7. 2D-збудлива тканина та спіральні хвилі

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

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

Для простої реалізації кожна клітина — окремий HH-нейрон. Зв'язок додається як додатковий струм: I_diff = D × Σ(V_neighbour − V_self). Сітка 100×100 з D ≈ 0,5 см²/с та dt = 0,05 мс комфортно працює в реальному часі на сучасному обладнанні.

8. Розширення та застосування

Фітцг'ю-Нагумо: спрощена модель

Для великомасштабних симуляцій тканини повна 4-змінна модель HH витратна. Модель Фітцг'ю-Нагумо зводить її до 2 змінних (швидкий активатор u, повільна змінна відновлення w), зберігаючи водночас збудливість і структуру біфуркацій. Вона працює у 4 рази швидше й схоплює спіральні хвилі, біжучі імпульси та коливання.

Багатокомпартментні нейрони

Реальні нейрони не є точковими об'єктами. Дендритні дерева складаються з багатьох циліндричних компартментів, зв'язаних осьовим опором. Кабельне рівняння (компартментна модель Ралла) розглядає кожен сегмент як HH-подібну ділянку, з'єднану із сусідами провідністю g_c = A/(R_a·Δx), даючи лінійну систему ЗДР, яку для стійкості можна розв'язувати неявними методами (Кранка-Ніколсон).

Стохастичні іонні канали

Окремі канали дискретні — вони або відкриваються, або закриваються — і їхня поведінка в малих нейронах імовірнісна. Стохастичний HH замінює ЗДР для m, h, n біноміальними оновленнями ланцюга Маркова, породжуючи канальний шум і спонтанне спайкування за підпорогових струмів. Це реалізується з переходами Binomial(N_channel, α·dt).