潘紅光，張 奇，米文毓，馬 彪

（1.西安科技大學(xué) 電氣與控制工程學(xué)院，陜西 西安710054；2.鄂爾多斯市神東工程設(shè)計(jì)有限公司，鄂爾多斯017000）

0 引 言

腦機(jī)接口（Brain-Machine Interface，BMI）技術(shù)近年來(lái)發(fā)展異常迅速［1-2］。BMI技術(shù)能夠通過(guò)采集人腦的電信號(hào)直接完成對(duì)輸出設(shè)備的操控，實(shí)現(xiàn)功能障礙者與外界的交流，從而提供了一種改善功能障礙者生活質(zhì)量的途徑［3］。其在醫(yī)療領(lǐng)域如殘疾人肢體功能恢復(fù)，生活領(lǐng)域如老年人生活輔助等方面應(yīng)用也越來(lái)越廣泛［4-5］。BMI系統(tǒng)主要組成部分為：解碼器、編碼器（本質(zhì)均為數(shù)學(xué)模型）和大腦皮層神經(jīng)元活動(dòng)的量測(cè)。其中，解碼器用于提取與任務(wù)相關(guān)的運(yùn)動(dòng)意圖；編碼器則將與運(yùn)動(dòng)相關(guān)的感知信息反饋到大腦。故BMI、大腦和假肢器官（如假臂）三者可構(gòu)成閉環(huán)的BMI系統(tǒng)，如圖1所示。

圖1 腦機(jī)接口系統(tǒng)Fig.1 System of brain-machine interface

目前針對(duì)BMI系統(tǒng)解碼器已有大量的研究。例如Goblu，孫京浩等人基于維納濾波和速度衰減卡爾曼濾波設(shè)計(jì)的解碼器，可以從連續(xù)的放電尖峰序列中解碼出速度、位置等信息［6-7］；Shanechi等人設(shè)計(jì)一個(gè)閉環(huán)的點(diǎn)過(guò)程過(guò)濾解碼器，可以研究反饋率和控制率對(duì)大腦控制連續(xù)BMI運(yùn)動(dòng)的影響（控制率表示從大腦向假體發(fā)送運(yùn)動(dòng)命令的頻率，反饋率表示向?qū)ο筇峁┘袤w的視覺(jué)反饋的頻率）［8］；潘家輝等人從多方面對(duì)P300成分進(jìn)行特征提取，并采用貝葉斯線性判別分析的方法進(jìn)行模式分類［9］。但是在上述解碼器中都未曾考慮腦電信號(hào)與時(shí)序的相關(guān)性。而腦電信號(hào)是一個(gè)動(dòng)態(tài)的時(shí)間序列并且其控制的外部設(shè)備（如假肢）的運(yùn)動(dòng)過(guò)程也是漸變的。而長(zhǎng)短期記憶網(wǎng)絡(luò)（Long Short-Term Memory，LSTM）適合處理與時(shí)間序列高度相關(guān)的問(wèn)題，并且在學(xué)習(xí)新信息的同時(shí)保持歷史信息持久化，因此基于LSTM網(wǎng)絡(luò)設(shè)計(jì)BMI系統(tǒng)解碼器，以期獲得更好的解碼效果。

近年來(lái)，針對(duì)BMI系統(tǒng)的研究已取得廣泛的進(jìn)展，例如Wodlinger等人同時(shí)從植入受試者左側(cè)運(yùn)動(dòng)皮層的2個(gè)96通道皮質(zhì)內(nèi)電極陣列中提取手形指令，實(shí)現(xiàn)被試者多自由度控制假肢［10］。Hong和Orsborn等人引入視覺(jué)反饋，證明閉環(huán)BMI系統(tǒng)較開(kāi)環(huán)BMI系統(tǒng)擁有更好的可靠性以及執(zhí)行效率［11-12］。Hortal和Abdullah等通過(guò)支持向量機(jī)和神經(jīng)網(wǎng)絡(luò)對(duì)腦電圖（electroencephalogram，EEG）信號(hào)進(jìn)行分類作為控制指令驅(qū)動(dòng)外部設(shè)備［13-14］。

基于Bullock等人提出的單關(guān)節(jié)信息傳輸（Single-Joint Information Transmission，SJIT）模型，設(shè)計(jì)基于LSTM網(wǎng)絡(luò)的解碼器來(lái)表征放電信號(hào)與肢體運(yùn)動(dòng)的非線性關(guān)系，并通過(guò)SJIT模型生成大量訓(xùn)練和測(cè)試數(shù)據(jù)，來(lái)訓(xùn)練、測(cè)試解碼器；此外，為了在感覺(jué)反饋通路缺失時(shí)，仍能準(zhǔn)確地恢復(fù)肢體運(yùn)動(dòng)功能，在設(shè)計(jì)長(zhǎng)短期記憶網(wǎng)絡(luò)解碼器的基礎(chǔ)上，結(jié)合數(shù)據(jù)驅(qū)動(dòng)的無(wú)模型控制（Model-Free Control，MFC）策略設(shè)計(jì)的輔助控制器構(gòu)造閉環(huán)BMI系統(tǒng)恢復(fù)關(guān)節(jié)活動(dòng)障礙者缺失的感覺(jué)反饋通路從而實(shí)現(xiàn)跟蹤期望軌跡的目的。

1 SJIT模型簡(jiǎn)介

Bullock等提出的SJIT模型［15］，如圖2所示。在模型中，GO為可變幅門信號(hào)，各神經(jīng)元組分別表示如下：DVV為期望速度向量（Desired Velocity Vector，DVV），OPV為輸出位置向量（Outflow Position Vector，OPV），OFPV為輸出力和位置向量（Outflow Force and Position Vector，OFPV），SFV為靜態(tài)力向量（Static Force Vector，SFV），IFV為慣性力向量（Inertial Force Vector，IFV），PPV為感知到的位置向量（Perceived Position Vector，PPV），DV為偏差向量（Difference Vector，DV），TPV為目標(biāo)位置向量（Target Position Vector，TPV），JPV為關(guān)節(jié)位置向量（Joint Position Vector，JPV），γD動(dòng)態(tài)γ神經(jīng)元（dynamic gamma motoneurons），γS靜態(tài)γ神經(jīng)元（static gamma motoneurons），α為α神 經(jīng) 元（alphamotoneuron），Ia為Ia型傳入纖維（type Ia afferent fibers），Ⅱ?yàn)棰蛐蛡魅肜w維（type Ia afferent fibers）。大腦區(qū)域包括區(qū)域4和區(qū)域5.

圖2 SJIT模型（“-”為抑制性反饋通路；無(wú)“-”為刺激性反饋）Fig.2 SJIT model

DV神經(jīng)元組計(jì)算TPV和PPV間的偏差。DV神經(jīng)元組平均放電活動(dòng)（average firing activity）ri由下式描述

式中0≤ri（t）≤1，下標(biāo)i對(duì)應(yīng)主動(dòng)肌（agonist muscle）。在后續(xù)介紹中，下標(biāo)i均對(duì)應(yīng)主動(dòng)??；而下標(biāo)j均對(duì)應(yīng)被動(dòng)肌（a ntagonistic muscle）。Ti為主動(dòng)肌i的TPV；xi（t）為PPV神經(jīng)元組平均放電率，這些神經(jīng)元可連續(xù)計(jì)算出主動(dòng)肌當(dāng)前位置；Br表示DV神經(jīng)元組基礎(chǔ)放電活動(dòng)（base firing activity）。

通過(guò)DV神經(jīng)元組連續(xù)計(jì)算偏差向量，進(jìn)而通過(guò)DVV神經(jīng)元組進(jìn)行比例放縮后得到DVV神經(jīng)元組平均放電活動(dòng)ui（t）

式中Bu為DVV神經(jīng)元組基礎(chǔ)放電活動(dòng)；g（t）為GO信號(hào)，并假定其來(lái)自于基底神經(jīng)節(jié)（basal ganglia）。DVV神經(jīng)元組的放電活動(dòng)僅僅在運(yùn)動(dòng)過(guò)程中進(jìn)行，且其平均放電活動(dòng)代表了相位移動(dòng)時(shí)間（phasic-movement time）。GO信號(hào)g（t）動(dòng)態(tài)為

式中 常量ò為慢積分率（slow integration rate）；g0為前腦決策中心的一個(gè)階躍輸入；常量C為GO信號(hào)飽和值。

OPV神經(jīng)元組可接收DVV和PPV神經(jīng)元組的信息，其平均放電活動(dòng)表示為

式中 η為比例因子。靜態(tài)運(yùn)動(dòng)神經(jīng)元組（static motoneurons）和動(dòng)態(tài)運(yùn)動(dòng)神經(jīng)元組（dynamic motoneurons）分別由γSi（t），γDi（t）表示，其平均放電活動(dòng)為

式中 ρ為放縮因子。Ia型和Ⅱ型肌肉傳輸纖維平均放電活動(dòng)由式（6）獲得

式中s1i（t）和s2i（t）分別為Ia型和Ⅱ肌肉傳輸纖維的平均放電活動(dòng)；pi為主動(dòng)肌i的位置；θ為靜態(tài)核袋纖維（static nuclear bag and chain fibers）的靈敏度；φ為動(dòng)態(tài)核袋纖維（dynamic nuclear bag fibers）的靈敏度。傳入纖維活動(dòng)飽和值可由式S（ω）=ω/（1+100ω2）給出。PPV神經(jīng)元組平均放電活動(dòng)xi（t）見(jiàn)式（7）

式中 τ為反饋延時(shí)，Θ為增益常量。IFV神經(jīng)元組平均放電活動(dòng)qi（t）見(jiàn)式（8）

式中 Λ為常數(shù)閾值。而SFV神經(jīng)元組平均放電活動(dòng)則由式（9）表示

式中 ψ為抑制性放縮參數(shù)（inhibitory scaling parameter）；h為常量增益，其用來(lái)控制外部負(fù)載補(bǔ)償（external load compensation）的大小和速度。OFPV神經(jīng)元組平均放電活動(dòng)ai（t）主要用來(lái)體現(xiàn)神經(jīng)元組的相位補(bǔ)償（phasic-tonic），其可由下式表示

α神經(jīng)元組平均放電活動(dòng)由式（11）表示

式中 δ為牽張反射（stretch reflex）增益。基于上述模型，臂部活動(dòng)可由下式描述

式中 主動(dòng)肌位置pi（t）始終位于主動(dòng)肌運(yùn)動(dòng)的近端和遠(yuǎn)端距離（origin-to-insertion distances）之間；類似地，pj（t）則表示被動(dòng)肌位置，且與主動(dòng)肌位置之間滿足pi（t）+pj（t）=1，此式是該模型的重要約束。Ei為施加于臂部關(guān)節(jié)的外部力。式M（ci（t），pi（t））=max｛ci（t）-pi（t），0｝為主動(dòng)肌i的生成力總和。M（ci（t）-pi（t））-M（cj（t）-pj（t））記作ΔM，ci（t）表示肌肉收縮活動(dòng)（contraction activity），其動(dòng)態(tài)表示為

式中v表征收縮率（contraction rate）的大小。

注意：該模型中，沒(méi)有考慮視覺(jué)反饋。在文中仿真中，均采用主動(dòng)肌位置pi（t）為肢體位置。

2 基于LSTM網(wǎng)絡(luò)設(shè)計(jì)的解碼器

對(duì)于健康的人而言，本身?yè)碛幸粋€(gè)完整的信息回路，如圖2所示，大腦區(qū)域4的DVV，OPV及OFPV神經(jīng)元組通過(guò)脊椎電路來(lái)傳遞肢體運(yùn)動(dòng)的指令，由感覺(jué)反饋部分進(jìn)行反饋。然而對(duì)于關(guān)節(jié)活動(dòng)障礙者而言，脊椎電路通路缺失，無(wú)法完成上述信息流閉環(huán)。因此則需要解碼器來(lái)恢復(fù)脊椎電路通路，即用解碼器來(lái)代替原有脊椎電路通路將大腦信號(hào)傳遞給假肢［16-18］。

2.1 解碼器設(shè)計(jì)

圖3 LSTM單元Fig.3 LSTM unit

LSTM網(wǎng)絡(luò)由Hochreiter S和Schmidhuber J在1997年提出，是為解決循環(huán)神經(jīng)網(wǎng)絡(luò)（Recurrent Neural Networks，RNN）“長(zhǎng)距離依賴”問(wèn)題的一種改進(jìn)算法［19-20］。LSTM的細(xì)胞單元結(jié)構(gòu)如圖3所示。

LSTM的細(xì)胞單元在RNN的基礎(chǔ)上加入了有3個(gè)門，分別是輸入門（Input Gate），遺忘門（Forget Gate），輸出門（Output Gate）。輸入門用來(lái)控制信息輸入，決定有多少將添加到記憶信息中；遺忘門用來(lái)控制細(xì)胞歷史狀態(tài)信息的保留，激活函數(shù)σf使得遺忘門的輸出值在［0，1］之間，當(dāng)遺忘門輸出接近0時(shí)的時(shí)候，表示將上一狀態(tài)的信息丟棄越多，接近1的時(shí)候，表示上一狀態(tài)的信息保留越多；輸出門用來(lái)控制信息輸出，決定多少記憶信息被用于下一階段的更新中。

上述過(guò)程可用以上公式表示，xt為輸入；ht為輸出；it為輸入門的輸出；ft為遺忘門的輸出；ct為當(dāng)前時(shí)刻t的細(xì)胞單元狀態(tài)；ot為輸出門的輸出，其中W和b為參數(shù)矩陣和偏置矩陣。

基于LSTM網(wǎng)絡(luò)設(shè)計(jì)解碼器來(lái)表征放電信號(hào)與肢體運(yùn)動(dòng)的非線性關(guān)系，選取DVV，OPV，OFPV神經(jīng)元組的放電活動(dòng)作為輸入，選取主動(dòng)肌與對(duì)抗肌的合力差ΔM作為L(zhǎng)STM網(wǎng)絡(luò)的輸出。

在輸入層，使輸入特征向量形式為x（t）=［x1（t），x2（t），x3（t），x4（t），x5（t），x6（t）］T并選取相對(duì)應(yīng)的放電活動(dòng)為：x1（t）=ai，x2（t）=aj，x3（t）=ui，x4（t）=uj，x5（t）=yi，x6（t）=yj.在輸出層，選取合力差作為輸出變量ht=ΔM，其中的隱層單元維數(shù)選取50，學(xué)習(xí)率選取為0.01.

2.2 數(shù)據(jù)生成及其離線測(cè)試

基于LSTM網(wǎng)絡(luò)設(shè)計(jì)的解碼器，其中的參數(shù)需要大量的數(shù)據(jù)訓(xùn)練來(lái)得到，所以此處利用SJIT模型生成大量的訓(xùn)練數(shù)據(jù)以及測(cè)試數(shù)據(jù)，如圖2所示，利用式（1）～（13）來(lái)構(gòu)建模型生成大腦皮層神經(jīng)元放電平均放電率、主動(dòng)肌與對(duì)抗肌的合力差、臂部的速度、位置以及關(guān)節(jié)轉(zhuǎn)矩等數(shù)據(jù)。生成的訓(xùn)練數(shù)據(jù)是對(duì)臂部關(guān)節(jié)的伸展任務(wù)進(jìn)行了共計(jì)1 600次模擬仿真，在每次臂部關(guān)節(jié)的伸展時(shí)，以采樣時(shí)間0.01 s對(duì)大腦區(qū)域4各神經(jīng)元組（DVV，OPV及OFPV）、主動(dòng)肌神經(jīng)元組以及被動(dòng)肌神經(jīng)元組的平均放電活動(dòng)進(jìn)行采樣，主動(dòng)肌與對(duì)抗肌的合力差ΔM（t）=M（ci（t），pi（t））-M（cj（t），pj（t））也被采樣記錄。因?yàn)槊看文M臂部關(guān)節(jié)的伸展耗時(shí)約為3s，所以經(jīng)過(guò)上述1 600次的模擬，可以得到共計(jì)480 000組的數(shù)據(jù)。文中，選取前470 000組數(shù)據(jù)進(jìn)行訓(xùn)練，后10 000組數(shù)據(jù)進(jìn)行測(cè)試。以下為具體變量和參數(shù)的設(shè)置：xi（0）=xj（0）=0.5，yi（0）=yj（0）=0.5，pi（0）=pj（0）=0.5，ui（0）=uj（0）=Bu，ri（0）=rj（0）=Br之外，其他變量的初始條件均設(shè)定為0.I=200，V=10，v=0.15，Br=0.1，Bu=0.01，Θ=0.5，θ=0.5，φ=1，η=0.7，ρ=0.04，λi=150，λj=10，Λ=0.001，δ=0.1，C=25，ò=0.05，ψ=4，h=0.01，τ=0.參數(shù)g0通過(guò)均值為0.75，方差為0.002 5的高斯分布隨機(jī)生成，生成后，其作為一組常數(shù)在仿真中存在。

圖4 解碼器離線測(cè)試Fig.4 Offline performance comparison

選取10 000組測(cè)試數(shù)據(jù)中［2 000，3 000］組的數(shù)據(jù)展示仿真結(jié)果，離線測(cè)試結(jié)果如圖4所示，由于測(cè)試數(shù)據(jù)通過(guò)采樣得到，是離散化數(shù)據(jù)，所以采用ΔM（k）表示，其中實(shí)紅線表示實(shí)際的ΔM（k），藍(lán)點(diǎn)畫線表示離線解碼的ΔM（k），并引入均方根誤差評(píng)價(jià)指標(biāo)（RMSE），作為解碼器的評(píng)價(jià)指標(biāo)。

式中 ΔˉM（k）為實(shí)際數(shù)據(jù)，即SJIT模型的采集的合力差，ΔM（k）為解碼數(shù)據(jù)。LSTM網(wǎng)絡(luò)設(shè)計(jì)的解碼器的均方根誤差值為4.3×10-2，表明LSTM網(wǎng)絡(luò)設(shè)計(jì)的解碼器良好的解碼效果。

3 閉環(huán)BMI系統(tǒng)構(gòu)建

感覺(jué)反饋通路缺失時(shí)，僅僅通過(guò)基于LSTM網(wǎng)絡(luò)設(shè)計(jì)解碼器難以準(zhǔn)確恢復(fù)肢體運(yùn)動(dòng)功能［15-17，22-23］。所以通過(guò)設(shè)計(jì)解碼器以及MFC輔助控制器，構(gòu)造閉環(huán)BMI系統(tǒng)跟蹤期望軌跡達(dá)到恢復(fù)肢體運(yùn)動(dòng)功能的目的。其中，期望軌跡為基于SJIT模型在加入解碼器情況下跟蹤目標(biāo)所產(chǎn)生的神經(jīng)元組PPV的平均放電率xi（t），如圖5所示。注意，文中人工反饋設(shè)計(jì)中沒(méi)有補(bǔ)償IFV和SFV神經(jīng)元組的反饋信息。

圖5 閉環(huán)腦機(jī)接口系統(tǒng)Fig.5 Closed-loop brain-machine interface system

3.1 輔助控制器設(shè)計(jì)

無(wú)模型控制（MFC）是一種只依賴于系統(tǒng)的輸入輸出數(shù)據(jù)，不依賴于被控系統(tǒng)模型的自適應(yīng)控制方法，并且不需要復(fù)雜的人工控制器參數(shù)整定［24-25］。選取控制輸入準(zhǔn)則函數(shù)如下

式中 φ（t）為被控系統(tǒng)的偽偏導(dǎo)數(shù)，它是一個(gè)時(shí)變參數(shù)；I（t）為被控系統(tǒng)輸入，即控制器輸出；ΔI（t）為I（t）的增量；λ為權(quán)重系數(shù)。在準(zhǔn)則函數(shù)中引入懲罰來(lái)限制輸入的變化。在腦機(jī)接口系統(tǒng)中，由于人類本體反饋的缺失，si和sj均等于0.

將公式（21）代入（20）中，對(duì)I（t）求導(dǎo)，并令其為0，可得

參數(shù)ρt的引入可加強(qiáng)公式（22）的通用性。參數(shù)λ的引入可避免出現(xiàn)分母為0的現(xiàn)象。由于偽偏導(dǎo)數(shù)φ（t）未知，不能直接用于式（22）可以利用相似的方法求出偽偏導(dǎo)的估計(jì)值

式中 ηt和μ與公式中的參數(shù)ρt和λ功能類似；即為求得的偽偏導(dǎo)數(shù)的估計(jì)值。將代入公式（22），便可得到控制器的輸出I（t）.根據(jù)設(shè)計(jì)好的腦機(jī)接口系統(tǒng)框架，用求得的I（t）代替單關(guān)節(jié)活動(dòng)信息傳輸電路模型公式（7-τ），即用公式（25）替換公式（7）。

4 仿真結(jié)果

首先為驗(yàn)證MFC輔助控制器的跟隨效果：基于上述閉環(huán)BMI系統(tǒng)框架，采用python進(jìn)行仿真實(shí)驗(yàn)，閉環(huán)BMI系統(tǒng)框架的各參數(shù)選取與2.2節(jié)參數(shù)相同。采樣時(shí)間選取0.01 s，g0=0.75，Ti（t）選?。撼跏嘉恢脼?.7，初始速度為-0.1/s，并以此速度移動(dòng)1 s后保持不變，MFC輔助控制器的相關(guān)參數(shù)選取ηt=1，μ=1，ρt=7，λ=2，^φ（0）=-2得到閉環(huán)BMI系統(tǒng)下PPV神經(jīng)元組的平均放電率跟蹤效果，如圖6所示，其中紅實(shí)線表示期望軌跡xi（t），藍(lán)點(diǎn)畫線表示MFC輔助控制器作用下的跟隨軌跡xi（t）.計(jì)算圖6中跟隨軌跡xi（t）的均方根誤差值RMSE以及相對(duì)誤差，均方根誤差為5.9×10-4，相對(duì)誤差如圖7所示，其最大值不超過(guò)1.2%.證明MFC輔助控制器的跟蹤效果良好。

圖6 PPV平均放電率對(duì)比Fig.6 Comparison of the average firing rate of PPV

圖7 PPV平均放電率相對(duì)誤差Fig.7 Relative error of average firing rate of PPV

驗(yàn)證構(gòu)建的閉環(huán)BMI系統(tǒng)有效性：采用構(gòu)建的閉環(huán)BMI系統(tǒng)與SJIT模型進(jìn)行比較。各參數(shù)選取同上。仿真結(jié)果如圖8所示，紅實(shí)線為SJIT模型的輸出軌跡，藍(lán)點(diǎn)畫線為文中構(gòu)建的閉環(huán)BMI系統(tǒng)的輸出軌跡，圖9為恢復(fù)軌跡的相對(duì)誤差。由仿真結(jié)果可知，恢復(fù)軌跡的相對(duì)誤差最大值不超過(guò)3%，說(shuō)明文中設(shè)計(jì)閉環(huán)BMI系統(tǒng)恢復(fù)了手臂的運(yùn)動(dòng)功能，證明構(gòu)建的閉環(huán)BMI系統(tǒng)的有效性。

圖8 手臂運(yùn)動(dòng)功能恢復(fù)效果Fig.8 Arm movement function recovery effect

圖9 恢復(fù)軌跡相對(duì)誤差Fig.9 Relative error of recovery trajectory

5 結(jié) 論

1）考慮腦電信號(hào)與時(shí)間序列存在相關(guān)性，提出基于LSTM網(wǎng)絡(luò)設(shè)計(jì)并訓(xùn)練解碼器，以此來(lái)恢復(fù)脊椎功能障礙者缺失的信息通路。

2）通過(guò)設(shè)計(jì)的解碼器并且結(jié)合基于MFC策略設(shè)計(jì)輔助控制器構(gòu)建的閉環(huán)BMI系統(tǒng)實(shí)現(xiàn)感覺(jué)反饋通路的恢復(fù)以及目標(biāo)軌跡的跟蹤。

3）仿真結(jié)果證明，基于LSTM網(wǎng)絡(luò)設(shè)計(jì)的解碼器離線、在線效果良好，基于MFC策略構(gòu)建的閉環(huán)BMI系統(tǒng)可實(shí)現(xiàn)對(duì)期望軌跡的跟蹤及缺失信息通路的恢復(fù)。