陳嘉卿，黃嘉嘉，許弢

（五邑大學(xué) 智能制造學(xué)部，廣東 江門 529020）

小腦模型是由Albus[1]于1975 年仿照小腦皮層結(jié)構(gòu)重建的關(guān)聯(lián)控制器經(jīng)多年發(fā)展形成的神經(jīng)網(wǎng)絡(luò)模型，該模型能對生物小腦學(xué)習(xí)行為進(jìn)行仿真，但其仿真結(jié)果與實(shí)際細(xì)胞電位變化還存在差異，其預(yù)測功能和記憶形成的機(jī)制尚未明確.大規(guī)模細(xì)胞的電生理測量難度非常大，生物實(shí)驗(yàn)需要花費(fèi)的時(shí)間和材料成本非常高，構(gòu)建大規(guī)模的仿生小腦模型可以解決這些難題.目前比較流行的小腦模型基礎(chǔ)架構(gòu)有Marr-Albus-Ito[1]的小腦模型、Hodgkin-Huxley-model 模型[2]等.基于這些基礎(chǔ)模型，科學(xué)家能很好地結(jié)合模型的特點(diǎn)去解決工程上的一些難題，比如對機(jī)器人控制的改進(jìn)[3-4]和利用小腦模型的生物行為特性進(jìn)行控制器開發(fā)[5-7]等，可見模型在工程上有很大的應(yīng)用價(jià)值.同時(shí)，科學(xué)家通過實(shí)時(shí)仿生小腦模型可以進(jìn)一步探討學(xué)習(xí)記憶在小腦中的形成與保存，對研究小腦的學(xué)習(xí)機(jī)制[8]有很大幫助.實(shí)時(shí)的模型一般以現(xiàn)場可編程門陣列（Field Programmable Gate Array，F(xiàn)PGA）[9]和圖形處理器（Graphics Processing Unit，GPU）[10]為基礎(chǔ)，值得注意的是隨著統(tǒng)一計(jì)算設(shè)備架構(gòu)（Compute Unified Device Architecture，CUDA）[11]技術(shù)的進(jìn)一步優(yōu)化，并行加速開始進(jìn)入了研究人員的視野，比如針對神經(jīng)網(wǎng)絡(luò)的并行計(jì)算改進(jìn)[12-13]，但GPU 的實(shí)時(shí)模型還有優(yōu)化的空間，在資源利用和時(shí)序上還能進(jìn)一步優(yōu)化.因此本文在Yamazaki 模型[14]的基礎(chǔ)上通過引入流的概念[15]優(yōu)化GPU 并行線程，為不同計(jì)算模塊和不同細(xì)胞分配相應(yīng)的線程并使用新的CUDA 特性（如區(qū)塊均衡和驅(qū)動內(nèi)存對齊等技術(shù)）加速網(wǎng)絡(luò)的計(jì)算，最后在不同細(xì)胞數(shù)量規(guī)模下觀察小腦模型的實(shí)時(shí)運(yùn)行規(guī)模和速度.

1 材料與方法

1.1 模型結(jié)構(gòu)

小腦皮層接受到刺激，進(jìn)行信息處理的信息流如圖1 的網(wǎng)絡(luò)結(jié)構(gòu)所示.以約10 萬個(gè)細(xì)胞規(guī)模的模型作為例子，模型由6 種細(xì)胞組成，分別為顆粒細(xì)胞（granule cells，GCs）、高爾基細(xì)胞（golgi cells，GOs）、藍(lán)狀細(xì)胞（basket cell，BS）、浦肯野細(xì)胞（purkinje cells，PC）、下橄欖核（inferior olive，IO）和小腦深部核團(tuán)（deep cerebellar nuclei，DCN）.還有3 種突觸纖維，分別為苔蘚纖維（mossy fibers，MF）、平行纖維（parallel fibers，PF）和攀延纖維（climbing fiber，CF）.整體信息流從MF流入顆粒細(xì)胞層和小腦核，信息通過顆粒層分別傳輸?shù)街虚g層和浦肯野細(xì)胞，最后浦肯野細(xì)胞抑制DCN，DCN 作為最后的輸出.指導(dǎo)信號由IO 發(fā)出，通過CF 傳遞給PC.

圖1 小腦模型網(wǎng)絡(luò)架構(gòu)圖

模型的規(guī)模包括102 400 個(gè)顆粒細(xì)胞，1 024 個(gè)高爾基細(xì)胞，16 個(gè)浦肯野細(xì)胞，16 個(gè)中間層細(xì)胞和1 個(gè)小腦核細(xì)胞.程序每次試驗(yàn)的仿真時(shí)長為1 666 ms，其他細(xì)節(jié)參考試驗(yàn)方案章節(jié).苔蘚纖維的突觸數(shù)目為100 條，平行纖維的突觸數(shù)目由顆粒細(xì)胞簇決定.

信號經(jīng)過顆粒層后由PF 輸入到PC 中，另外一路信號通過中間層的星狀細(xì)胞和BS 處理后發(fā)出抑制信號給PC.PC 同時(shí)受到了顆粒層的刺激和中間層的抑制.CF 的信號只作為突觸權(quán)重調(diào)節(jié)的誘導(dǎo)，誘導(dǎo)LTP（long-term potentiation）和LTD（long-term depression）[17]的產(chǎn)生.

模型的輸出是DCN，它同時(shí)受到浦肯野細(xì)胞抑制輸入和苔蘚纖維的刺激輸入.類似于輸出的時(shí)候加入一個(gè)正則項(xiàng)，最后通過調(diào)節(jié)平行纖維的權(quán)重達(dá)到學(xué)習(xí)的效果.

圖1 中箭頭表示刺激信號傳遞方向，圓頭表示抑制信號的傳遞方向.顆粒層中的圓圈表示顆粒細(xì)胞，方框?yàn)楦郀柣?xì)胞.纖維上的數(shù)字代表其規(guī)模.

1.2 計(jì)算機(jī)制

每個(gè)細(xì)胞，通過建立積分電導(dǎo)模型[18]式（1）以及概率電導(dǎo)模型式（2）來計(jì)算模型中的單個(gè)神經(jīng)元以及一簇神經(jīng)元隨時(shí)間變化的活性.最后神經(jīng)元的膜電位變化按式（3）計(jì)算.

其中，g是電導(dǎo)，是最大電導(dǎo)，用下標(biāo)c表示不同離子通道；wj是突觸權(quán)重；α(t-s)是電勢衰減函數(shù)；δj(s)是脈沖發(fā)射情況，用1 和0 分別代表有或沒有觸發(fā)動作電位.

其中，gc代表相關(guān)離子通道的導(dǎo)通概率，m是化學(xué)能推動離子通道打開的概率，n次方代表有多少個(gè)這樣的離子通道.

其中，V(t) 是t時(shí)刻的神經(jīng)元膜電位；Ec是觸發(fā)電勢；C是膜電容.該公式是每個(gè)細(xì)胞總電勢的積分公式.這里以浦肯野細(xì)胞作為例子，刺激信號的正向模型中預(yù)測計(jì)算為：

其中，PC(t)是t時(shí)刻MF 輸入經(jīng)過PF 后的結(jié)果，為第i個(gè)苔蘚纖維鏈接的突觸權(quán)重，MFi(t)在t時(shí)刻的時(shí)候刺激信號.

按照主題內(nèi)容，可將本次統(tǒng)計(jì)的189篇論文劃分為9個(gè)類別（詳見表1）。其中，研究范疇主要圍繞MOOC環(huán)境下圖書館及館員的角色定位、MOOC環(huán)境下圖書館服務(wù)創(chuàng)新展開，研究切入點(diǎn)主要包括信息素養(yǎng)教育、版權(quán)服務(wù)、MOOC平臺技術(shù)支持與優(yōu)化，研究基礎(chǔ)側(cè)重于對國外圖書館MOOC實(shí)踐經(jīng)驗(yàn)的借鑒。

當(dāng)IO 被激發(fā)的時(shí)候PC 也同時(shí)被激發(fā)，計(jì)算在固定窗口下GR 的觸發(fā)數(shù)量作為減少權(quán)重的指數(shù)，公式如下：

最終結(jié)果是由苔蘚纖維傳輸?shù)拇碳ば盘柡推挚弦凹?xì)胞的抑制共同作用輸出的：

其中，第二項(xiàng)主要為負(fù)數(shù)，對CN 起抑制作用.

最后根據(jù)CN 是否激活來判斷輸出：

其中，OUT 作為最終輸出，可以根據(jù)不同任務(wù)改變輸出的形式，比如本文采用激活率作為標(biāo)準(zhǔn)，在一段時(shí)間內(nèi)激活率達(dá)到一定程度作為模型觸發(fā)的輸出標(biāo)準(zhǔn).

1.3 硬件實(shí)現(xiàn)

本文代碼使用C++語言，并使用CUDA 框架提供的API 在內(nèi)存和顯存之間傳輸數(shù)據(jù).先構(gòu)建數(shù)據(jù)結(jié)構(gòu)，每個(gè)種類的細(xì)胞包含細(xì)胞數(shù)目、電導(dǎo)和觸發(fā)電勢等各項(xiàng)參數(shù).自動平衡CUDA 分配的區(qū)塊和塊內(nèi)線程數(shù)，盡可能利用顯卡的所有資源，比傳統(tǒng)的固定區(qū)塊數(shù)計(jì)算要更加有效.其次，在多流處理上采用原子操作（類似數(shù)據(jù)的加鎖和解鎖），防止在異步操作時(shí)候取舊數(shù)據(jù)或者覆蓋了不應(yīng)該覆蓋的數(shù)據(jù).同時(shí)，優(yōu)化對比數(shù)據(jù)大小的機(jī)制，采用類似減法器的對比方案.

如圖2 所示，不同的工作流處理指定的運(yùn)算，因?yàn)檩斎胄盘柺歉怕誓Ｐ?，所以由梅森旋轉(zhuǎn)法[19]在顯存中生成隨機(jī)數(shù)并作為刺激信號的發(fā)生器.Aph 是每個(gè)細(xì)胞超極化后電位恢復(fù)的處理過程，超極化指前一個(gè)時(shí)刻被激活后，則當(dāng)前時(shí)刻膜電位就變成衰減的膜電位.英文縮寫（如GR）表示該細(xì)胞（顆粒細(xì)胞）的內(nèi)部膜電位變化處理，一般是針對刺激信號進(jìn)行膜電位積分.而帶箭頭的英文縮寫（如MF→GR）表示信號從輸出纖維（苔蘚纖維）到目標(biāo)細(xì)胞（顆粒細(xì)胞）的過程.圖中每個(gè)塊的高度代表運(yùn)算規(guī)模，最長的運(yùn)算模塊統(tǒng)一為顆粒細(xì)胞數(shù)目與2 的次方對齊，在10 萬規(guī)模細(xì)胞下，單個(gè)方框的規(guī)模是102 400.一共分為100 塊，每個(gè)塊內(nèi)有1 024 個(gè)線程數(shù).每次任務(wù)的塊大小相同能增加運(yùn)算效率.

圖2 GPU 計(jì)算流架構(gòu)圖

2 實(shí)驗(yàn)結(jié)果與分析

為了測試模型的性能，本文使用8 Gb 顯存的NVIDIA 顯卡，單精度浮點(diǎn)性能不低于8.1 TFlops.通過梅森旋轉(zhuǎn)算法生成偽隨機(jī)序列模擬細(xì)胞在小腦中的電生理活動.最后通過運(yùn)行不同規(guī)模的小腦網(wǎng)絡(luò)對比不同方案下的耗時(shí).

2.1 實(shí)驗(yàn)方案

本文進(jìn)行了小鼠的延遲眨眼條件反射dEBC（delay eyeblink conditioning）[20]仿真試驗(yàn)，觀察小腦對于延遲眨眼條件反射記憶的形成.刺激信號如圖3 所示，總仿真時(shí)長1 666 ms.200 ms 之前作為穩(wěn)定時(shí)間，從205 ms 開始發(fā)送有條件刺激，在450 ms 的時(shí)候發(fā)送無條件刺激并持續(xù)30 ms，最后有條件刺激信號在1 200 ms 時(shí)結(jié)束.這樣一個(gè)流程作為一次試驗(yàn)，重復(fù)100 次以達(dá)到現(xiàn)實(shí)小鼠學(xué)習(xí)到眨眼條件反射的最低次數(shù)，觀察過程中各個(gè)細(xì)胞的電生理變化.

圖3 延遲眨眼刺激信號

圖3 中藍(lán)色部分為條件刺激信號，以均勻分布的概率形式刺激顆粒細(xì)胞神經(jīng)元的觸發(fā)，黃色條狀作為無條件刺激信號輸入到下橄欖體（IO）.而圖中上部作為普通訓(xùn)練時(shí)候的刺激信號，下部作為測試時(shí)候的刺激信號.因?yàn)閷W(xué)會延遲眨眼后小鼠會在有條件刺激下產(chǎn)生無條件刺激的反應(yīng).

2.2 實(shí)驗(yàn)分析

模型仿真延遲眨眼條件反射試驗(yàn)結(jié)果如下：

如圖4 所示，上部分為不同訓(xùn)練次數(shù)下各個(gè)細(xì)胞的神經(jīng)元尖峰放電情況，下部分為學(xué)習(xí)過程中各個(gè)細(xì)胞在各個(gè)階段的電生理仿真圖像.觀察到PC 和CN 隨著訓(xùn)練次數(shù)變化明顯，CN 在IO 激活位置的激活率逐漸增加.其中trial 1 是網(wǎng)絡(luò)剛剛開始學(xué)習(xí)的時(shí)候，明顯發(fā)現(xiàn)這個(gè)時(shí)候的PC 的激活率很高，但隨著學(xué)習(xí)的進(jìn)行，激活率慢慢下降，到了trial 30 已經(jīng)有明顯的學(xué)習(xí)跡象.最后在trial 50的時(shí)候習(xí)得了延遲眨眼，判斷標(biāo)準(zhǔn)是CN 的激活率在IO 給出信號的位置激活率顯著上升，表明CN在這個(gè)時(shí)候的輸出.再隨著學(xué)習(xí)的進(jìn)行，在trial 100 的時(shí)候已經(jīng)會提前預(yù)判眨眼的動作，觀察CN最開始的激活位置比IO 的指導(dǎo)信號還提前了，到此該模型已經(jīng)學(xué)會了延遲眨眼條件反射試驗(yàn).

圖4 電生理仿真圖

觀察發(fā)現(xiàn)GR 和GO 在學(xué)習(xí)的過程中激活情況變化不大，而PC 的變化很大，可以看出學(xué)習(xí)的位點(diǎn)并不在前面信號預(yù)處理的細(xì)胞上而是在浦肯野細(xì)胞上.而浦肯野細(xì)胞在IO 的指導(dǎo)信號下，受到LTP 和LTD 的機(jī)制影響調(diào)節(jié)平行纖維（PF）的權(quán)重.

圖5 是顆粒細(xì)胞（GR）和高爾基細(xì)胞（GO）在一次仿真過程中的神經(jīng)元激活情況，圖中每個(gè)黑點(diǎn)代表該序列的神經(jīng)元細(xì)胞在該時(shí)間點(diǎn)激活，縱軸是每個(gè)不同的細(xì)胞，并按一定順序排序.可以看出不同簇的細(xì)胞對刺激信號產(chǎn)生不同的反應(yīng)，在收到刺激信號的時(shí)候高爾基細(xì)胞的平均激活率顯著提高，顆粒細(xì)胞則是表現(xiàn)出對刺激信號的不同特性，有的反而激活率下降.這樣有利于信號的模糊化和預(yù)處理特性，使信號的隱藏特征通過網(wǎng)絡(luò)的隨機(jī)連接關(guān)系展現(xiàn)出來.

圖5 細(xì)胞激活散點(diǎn)圖

3 討論

針對實(shí)時(shí)仿生小腦模型的加速和擴(kuò)大運(yùn)算規(guī)模的問題，本文提出了一種基于GPU 的流式處理優(yōu)化方法.通過試驗(yàn)，仿真小腦模型能學(xué)習(xí)延遲眨眼條件反射，并仿真其學(xué)習(xí)過程中神經(jīng)元的激活狀態(tài)變化，而學(xué)習(xí)記憶離不開突觸可塑性的形成，平行纖維是學(xué)習(xí)延遲眨眼條件反射的關(guān)鍵位點(diǎn).GR和GO 細(xì)胞則對數(shù)據(jù)起著清洗作用，挖掘輸入信號的潛在特征同時(shí)也給神經(jīng)網(wǎng)絡(luò)提供了更多的輸入模板.在不同細(xì)胞規(guī)模下的模型對比仿真表明，使用流方案的模型運(yùn)行速度更快，這個(gè)加速效果除了是流技術(shù)帶來的高效并行化[12]，還離不開優(yōu)化運(yùn)算單元分配[13].保證每次進(jìn)入顯卡的任務(wù)序列的運(yùn)算單元數(shù)目盡可能相同.

本模型相比于Yamazaki[15]的模型，使用流方案并行多模塊運(yùn)行，同時(shí)將所有數(shù)據(jù)都存放在顯存當(dāng)中，減少數(shù)據(jù)拷貝時(shí)間，加快了模型運(yùn)行速度.而在硬件方面對比FPGA 有著更大規(guī)模的神經(jīng)元，不受限于硬件上資源的緊缺.本文的關(guān)鍵點(diǎn)在于將不同運(yùn)算量的細(xì)胞分配到不同的流中，使并行程序中的有效并行數(shù)目得到充分的發(fā)揮.同時(shí)使用顯存存儲數(shù)據(jù)，減少了數(shù)據(jù)轉(zhuǎn)移帶來的性能損失，并且利用CUDA 的新特性加速計(jì)算單元.

圖6 是原方案和使用流方案的模型在不同規(guī)模神經(jīng)元下的仿真耗時(shí)對比，最大取40 萬規(guī)模進(jìn)行比較是因?yàn)楸灸Ｐ蛯?shí)時(shí)模擬1 666 ms 在40 萬規(guī)模左右接近極限.從圖中可以看出使用流方案的模型耗時(shí)更少.在小規(guī)模下兩者差距不明顯，但在大規(guī)模下兩者差距逐漸拉開，而且流方案下耗時(shí)方差更小更加穩(wěn)定，但隨著規(guī)模變大方差也逐漸變大.在40 萬規(guī)模下，原方案運(yùn)行時(shí)間已超過仿真總時(shí)長，不能滿足實(shí)時(shí)要求.

圖6 仿真耗時(shí)圖

4 結(jié)語

本論文主要利用 CUDA 新特性以及均衡各個(gè)線程的計(jì)算量來優(yōu)化模型，提高了模型的運(yùn)行速度，并且在相同規(guī)模下模型耗時(shí)更短.通過延遲眨眼條件反射實(shí)驗(yàn)的仿真，驗(yàn)證了模型的基本學(xué)習(xí)功能.本模型使用C++語言直接編寫CUDA 腳本.目前模型中的數(shù)值運(yùn)算采用矩陣運(yùn)算，但是有很多細(xì)胞不參與某些時(shí)刻的運(yùn)算，是否可以采用稀疏矩陣[21]或者優(yōu)化算法本身加速運(yùn)算是未來要研究的方向.