神經(jīng)核團異常放電的數(shù)值模擬與FPGA實現(xiàn)

王超，逯邁

蘭州交通大學(xué)光電技術(shù)與智能控制教育部重點實驗室，甘肅蘭州730070

前言

功能相似的神經(jīng)細胞聚集成界限明顯的群落，稱為神經(jīng)核團或神經(jīng)節(jié)?；咨窠?jīng)節(jié)，簡稱基底節(jié)，是大腦皮層下一系列神經(jīng)核團的總稱，是神經(jīng)系統(tǒng)中調(diào)節(jié)機體運動的重要結(jié)構(gòu)，與隨意運動的產(chǎn)生、肌肉緊張和感覺信息的處理等直接相關(guān)，同時還被認為與腦的感知、認知、記憶、學(xué)習(xí)等機制有著密切的聯(lián)系［1］。基底節(jié)主要包括紋狀體、蒼白球外側(cè)（the External Segment of the Globus Pallidus, GPe）、蒼白球內(nèi)側(cè)（the Internal Segment of the Globus Pallidus,GPi）、黑質(zhì)（the Substantia Nigra, SN）、丘腦底核（Subthalamic Nucleus, STN）等核團?；咨窠?jīng)節(jié)的病變會導(dǎo)致帕金森病、亨廷頓舞蹈病、癲癇等運動紊亂性疾病，是腦神經(jīng)科學(xué)研究的重點［2］。

建立基底節(jié)的數(shù)學(xué)模型并利用計算機軟件進行數(shù)值模擬有助于模擬基底節(jié)中各神經(jīng)核團的放電特性及其運動控制功能等［3］。1989年，Albin 等［4］提出宏觀的基底節(jié)直接-間接通路模型，研究帕金森病的發(fā)病機制；2001年，Gurney等［5］提出基底節(jié)的選擇-控制模型，模擬基底節(jié)的選擇-控制功能；2004年，Rubin等［6］提出一個精細的基底節(jié)神經(jīng)核團網(wǎng)絡(luò)結(jié)構(gòu)，并研究深部腦刺激治療帕金森病的作用機理；2001年，Humphries 等［7］利用Integrate-Fire（IF）神經(jīng)元模型搭建基底節(jié)核團網(wǎng)絡(luò)，模擬基底節(jié)的動作選擇功能和振蕩放電；2013年，Thibeault等［8］利用Izhikevich神經(jīng)元模型搭建基底節(jié)核團網(wǎng)絡(luò)并模擬帕金森病不同放電狀態(tài)。

近年來利用硬件模擬神經(jīng)元及神經(jīng)網(wǎng)絡(luò)的放電特性也得到越來越多的關(guān)注。2004年，Grass 等［9］利用Xilinx 公司的FPGA 芯片實現(xiàn)了Hodgkin-Huxley（HH）神經(jīng)元模型；2008年，Cassidy 等［10］提出利用Izhikevich神經(jīng)元模型硬件實現(xiàn)大規(guī)模神經(jīng)網(wǎng)絡(luò)的方法；2013年，Paolo 等［11］將Izhikevich 模型搭建的神經(jīng)網(wǎng)絡(luò)應(yīng)用于腦機接口；2015年，王金龍等［12］基于FPGA 芯片實現(xiàn)HH 模型，為硬件模擬神經(jīng)元提供了新的參考。但是大部分研究集中于神經(jīng)元及神經(jīng)元網(wǎng)絡(luò)的硬件實現(xiàn)，對于神經(jīng)核團網(wǎng)絡(luò)的硬件實現(xiàn)只有很少一部分。2006年，Prescott 等［13］設(shè)計制造基于基底節(jié)選擇-控制模型的機器人，模仿基底節(jié)的選擇與控制功能；2015年，Yang等［14］利用簡化的突觸模型在FPGA芯片實現(xiàn)一種基底節(jié)核團網(wǎng)絡(luò)模型。

模擬神經(jīng)核團的不同放電狀態(tài)可以深入研究相關(guān)神經(jīng)性疾病的發(fā)病機制，探究治療方式的作用機理，如深部腦刺激、經(jīng)顱磁刺激等。采用價格低廉、運行速度快的FPGA 芯片實現(xiàn)的神經(jīng)核團網(wǎng)絡(luò)具有重要的應(yīng)用價值，如替代動物活體實驗，降低醫(yī)學(xué)實驗成本；應(yīng)用于醫(yī)療康復(fù)工程，替代受損的神經(jīng)核團［15］；腦機接口技術(shù)，輔助神經(jīng)性疾病患者的生活，為研究和治療神經(jīng)性疾病提供新的思路［16］；嵌入智能機器人，模擬人類簡單的選擇運動功能，實現(xiàn)類腦式的機器人操控，為智能機器人的設(shè)計提供新的方向。

本研究依據(jù)Rubin 等［6］提出的基底節(jié)神經(jīng)核團網(wǎng)絡(luò)結(jié)構(gòu)，以Izhikevich 神經(jīng)元模型為神經(jīng)核團網(wǎng)絡(luò)的節(jié)點，選用化學(xué)突觸連接各神經(jīng)元，運用Simulink軟件和DSP Builder軟件分別對神經(jīng)核團網(wǎng)絡(luò)進行建模。對比分析兩種軟件仿真結(jié)果的相對標準誤差，驗證建模的一致性。通過DSP Builder 軟件與Quartus II 軟件聯(lián)合編譯的方法，對神經(jīng)核團網(wǎng)絡(luò)進行FPGA 硬件實現(xiàn)，再現(xiàn)基底節(jié)神經(jīng)核團正常放電和異常放電兩種狀態(tài)。

1 基底節(jié)神經(jīng)核團網(wǎng)絡(luò)模型

2004年，Rubin 等［6］提出一種基底節(jié)神經(jīng)核團網(wǎng)絡(luò)模型，研究各神經(jīng)核團的放電波形和深部腦刺激作用機理，其提出的網(wǎng)絡(luò)模型包含GPe、STN、GPi、丘腦皮層（Thalamic Cortex,TC）這4 種神經(jīng)核團。如圖1所示，紅色箭頭表示化學(xué)突觸興奮性輸入，黑色箭頭表示化學(xué)突觸抑制性輸入，綠色箭頭表示輸入到TC 神經(jīng)元的皮層運動信號，3 個STN 神經(jīng)元興奮性輸入到1 個GPe 神經(jīng)元，2 個GPe 神經(jīng)元抑制性輸入到1 個STN 神經(jīng)元，1 個GPe 神經(jīng)元抑制性輸入到1個GPi 和1 個STN 神經(jīng)元興奮性輸入到同一個GPi神經(jīng)元，最后由GPi 神經(jīng)元匯總輸入到TC 神經(jīng)元［6］。在異常放電時，來自GPi節(jié)律性的抑制性突觸輸入減弱了TC 中繼神經(jīng)元對來自皮層的運動信號的響應(yīng)能力，導(dǎo)致帕金森癥的運動障礙。

圖1 基底節(jié)神經(jīng)核團網(wǎng)絡(luò)Fig.1 Network of basal ganglia nuclei

2 Simulink建模仿真

2.1 神經(jīng)元模型

搭建基底節(jié)神經(jīng)核團網(wǎng)絡(luò)首先要選擇一個合適的神經(jīng)元模型。HH 模型過于復(fù)雜，消耗硬件資源多；IF 模型放電模式單一，不能表現(xiàn)神經(jīng)核團的多種放電模式［17-18］。本研究選用的Izhikevich 神經(jīng)元模型可以表現(xiàn)20 種放電模式，其中包括神經(jīng)核團的簇放電和峰放電［19］。Izhikevich神經(jīng)元數(shù)學(xué)模型如下：

其中，V是神經(jīng)元膜電位，單位為mV；u是膜電位恢復(fù)變量，表示K+通道的激活和Na+通道的失活，膜電位恢復(fù)變量u為膜電位V提供負反饋；a、b、c、d為無量綱參數(shù)；Iapp為外部輸入電流；Isyn為輸入該神經(jīng)元所有突觸電流的總和。

2.2 化學(xué)突觸模型

化學(xué)突觸采用文獻［20］中的模型，該模型可以實現(xiàn)神經(jīng)核團網(wǎng)絡(luò)中化學(xué)突觸的抑制性輸入和興奮性輸入，其表達式為：

其中，Isyn表示突觸后電流；Gsyn表示耦合強度，單位為mS/cm2；S表示突觸后膜上化學(xué)門控離子通道開放比例；ES是突觸可逆電位，單位為mV；F(Vpre)表示突觸前神經(jīng)元膜電位的S型函數(shù)；α和β分別表示離子通道打開和關(guān)閉的速率；Vpre是突觸前神經(jīng)元的膜電位，單位為mV；θsyn是突觸前可逆電位，單位為mV［21］。當ES=0 mV時，表示化學(xué)突觸興奮性輸入；當ES=-75 mV，表示化學(xué)突觸抑制性輸入。表1為各神經(jīng)元間突觸系數(shù)Gsyn的值［14］。

表1 突觸系數(shù)Gsyn（mS/cm2）Tab.1 Synaptic coefficient Gsyn(mS/cm2)

2.3 皮層運動信號

根據(jù)Rubin 等［6］和Yang 等［14］的研究，輸入到TC核團的皮層運動信號ISM用式（3）模擬：

其中，H為階躍函數(shù)；皮層電流最大值imax=30 pA；脈沖間隔時間參數(shù)分別為ρ=25 ms，δ=3 ms。波形如圖2所示。

圖2 皮層運動信號ISMFig.2 Cortical motor signal ISM

根據(jù)Thibeault 等［8］和Yang 等［14］的研究，對a、b、c、d參數(shù)進行不同設(shè)置，可以表現(xiàn)出GPe、STN、GPi、TC神經(jīng)元不同放電特性。具體參數(shù)設(shè)置見表2和3。

表2 神經(jīng)元正常放電參數(shù)Tab.2 Parameters for neurons in normal discharge

表3 神經(jīng)元異常放電參數(shù)Tab.3 Parameters for neurons in abnormal discharge

2.4 Simulink仿真結(jié)果

圖3為Simulink 軟件仿真結(jié)果，分別為GPe、STN、GPi 和TC 神經(jīng)元正常放電和異常放電的數(shù)值模擬波形。可以觀察到GPe 和GPi 神經(jīng)元異常放電時，放電波形從峰放電變?yōu)榇胤烹?。當STN 神經(jīng)元異常放電時，放電頻率相對于正常放電更快。當TC神經(jīng)元異常放電時，不能對運動信號做出正常響應(yīng)，出現(xiàn)信號的缺失。

圖3 GPe、STN、GPi和TC神經(jīng)元正常放電和異常放電的Simulink模擬結(jié)果Fig.3 Simulink simulation results of GPe,STN,GPi and TC neurons in normal and abnormal discharges

3 DSP Builder建模仿真

DSP Builder 軟件是Altera 公司提供的一種系統(tǒng)級開發(fā)工具，將建模仿真和硬件實現(xiàn)結(jié)合起來，在Matlab/Simulink 環(huán)境下以庫文件的形式進行圖形化的設(shè)計和仿真，并將Quartus II 軟件作為底層設(shè)計工具置于后臺。將建好的DSP Builder 設(shè)計文件（.mdl）通過Signal Compiler 轉(zhuǎn)化成硬件描述語言HDL 文件（.vhd），進而可以在Quartus II 環(huán)境下進行編譯和FPGA硬件配置。

3.1 神經(jīng)元及突觸建模

在硬件實現(xiàn)時，Izhikevich 原模型資源消耗相對過大，不利于神經(jīng)核團網(wǎng)絡(luò)的實現(xiàn)。2012年，Soleimani 等［22］利用四階線性近似的方法簡化了Izhikevich模型，模型如下：

根據(jù)Soleimani 等［22］提出的選擇算法，設(shè)定m1=0.75，m2=0.375，m3=11，m4=62.5。

對神經(jīng)元和化學(xué)突觸進行建模首先要選擇合適的離散方法，考慮硬件資源的消耗，采用歐拉法對神經(jīng)元模型和化學(xué)突觸進行離散化處理。用歐拉法對Izhikevich模型進行處理，得到差分方程：

圖4為神經(jīng)元的建模過程，采用移位器和加法器替代乘法器，減少乘法器的使用，節(jié)約硬件資源［23］。移位器可以根據(jù)命令向左或向右移動數(shù)據(jù)位，移位操作為算術(shù)移位。用移位器和加法器的運算值趨近于模型中m1、m2、a、b和△T的值，以此替代乘法運算，如單位時間△T=0.02≈1/64+1/256。

圖4 Izhikevich模型建模流程Fig.4 Establishment of Izhikevich model

運用DSP Builder軟件對神經(jīng)元模型和突觸模型分別建模，然后用化學(xué)突觸將各個神經(jīng)元按照神經(jīng)核團網(wǎng)絡(luò)結(jié)構(gòu)依次連接，構(gòu)成基底節(jié)神經(jīng)核團網(wǎng)絡(luò)。在DSP Builder 建模時需要考慮時序，可以更準確地將突觸的延時表現(xiàn)出來，取神經(jīng)元模型計算1次的時間作為突觸傳導(dǎo)的延時。

3.2 DSP Builder軟件仿真結(jié)果

考慮硬件資源的限制，本研究實現(xiàn)的神經(jīng)核團網(wǎng)絡(luò)由3 個GPe 神經(jīng)元、3 個STN 神經(jīng)元、3 個GPi 神經(jīng)元、1 個TC 神經(jīng)元和30 個化學(xué)突觸構(gòu)成，最后由GPi 神經(jīng)元匯總輸入到1 個TC 神經(jīng)元。圖5為DSP Builder軟件的仿真結(jié)果，分別為GPe、STN、GPi和TC神經(jīng)元正常放電和異常放電的仿真結(jié)果，與Simulink仿真結(jié)果基本一致。

圖5 GPe、STN、GPi和TC神經(jīng)元正常放電和異常放電的DSP Builder模擬結(jié)果Fig.5 DSP Builder simulation results of GPe,STN,GPi and TC neurons in normal and abnormal discharges

3.3 波形誤差分析

為了量化評估Izhikevich線性簡化模型搭建的神經(jīng)核團網(wǎng)絡(luò)的精確性，采用標準誤差公式進行計算分析［22］。選取一個周期內(nèi)兩種軟件的模擬結(jié)果，即0～150 ms 樣本值進行對比分析。因為Simulink 軟件沒有時序控制模塊，為了更好地對比兩款軟件的仿真結(jié)果，DSP Builder 軟件仿真中不考慮突觸延時，并對仿真時間進行調(diào)整。標準誤差（Root Mean Square Error，RMSE）計算見式（6）：

其中，Vorigin為原模型電位，VPWL為線性簡化模型電位。計算結(jié)果如表4所示，可以發(fā)現(xiàn)各神經(jīng)元在不同狀態(tài)的標準誤差均小于0.1，其中TC 神經(jīng)元的誤差最大。造成誤差的主要原因是線性簡化模型中利用移位器代替了乘法器，a、b、c、d等值有略微改變，此外硬件實現(xiàn)中考慮資源的消耗，對輸出波形的位數(shù)進行了限制。

表4 各神經(jīng)核團的仿真結(jié)果標準誤差Tab.4 Root mean squared error of simulation results of nuclei each neuron

4 神經(jīng)核團的硬件實現(xiàn)

本研究所用FPGA 開發(fā)板為Altera 公司CycloneIV 系列EP4CE115F23I7，數(shù)模轉(zhuǎn)換模塊為DAC900E，采用架構(gòu)在MATLAB/Simulink 上的DSP Builder 輔助設(shè)計工具箱對神經(jīng)元、化學(xué)突觸和神經(jīng)核團網(wǎng)絡(luò)進行建模，VerilogHDL 語言描述分頻文件。所用各軟件版本為Quartus II 13.0、DSP Builder 13.0以及MATLAB 2013。

將DSP Builder軟件中建模完成的神經(jīng)核團網(wǎng)絡(luò)進行編譯，然后導(dǎo)入Quartus II 軟件加入分頻文件并配置管腳，運行成功后導(dǎo)入FPGA 芯片，通過DA 芯片輸出到示波器，觀察到各神經(jīng)元放電波形。圖6為GPe、STN、GPi 和TC 神經(jīng)元正常放電和異常放電時的硬件實現(xiàn)波形。圖中的毛刺是外部干擾的原因，并不影響結(jié)果的正確性。圖7為實現(xiàn)的神經(jīng)核團網(wǎng)絡(luò)實物圖。

圖6 GPe、STN、GPi和TC神經(jīng)元正常放電和異常放電的硬件結(jié)果Fig.6 Hardware implementation results of GPe,STN,GPi and TC neurons in normal and abnormal discharges

圖7 神經(jīng)核團網(wǎng)絡(luò)硬件模型Fig.7 Hardware model of the network of nuclei

5 討論與分析

帕金森病是中老年人常見的一種神經(jīng)退行性疾病。目前帕金森病的致病原因尚不清楚，但是研究表明與環(huán)境、遺傳和衰老等多種因素共同作用有關(guān)。典型的臨床癥狀有肌肉僵直，靜止性震顫和運動遲緩等。其中，TC 神經(jīng)元中繼能力的喪失是帕金森病的主要癥狀之一［24］。Rubin 等［6］和Yang 等［14］在研究中提出一種丘腦神經(jīng)元對皮層運動信號的中繼可靠性指標（Reliability Index, RI），是評價帕金森狀態(tài)的一個重要指標。

其中，Nerrors表示丘腦元沒有準確響應(yīng)皮層運動信號脈沖的個數(shù)，NSM表示皮層運動信號脈沖的總個數(shù)。當RI＞0.9 時，表示正常狀態(tài)；RI＜0.6 時，表示帕金森狀態(tài)；其余為中間狀態(tài)。

圖8、圖9為實現(xiàn)的神經(jīng)核團網(wǎng)絡(luò)在異常放電和正常放電時TC 神經(jīng)元對皮層運動信號的響應(yīng)波形。上欄波形為皮層運動信號，下欄波形為TC 對運動信號的響應(yīng)。

圖8 丘腦皮層神經(jīng)元正常放電Fig.8 TC neuron in normal discharge

選取其中20 個運動信號脈沖個數(shù)進行計算（如圖9紅框）。正常放電時，RI=1.0＞0.9，TC 神經(jīng)元沒有出現(xiàn)信號缺失；異常放電時，RI=0.3＜0.6，TC 神經(jīng)元出現(xiàn)信號的缺失，證明實現(xiàn)的基底節(jié)神經(jīng)核團表現(xiàn)為帕金森狀態(tài)。

圖9 丘腦皮層神經(jīng)元異常放電Fig.9 TC neuron in abnormal discharge

RI 可以驗證硬件實現(xiàn)的基底節(jié)神經(jīng)核團網(wǎng)絡(luò)可以展現(xiàn)帕金森病的一種放電狀態(tài)。但是由于RI過于單一，僅限于理論上的驗證，還要結(jié)合各神經(jīng)核團的放電波形等其他因素進一步驗證。

帕金森病與蒼白球核團和STN 核團的異常放電有著明顯聯(lián)系。在帕金森病相關(guān)臨床研究和動物模型生理實驗中發(fā)現(xiàn)，GPe 和GPi神經(jīng)元出現(xiàn)明顯的爆發(fā)式簇放電［25］，并且GPi低頻爆發(fā)式簇放電與帕金森癥的震顫有關(guān)［26］。從模擬的神經(jīng)核團放電波形可以觀察到，在異常狀態(tài)時，GPe 和GPi 神經(jīng)元再現(xiàn)了帕金森病理模型的簇放電。同時可以觀察到STN 神經(jīng)元的放電波形出現(xiàn)異常，放電頻率加快。這也符合病理模型中對STN 核團的研究，STN 神經(jīng)元放電活動的改變引發(fā)明顯的行為異常?；坠?jié)SN 致密部多巴胺能分泌的減少影響GPe 核團對STN 核團的抑制，使STN 神經(jīng)元連續(xù)非正常放電，并且放電頻率增加，導(dǎo)致輸出核團GPi 的興奮，過度抑制皮層的相關(guān)運動區(qū)域，進而導(dǎo)致人體產(chǎn)生運動遲緩，動作保持困難和震顫等運動功能障礙［27］。

現(xiàn)在大多數(shù)研究集中于實現(xiàn)單個神經(jīng)元和鏈式神經(jīng)網(wǎng)絡(luò)等簡單的神經(jīng)網(wǎng)絡(luò)或者稍微復(fù)雜的小世界神經(jīng)網(wǎng)絡(luò)［28］，但是對于功能性的神經(jīng)核團網(wǎng)絡(luò)研究較少，并且更多的研究限于數(shù)值模擬仿真。如Thibeault 等［8］通過數(shù)值模擬研究基底節(jié)神經(jīng)核團放電特性并研究腦深部電刺激對神經(jīng)核團網(wǎng)絡(luò)的影響，但是沒有涉及硬件實現(xiàn)。Humphries 等［7］、Prescott等［13］選用IF 模型或者人工神經(jīng)元等相對簡單的神經(jīng)元模型硬件實現(xiàn)神經(jīng)核團網(wǎng)絡(luò)，但是由于神經(jīng)元模型簡單，不能表現(xiàn)神經(jīng)核團多種放電模式。本研究則選用動力特性較為復(fù)雜的Izhikevich模型搭建神經(jīng)核團網(wǎng)絡(luò)，不僅通過軟件模擬神經(jīng)核團的兩種放電狀態(tài)，而且結(jié)合FPGA進行了硬件實現(xiàn)。

6 結(jié)論

多種神經(jīng)性疾病的產(chǎn)生與神經(jīng)核團的異常放電有著重要聯(lián)系，數(shù)值模擬和FPGA 硬件模擬神經(jīng)核團異常放電為研究神經(jīng)核團放電狀態(tài)與神經(jīng)性疾病間的聯(lián)系提供了更多的研究方式。本研究在Simulink、DSP Builder 和Quartus II 聯(lián)合仿真編譯的基礎(chǔ)上，通過單片大資源FPGA 模擬了基底節(jié)神經(jīng)核團網(wǎng)絡(luò)異常放電和正常放電兩種放電模式。采用相對標準誤差公式量化分析了線性簡化Izhikevich 模型的誤差，驗證了兩種軟件建模的一致性。在神經(jīng)核團網(wǎng)絡(luò)正常狀態(tài)時，RI=1.0，表明TC神經(jīng)元如實地中轉(zhuǎn)了傳入的每一個運動信號；當核團網(wǎng)絡(luò)處于異常狀態(tài)時，RI=0.3，小于文獻［14］所實現(xiàn)的神經(jīng)核團帕金森狀態(tài)的RI 值（0.55），表明TC 神經(jīng)元對運動信號的中轉(zhuǎn)出現(xiàn)明顯的缺失。此外，結(jié)合相關(guān)帕金森病理實驗記錄進行分析，證明實現(xiàn)的神經(jīng)核團網(wǎng)絡(luò)可以模擬帕金森病的一種放電狀態(tài)。

由于FPGA 芯片資源的限制，本研究只實現(xiàn)了一個小規(guī)模的神經(jīng)核團網(wǎng)絡(luò)，并且實現(xiàn)的神經(jīng)核團放電模式較少，僅有峰放電和簇放電。后續(xù)的研究工作可以在這個基礎(chǔ)上完成大規(guī)模、多種放電功能的神經(jīng)核團網(wǎng)絡(luò)的硬件實現(xiàn)。用移位器和加法器代替乘法器，可以大大節(jié)約FPGA 資源的消耗，為以后大規(guī)模神經(jīng)元網(wǎng)絡(luò)的實現(xiàn)提供新的參考。此外，硬件實現(xiàn)的基底神經(jīng)核團可以應(yīng)用于替代動物活體實驗、腦機接口和醫(yī)療康復(fù)工程的研究中。