第三章 研究方法
多數的神經元當接收到外界刺激後，會產生一系列的電流脈衝，稱之為動作 電位。通常這些動作電位在傳遞時，會保持一個穩定的速度與振幅。動作電位在 神經元之間的溝通與訊息的傳遞上，扮演著很重要的角色，所以瞭解動作電位產 生的機制與如何進行傳導，細胞膜周遭的離子濃度變化與動作電位的相關聯性…
等，一直是眾人研究的目標。
3-1 The Hodgkin-Huxley Model
在 1952 年，Alan Hodgkin 和 Andrew Huxley 首次說明了神經組織產生動作 電位與傳遞動作電位的整個離子機制。並且提出了一個簡單的模型，明確的描述 動作電位與離子電流之間的關係，稱為Hodgkin-Huxley model ( HH model )。[3]
3-1-1 HH model
dt t C dV t I t
Im ionic m ( )
) ( )
( = + (1)
Im ： the total membrane current density Iionic ： the total ionic current density
V ： the displacement of the membrane potential from its resting value Cm ： the membrane capacity per unit area
HH model 是一個屬於細胞尺度的離子系統模型，單純的考慮流過細胞膜表 面的電流與其周圍的離子電流變化；其中離子電流包含了兩項主要的元素，分別 是鈉離子電流(INa)與鉀離子電流(IK)，至於其他效應較微弱的離子電流，則採用 一組等效電流(Ileak)來表示。[2,3,17]
Iionic =INa +IK +Ileak (2) INa ： the sodium current
IK ： the potassium current
Ileak ： the effective current of other ions
在HH model中，細胞膜上的動作電位與離子電流的關係，可以由Fig. 15 簡 單的電路示意圖來表示。由圖可以看出整個HH model的基本概念，簡單的透過 四組並聯的電流路徑，運用離子電流的改變來調控動作電位的變化，其中INa及IK 兩組離子電流與時間、膜電位有關，是整個模型最主要的調變項，另兩組則是屬 於被動項。
Ileak
Cm
Fig.15 : HH model 的基本電路示意圖，將細胞膜兩側的膜電位與離子電流的關 係，用四個並聯的電路來表示，其中INa及IK 是主要的調控項，其他兩組為被動 項。[2]
3-1-2 離子電流：
離子電流在HH model 中扮演著相當重要的角色，尤其是INa及IK 這兩項，是 控制整體膜電位變化的關鍵。而對於每一個單獨的離子電流來說，都必須遵守 Ohm’s law；所以具有共同的形式：
Ii(t)= g i(V(t),t)(V(t)−Ei) , i = K , Na , leak ； (3) gi ： the voltage-dependent ionic conductance
Ei ： the ionic reversal potential
為了描述離子電流中gi 的動態變化機制，Hodgkin和Huxley引用了一個因 子，稱為gating particle，他們利用這個因子來調控離子通道電導率的變化，藉此 來 描 述 離 子 通 道 的 活 性 與 狀 態 的 變 化 。 其 中gating particle 又 可 以 區 分 為
activating 和 inactivating 兩種。以下我們將分別對鉀離子電流與鈉離子電流做 介紹。
1.鉀離子電流 (The potassium current IK )
在 HH model 中，鉀離子電流的變化調控較為單純，只有一個 activating gating particle --- n 來調控其變化。根據 Hodgkin 和 Huxley 假設：
n4 K K
g
g = (4)
n ： the activating gating particle
K
g ： the maximal potassium conductance
)
4(
K K
K n V E
I = g −
(5) EK ： the potassium battery
假設只有兩種狀態存在(open 和 close)，則若 open 狀態的機率為 n，close 狀 態的機率則為(1-n)。那麼鉀離子通道在這兩種狀態的變化可以簡單的表示為：
βn n 1- n
αn
αn , βn ： the voltage-dependent rate constants
我們可以將上述的變化關係式數學化，用一個一階 ODE 來表示：
V n V n dt
dn
n
n( )(1 ) β ( )
α − −
= (6)
解 eq. (6)，我們可以得到：
) ( 0
)
(
nt
e n n n
n
τ−
−
−
=
∞ ∞ (7)
n n
n n
β α
α
= +
∞ (8)
n n
n α β
τ = +1
(9)
) 0 ( ) 0 (
) 0 (
0
n n
n n
β α
α
= + (10)
根據 Hodgkin 和 Huxley 透過實驗的量測，分析 rate constants 對電位差的關 係後，取近似值所得到的結果：
] 1 [
100 ) 10 (
10 ) (10
−
= −−V
n
e V V
α (11)
( ) 0.125 (80)
V
e V
−
βn = (12)
2.鈉離子電流( The sodium current INa) Hodgkin 和 Huxley 假設：
Na Nam3h
g
g = (13)
m ： the activating gating particle h ： the inactivating gating particle
Na
g ： the maximal sodium conductance
同理：
m V m
dt V dm
m
m( )(1 ) β ( )
α − −
= (14)
h V h
dt V dh
h
h( )(1 ) β ( )
α − −
= (15)
解 eq. (14)、(15)，可以得到
) ( 0
)
(
mt
e m m
m
m
τ−
−
−
=
∞ ∞ (16)) ( 0
)
(
ht
e h h h
h
τ−
−
−
=
∞ ∞ (17)
m m
m m
β α
α
= +
∞ (18)
h h
h h
β α
α
= +
∞ (19)
m m
m α β
τ = +1
(20)
h h
h α β
τ = +1
(21)
) 0 ( ) 0 (
) 0 (
0
m m
m m
β α
α
= + (22)
) 0 ( ) 0 (
) 0 (
0
h h
h h
β α
α
= + (23)
根據實驗量測結果，分析得：
] 1 [
10 ) 25
( )
10 (25
−
= −−V
m
e V V
α (24) ( ) 0.07 (20)
V
h V e
−
α = (25)
18)
4
() (
V
m
V e
−
β =
(26)1 ) 1
(
10 ) (30
+
= −V
h
e
β V (27)
所以整個 HH model 可以寫成：
inj( ) m Kn4(V EK) Nam3h(V ENa) L(V Vrest) dt
C dV t
I = + g − + g − + g −
(28)
參數 數值 參數 數值
Cm 1.0 (μF/cm2) g K 36 (mS/cm2) 120 (mS/cm2) 0.3 (mS/cm2) EK -12 (mV) ENa 115 (mV) Vrest 10.613 (mV) n0 0.3176
m0 0.0530 h0 0.5961
T 6.3 (℃)
Na
g
L
g
Table1. Hodgkin-Huxley model 所用的參數。[2,3,17]
3-1-3 數值分析法
由 eq. (28) 可以知道 HH model 是一個非線性的 ODE，沒有辦法得到精確解。
所以我們使用以下的數值方法來計算：
Euler Method
對任意一個一階 ODE.可以做以下的近似：
f x t f
dt
dx approximate
Δ •
= Δ
⎯
⎯
⎯
⎯
⎯ →
⎯
=
所以，eq. (28) 可以改寫成以下的近似式：
{ [ ( )4( K) Na( t)3 t( mt Na) L( mt rest)]}
t m t K t inj m t m t t
m I n V E m h V E V V
C V t
V+Δ = + Δ − g − + g − + g −
(29)
透過eq.(29)，結合eq.(6)、(14)、(15)，只要我們給定時間為t時的Vm、n、m、
h ，就可以計算出時間為t+Δt 時的Vm值，得到的結果如Fig.16 所示。當然這個 方法的精確度取決於我們選定的Δt之值，越小則近似值越準確[31]。在我們後續 HH model的模擬中，我們均採用Δt = 0.01 msec。
a b
Fig. 16:我們輸入一 0.5 msec , 0.4 nA 的電流脈衝至 HH model，得到(a) 膜電位隨 時間變化關係。(b) gating particle 隨時間變化關係圖。
a b
c d
Fig. 17: 輸入不同的 DC 電流至 HH model 中，觀察膜電位的變化情況。(a) 0.06 nA；(b) 0.065 nA；(c) 0.18 nA；(d) 0.25 nA。
有關DC電流與動作電位的產生情況，由Fig.17(a)可以知道，當輸入的電流值
小於Ithreshold 的時候，不會有動作電位的產生；當輸入電流大於Ithreshold 的時候，
如Fig.17(b)，會開始產生動作電位。若持續增強輸入電流至 0.18 nA，會開始產 生週期性的動作電位，如Fig.17(c)。若再增強電流值可以發現動作電位的產生週 期會逐漸變小，如Fig.17(d)。
3-1-4 HH model 的優缺點
現行神經元模擬被採用的模型很多，但是要研究與探討神經網路的動態變 化，必須要兼具幾項特點。首先，要能夠描述每一個神經元產生興奮時的動態變 化。再來就是要能夠知道神經元之間是如何連接成網路，彼此間如何相互聯繫，
如何傳遞與溝通訊息。
在 2004 年，Izhikevich 對於各種常用的神經元模型做了一系列的整理與比 較。比較的項目主要集中在下面幾個特點：1) 是否具有完整的生物物理意涵。 2) 是否可以完整的表現出神經元應有的生理現象。 3) 能發揮的功效與所需的運算 步驟。[8]其整理的圖表如 Fig.16 所示。
以最常被使用的integrate-and-firing model(IF model) 與 HH model 來說，雖 然IF model 較為簡單，也常常被拿來研究神經網路的性質與現象，但是他在生 物物理上的意義是缺乏的，只能簡單的表現出神經元產生神經衝動的模式，膜電 位超過threshold 以上的現象與資訊完全無法得知，更無法進行相關性質的探討；
不過他的運算極為簡易，所需要的運算步驟也最少，若不著重於神經元的基本特 性探討，適合用來做數目龐大的神經網路性質的模擬運算。
相較於IF model 來說，HH model 這個最早被提出來的模型，大概是計算神 經科學上最重要的模型。透過許多實驗的驗證支持，真實將膜電位與離子電流的 關係，用簡單的四道方程式與十個參數表示出來，可以確實的描述各項神經元表 現出來的行為與生理狀況，詳細探討各種神經元產生的生理反應，其所具有的生 物物理意涵十分豐富。不過，他也有侷限性存在，就是模擬所需要的運算步驟過 多，所以無法做較多數目的神經元網路運算；另外，由於時間的解析度高，故難 以進行長時間性質變化的模擬研究。[19]
a
b
Fig. 18: 各種神經元模型的比較表。(a) 各模型在生物物理上的意義與所需要的 運算步驟。(b) 各神經元模型可以表現出來的神經元生理現象一覽表。[8]
3-2 腦神經網路
我們考慮一個 2-D 的腦神經網路系統，將 N 個神經元隨機的分佈在大小為 L×L 的基底之上；隨著時間的演化，神經元之間會隨機的形成突觸，彼此相互連 結，我們透過刺激其中部分的神經元，讓整個神經網路透過突觸進行溝通聯絡以 及訊息的傳遞，經由學習與記憶效應，以此建構出整個腦神經網路的運轉機制與 生理效應。
3-2-1 突觸連結的形成(connection)
由於神經元生成後，軸突的生長是隨著生長錐的引導，透過環境中生長因子 的誘導，隨機的向外生長，直到尋找到樹突上的刺突，相互接觸後，經由學習與 記憶效應，一連串複雜的生理作用，才能形成穩定而強健的突觸聯結。所以，在 我們後續的模擬中，我們會採用一個簡單的隨機尋找機制，讓神經元彼此之間採 用隨機尋找的模式，產生突觸連結。
首先，我們先給定一個突觸聯結的形成機率，讓神經元透過隨機尋找的模 式，然後以形成機率決定連結的建立與否。此連結形成的機率將依據不同的尋找 機制，與距離的α次方有關。
我們設定任意兩神經元間形成突觸連結的機率：
α
r
ijP = K
， K = 0.005。 (34)rij：任兩神經元間的距離。
α：隨機尋找的機制參數。
a b
c d
Fig. 19:在 N=50, α=0.5，K=0.005 的情況下，突觸聯結形成的情況。(a)所有神經 元的分佈情況；(b)Nc=5；(c)Nc=10；(d)Nc=20。(紅色代表起點，藍色代表終點)。
a b
c d
Fig. 20:在 N=50, α=1，K=0.005 的情況下，突觸聯結形成的情況。(a)所有神經 元的分佈情況；(b)Nc=5；(c)Nc=10；(d)Nc=20。(紅色代表起點，藍色代表終點)。
a b
c d
Fig. 21:在 N=50, α=2，K=0.005 的情況下，突觸聯結形成的情況。(a)所有神經 元的分佈情況；(b)Nc=5；(c)Nc=10；(d)Nc=20。(紅色代表起點，藍色代表終點)。
α是主要聯結機制的調控參數，決定神經元是採用何種機制尋找與建立突觸 聯結。當α=1 時，表示整個尋找機制是隨機的進行聯結，以自己為中心向外尋 找並隨機建立聯結；當α<1 時，則代表整個網路的突觸建立傾向於神經元之間 距離較遠的會產生聯結；而α>1 時，則反過來對於距離較遠的神經元傾向於不 形成聯結。另外，突觸聯結形成的係數 K 跟網路的活躍程度與分佈密度有關，
此處為了簡化與方便，我們將其設為常數。[15]
3-2-2 突觸可塑性的演化(evolution of synaptic plasticity )
Hebbian Learning Rule
“ When an axon of cell A is near enough to excite cell B or repeatedly or consistently takes part in firing it , some growth or metabolic change takes place in one or both cells such that A’s efficiency , as one of the cells firing B , is increase ” [7,16]
Hebbian Learning Rule 是由 Donald Hebb 在 1949 年所提出的，主要在講述 突觸前與突觸後神經元間的興奮時間對於突觸的性質變化影響之關係；根據他當 時提出的假設，決定此一性質的關鍵為神經興奮的時間演化順序關係，而主要的 內涵為一個細胞尺度兼具時空對稱性質的變化機制。其主要的精神為：當突觸前 神經元 i 聯結到突觸後神經元 j，兩者間的突觸聯結其強度變化只與兩神經元產 生神經脈衝的相對時間有關；如果突觸前神經元 i 的神經脈衝產生時間比突觸後 神經元 j 要早，那由 i 連結到 j 的突觸強度就會被增強，反之則會被減弱。[7]
不過在近幾年來，透過相關的實驗證據顯示，Hebbian Learning Rule 需要做 適當的修正。實驗的結果發現突觸可塑性的變化在時間與空間上並不像 Donald Hebb 所提出的情況，實際上量測到的結果在時空上是不對稱的。[21]經由相關的 研究，突觸前與突觸後神經元產生的神經脈衝所誘發之 long-term potentiation (LTP) / long-term depression (LTD)變化的時間視窗「critical windows of spiking timing」被提出，由其中的變化情況，可以系統的看出兩者間的相互關連。
Fig. 22: Critical window for synaptic modulation。 LTP/LTD 的產生與突觸前神經 元跟突觸後神經元兩者的興奮時間是相關聯的。由圖中可以清楚的看出突觸後電 流相對於兩者興奮產生的時間差之變化情況，這種相關聯的變化在1998 年被 Bi 和Poo 提出，並且稱為 Spike-timing-dependent synaptic plasticity (STDP)。[4-7,12]
故根據近年來的實驗結果，Bi 和 Poo 在 1998 年，整理其研究結果，根據從 實驗得到的兩神經元產生神經脈衝的時間差對興奮性突觸後電流(EPSC)的影響 關係，提出了spike-timing-dependent plasticity (STDP) rule ，針對突觸可塑性變 化的時間不對稱方面，修正了Hebbian Learning Rule，以此來作為突觸可塑性變 化的模型。[12]
STDP synaptic modulation model
(35)
⎩⎨
⎧
<
Δ Δ
−
>
Δ Δ
= − Δ Δ
−
−
+ +
0 )
/ exp(
0 )
/ ) exp(
( A t t
t t
t A wij
τ τ 其中，我們採用
A+ = 0.005 ， A- = 0.00525 ， τ+ = τ- = 20 (msec) ，Δt = tpost - tpre [4]
![Fig. 22: Critical window for synaptic modulation。 LTP/LTD 的產生與突觸前神經 元跟突觸後神經元兩者的興奮時間是相關聯的。由圖中可以清楚的看出突觸後電 流相對於兩者興奮產生的時間差之變化情況，這種相關聯的變化在 1998 年被 Bi 和 Poo 提出，並且稱為 Spike-timing-dependent synaptic plasticity (STDP)。[4-7,12]](https://thumb-ap.123doks.com/thumbv2/9libinfo/7113886.34110/12.892.270.628.110.397/Fig中可楚的看出突觸後電流相兩者興奮產生時間差之變化情況這種.webp)
