基于外側腓腸肌sEMG的帕金森神經調控步態量化測評方法

孫玉波，李海濤，舒智林，于洋，韓建達，梁思泉，于寧波

（1.南開大學 人工智能學院, 天津 300350; 2.南開大學 天津市智能機器人技術重點實驗室, 天津 300350; 3.天津市環湖醫院 神經外科, 天津 300350; 4.天津市環湖醫院 康復醫學科, 天津 300350）

帕金森病(Parkinson’s disease, PD)是一種常見的神經退行性疾病[1]。我國是世界帕金森第一大國，且病人數量隨人口老齡化持續上升。由于腦部黑質多巴胺能神經元的進行性缺失，導致患者出現運動癥狀如步態障礙和震顫，非運動癥狀如抑郁和焦慮等，嚴重影響患者的生活和自理能力[2-5]。

腦深部電刺激(deep brain stimulation, DBS)已被證明是一種有效的帕金森病治療手段[6-8]。術后DBS神經調控參數的調整需要量化、持續、精確的評估，對患者運動癥狀的動態測量與日常生活中的跟蹤測量也有助于評估患者的病情進展情況[9-10]。但是，目前這主要依賴于醫生的經驗和判斷，迫切需要客觀、量化的測評手段和分析方法。

在人體步行過程中足跟擊地時，腓腸肌的激活是保證腳踝穩定的重要因素，所以人體步行時腓腸肌表面肌電信號(surface electromyogram signal, sEMG)具有周期性特征[11]。根據腓腸激活時刻的檢測，可以對足跟擊地時刻進行識別，并提取步態周期時間信息。步態障礙是帕金森患者的典型運動癥狀[12]，步態量化測評對于帕金森病篩查、鑒別診斷和治療評價都具有非常重要的意義。有研究表明患者下肢腓腸肌肌電信號變化與患病程度密切相關[13-14]，因此帕金森患者的腓腸肌表面肌電信號分析對患者步態量化測評也具有重要意義[15-16]。

DBS神經調控治療，是對患者大腦直接施加電刺激，這可能對與步態直接相關的腓腸肌的神經指令造成影響。在步態研究的前向動力學方法中，神經指令作為輸入，直接確定了肌肉激活的幅度。這可以通過從肌電信號中提取的神經激活度與肌肉激活度進行量化[17]。

帕金森患者腦部黑質多巴胺能神經元的進行性缺失，會導致基底神經節輸出過度。基底神經節下行調節腦干活動，腦干是負責調節人體運動節律的重要部分，因此，運動節律失調是帕金森患者步態障礙的一個主要特征[18]。為了調整基底神經節的輸出，大部分針對帕金森患者的DBS電極刺激靶點為基底神經節中的丘腦底核[8]。針對帕金森的病理特征和DBS神經調控治療的特點，本文提出與步態節律密切相關的腓腸肌神經激活度與肌肉激活度變異性指標。

足底壓力可以充分反映步態。地面反作用力是人在行走過程中受到的唯一外力。在支撐相，足底與地面相互作用，地面反作用力對于步態有著根本性的作用，該信息可以通過力傳感器采集[19-21]。本文通過足底壓力進行步態在時域的準確識別，提取在臨床廣泛應用的步態周期時間變異性指標，對所提基于表面肌電信號的變異性指標進行驗證。

本文針對帕金森患者腦深部電刺激神經調控治療步態測評需求，分析步態與帕金森病理特點，提出基于外側腓腸肌(lateral gastrocnemius,LG)表面肌電信號的步態分析方法和量化測評指標，建立一個肌電和足底壓力信號同步采集與無線傳輸系統，并開展注冊臨床研究進行有效性驗證，為帕金森患者提供了一種輕量化的、基于外側腓腸肌sEMG的步態量化測評方案。本文的工作特色和創新具體如下：

1) 針對步態與帕金森病理特點，本文提出了基于外側腓腸肌表面肌電信號的步相識別方法與步態變異性指標提取方法。采用平滑先驗算法去除肌電信號中的低頻趨勢并保留帕金森臨床癥狀相關信息，通過小波變換與激活程度識別，檢測足跟擊地時刻、進行步相識別，進而對肌電信號進行劃分，計算神經激活度、肌肉激活度、及其變異性指標。

2) 帕金森患者的運動障礙表現極易受到儀器和場景等帶來的心理暗示的影響，為客觀測評帶來挑戰。為此，本文采用無線通訊、穿戴式、輕量化的肌電和足底壓力采集設備，建立數據同步采集方法，設計基于在臨床廣泛應用的10 m步行測試的實驗方案，構建面向帕金森患者的步態量化測評系統和范式。

3) 取得倫理許可，針對接受DBS治療的帕金森患者開展注冊臨床研究，基于時間同步的足底壓力信號對所提方法和指標的有效性進行了驗證。證明了本文所提步態量化測評方法的有效性和適用于帕金森患者的步態量化測評系統的可行性，為帕金森神經調控治療提供重要臨床手段。

1 系統設計與實驗范式

1.1 系統設計

人體在步行過程中受到的唯一外力是地面反作用力，而人的運動由肌肉驅動，因此肌電信號和足底壓力中包含著步態的關鍵信息。基于帕金森的病理特點，病人步態障礙極易受到聲音、視覺等提示的影響。為了減小對病人的心理暗示和干擾，實現對患者自然真實、客觀量化的步態測評，本文構建了一個采用無線通信、穿戴式、輕量化的肌電和足底壓力信號同步采集系統。系統設計示意圖如圖1所示。

圖1 系統示意圖Fig.1 System diagram

本文采用Delsys trigno avanti設備(delsys incorporated, USA)記錄下肢外側腓腸肌表面肌電信號，采樣頻率為2 000 Hz。采用無線壓力鞋墊(ZIGUN AVOIN company, Finland)采集患者在行走過程中的足底壓力，采樣頻率為100 Hz，最大測量值為1000 N, 分辨率為0.1 N/cm，壓力鞋墊質地柔軟且厚度不超過1.2 mm，可以方便、舒適地放到病人的鞋中。

兩種信號采樣頻率不同，并且是由兩臺獨立的計算機分別進行采集。為此，需設計一個能夠實現信號同步的打標軟件模塊，如圖1所示。圖1中左側虛線箭頭表示的均為無線連接通信，而中間和右側實線箭頭表示的均為有線連接通信。

同步軟件單獨在一臺計算機中運行，程序為Windows Form窗體程序，在Microsoft Visual Studio 2019中基于C#語言開發。通過串口分別連接到表面肌電和足底壓力設備對應的計算機。按對應的波特率（肌電：110；壓力鞋墊：115 200）和數據位（肌電：5；壓力鞋墊：8）向對應設備發送一個字節數據，作為同步信號。

1.2 實驗范式

基于在臨床廣泛應用的10 m步行測試，設計測評實驗范式。受試者在起點處靜止站立，在聽到開始指令后向前自然行走，抵達終點后停止。通過無線足底壓力測量單元和表面肌電測量單元對受試者的步態數據進行采集。

步行過程可分為3個部分：啟動、步行和停止。對于神經功能障礙患者，其啟動和停止部分的步態易受到多種因素影響。因此，本文將每次測試的前兩個步態周期（啟動）和最后兩個步態周期（停止）的數據剔除，僅采用中間步行部分的數據進行分析，如圖2所示，從左到右，步行過程劃分為3個部分：啟動、步行和停止。

圖2 基于標準臨床10 m步行測試的實驗范式Fig.2 Experimental paradigm based on the standard clinical 10-meter walking test

2 信號處理與特征提取

通過建立的肌電和足底壓力信號同步采集系統，在臨床10 m步行測試實驗過程中同步采集數據，通過無線傳輸至計算機，隨后分別對兩種信號進行分析與特征提取，流程如圖3所示。

圖3 步態變異性指標提取流程Fig.3 Diagram of gait variability features extraction

對時間對齊后的表面肌電信號進行濾波等預處理，然后基于小波變換進行足跟擊地時刻識別，提取步態周期時間，隨后對肌電數據進行劃分并提取神經激活度、肌肉激活度及其變異性指標。基于足底壓力數據,在時域對步態進行準確識別，進而計算步態周期的變異性系數，驗證肌電信號變異性指標的有效性。

2.1 外側腓腸肌sEMG信號分析與特征提取

在表面肌電信號預處理過程中，需要去除肌電信號的運動偽跡以及由電極與皮膚的不穩定接觸造成的低頻趨勢。通常做法是對表面肌電信號進行截止頻率不超過20 Hz的高通濾波。然而，對帕金森患者來說，肌電信號低頻成分中存在著與其臨床特征相關的重要信息。因此，在去除低頻趨勢的同時，要避免相關信息的丟失。為此，本文采用平滑先驗算法。

對時間對齊后的表面肌電信號的預處理流程如圖4所示，分為3個步驟：通過平滑先驗算法去趨勢、對肌電信號進行幅值歸一化處理、全波整流。

圖4 表面肌電信號預處理流程Fig.4 Preprocessing workflow for sEMG signals

假設表面肌電信號的時間序列由有效信號和低頻趨勢成分兩部分組成：

式中：zs是去掉趨勢項之后的表面肌電信號序列；zt是表面肌電信號中的低頻趨勢成分。

對有效信號的估計可以表示為

式中 α為正則化參數。D2的取值為

通過平滑先驗算法去趨勢后，進行幅值歸一化處理。歸一化尺度一種可行的選擇是表面肌(surface electromyogram, sEMG)信號絕對值中的最大值[22]。但為了排除異常值和采集過程中由于運動而產生的高幅值異常值對激活時刻識別的影響，歸一化尺度參數選擇為肌電信號50個最大絕對值中的最小值[23]。經過全波整流后，得到去除了基線漂移、低頻趨勢和高頻噪聲的肌電信號；再通過低通濾波器得到線性包絡。

在人體步行過程中足跟擊地時，腓腸肌的激活是保證腳踝穩定的重要因素。外側腓腸肌表面肌電信號具有與足底壓力信號相同的周期性特征，由傳感器采集到的原始足底壓力信號與同側下肢的外側腓腸肌表面肌電信號如圖5所示。本文提出了一種基于小波變換和閾值設定的步相識別方法。通過對外側腓腸肌的激活狀態識別，確認足跟擊地時刻，進而識別步相。

圖5 一名帕金森患者的時間對齊的足底壓力原始信號與外側腓腸肌表面肌電原始信號Fig.5 Synchronized raw data of plantar pressure signal and LG sEMG signal in a PD patient

小波變換可以同時保留時域和頻域的信息，有研究表明小波變換在鑒別肌電信號模式方面優于其他時頻方法，并且相比傳統方法提供了更好的時頻分辨率，更加適用于不穩定的肌電信號分析[24-25]。

使用離散小波變換代替均方根方法，對表面肌電信號進行處理，可以更加準確地識別肌肉激活時刻，去除噪聲影響。經過離散小波變換，高頻噪聲分量主要集中在細節系數中，而近似系數主要表示信號的低頻信息。在識別激活時刻過程中，主要關注的是肌電信號幅值的低頻變化，關注低頻信息。在肌肉激活程度識別過程中，離散小波變換方法具有更高的準確度[26]。

基于小波的方法是分析非平穩信號如肌電信號的實用技術。離散小波變換(discrete wavelet transform, DWT)對母小波進行縮放和移位，并將離散時間信號分解成一組信號。DWT使用高通和低通濾波器與降采樣方法，將信號分解成多分辨率系數。原始肌電信號S經低通濾波器和高通濾波器，獲得第一級的近似系數(cA1)和細節系數(cD1)。為了得到多分辨率子集，對上一級的近似系數進行降采樣和濾波，并重復此操作，獲得近似系數 cAn，細節系數 cD1，cD2，…，cDn。離散小波變換分解流程如圖6所示。

圖6 表面肌電信號分解樹Fig.6 Tree decomposition of the sEMG signal

相關研究表明基于Haar小波的DWT可以更加有效地對不同動作進行識別，對肌肉活動檢測準確度較高[26]。并且基于Haar小波的小波變換方法可以更有效地去除肌電信號中的噪聲[27]。因此本文中，對預處理后的外側腓腸肌表面肌電信號進行離散小波變換，采用Haar小波作為母小波，進行5級分解，獲得近似系數cA5。隨后通過3 Hz的巴特沃斯低通濾波器，獲得近似系數的線性包絡，再通過閾值方法，檢測激活狀態，從而實現對足跟擊地時刻的識別。通過外側腓腸肌表面肌電信號進行足跟擊地時刻識別的具體流程與對應信號示例如圖7所示。

圖7 帕金森患者外側腓腸肌sEMG信號激活識別處理過程Fig.7 Schematic diagram of LG sEMG activation recognition process in PD patients

本文中，使用了已經得到廣泛認可的肌肉活動檢測閾值設定方法[28]，閾值T由式(4)所得。其中，本文取患者在靜止站立時長度為2 s的外側腓腸肌表面肌電信號作為基線，預處理后提取其平均值M和標準差STD。將近似參數包絡線中超過閾值T的時刻定義為激活時刻。從非激活向激活過渡的時刻定義為足跟擊地時刻。

隨后，將經過預處理與低通濾波求出的線性包絡e(t)，按照足跟擊地時刻進行劃分。再通過神經激活模型計算每段信號的肌肉神經激活度u(t)。

2.2 神經激活度與肌肉激活度

肌電信號的線性包絡e(t)轉化為神經激活度u(t)的過程稱為激活動力學。當肌肉纖維被單個動作電位激活時，肌肉就會產生抽動反應。這個響應可以用臨界阻尼線性二階微分系統很好地表示，如式(5)所示。這種類型的響應是根據肌電信號的輸入e(t)來確定神經激活度u(t)的基礎。

式中M、B和K是定義二階動力學系統的常數。在實驗室中采集到的數據是離散的時間序列信號，離散方程表示為

式中：d是采集設備的延遲，由于所采用設備延遲為48 ms，采樣周期0.5 ms，所以d取96； β1和 β2是定義二階動力學系統的參數。這些參數將e(t)映射到神經激活度u(t)。參數 α、 β1和 β2滿足以下條件：

為了表征神經激活度與肌肉激活度之間的非線性關系，引入一個函數：

式中a(t)是計算所得的肌肉激活度。

2.3 足底壓力數據處理與步態識別

將時間對齊后的足底壓力信號通過一個截止頻率為20 Hz的低通濾波器，消除高頻噪聲。基于在足底壓力信號中識別出的周期性特征進行步態周期識別。

在進行步態識別和劃分時，考慮到患者體重的特異性以及足底壓力信號在擺動相過程中具有低頻噪聲，設定區分擺動相和支撐相的閾值為體重的5%。記錄每個步態周期擺動相的結束時刻，并作為下一個步態周期支撐相的開始時刻，如圖8所示，為一名PD患者的足底壓力曲線與步態識別示意圖。隨后可提取出每個步態周期的時間長度。

圖8 基于足底壓力信號的步態識別Fig.8 Gait recognition based on plantar pressure

在本文中，根據足底壓力信號劃分步態周期，計算得到步態周期時間變異性指標，與表面肌電信號變異性指標進行交叉驗證。

2.4 步態變異性指標

運動節律失調是帕金森患者步態障礙的一個主要特征，并且DBS神經調控可能會對與步態相關的外側腓腸肌的神經指令造成影響。因此，確定步態變異性作為主要測評指標。

步態周期時間的變異性是公認可以反映帕金森患者運動障礙的評估指標。具體計算方法如式(9)所示，由步態周期時間的平均值和標準差計算得到：

計算肌電信號激活度變異性時，令步態周期個數為n。求得每個步態周期相應的平均肌肉神經激活度u(k)，k∈[1,n]。

式中：t∈[1,mk]；mk為第k個步態周期的樣本長度。計算神經激活度變異指標：

對肌肉激活度a(t)計算變異指標的過程與u(t)相同：

3 臨床實驗與結果分析

3.1 倫理與受試者招募

本研究獲得天津市環湖醫院倫理委員會批準(2019—35)，并在中國臨床試驗注冊中心注冊(ChiCTR1900022715)。入組標準為：

1） 需要重新調整DBS治療參數；

2） 年齡40歲及以上；

3） 學歷小學及以上；

4） 能夠自主站立和行走；

5） 無其他會對步態造成影響的肢體或肌肉相關疾病。

本研究中，5名接受了腦深部電刺激手術、并符合入組標準的帕金森患者被招募，其中包括2名男性和3名女性，基本信息如表1所示。所招募的受試者平均年齡61.8歲，平均身高165.0 cm，平均體重67.6 kg。

表1 招募病人的基本信息Table 1 Information of PD patients

實驗開始前，每位帕金森病患者都被充分告知實驗目的和步驟，并簽署了同意書。患者的臨床特征、所使用的DBS儀器種類、初始DBS刺激參數以及神經調控過程中所使用的刺激頻率等詳細信息見表2。

表2 病人DBS神經調控初始參數Table 2 Patients’ initial DBS programming parameters

3.2 步態量化測評結果與分析

本文對5名帕金森患者在接受DBS神經調控時的步態變異性進行了量化測評。根據表面肌電信號，計算出外側腓腸肌肌電信號變異性指標。根據足底壓力信號，識別步態周期，計算步態周期時間的變異性指標。在不同DBS神經調控頻率下，患者的步態變異性指標有顯著不同，具體結果如表3所示。

表3 患者神經調控測試中提取的步態周期時間變異性指標與外側腓腸肌肌電信號變異性指標Table 3 CV of step time and GL sEMG under different DBS parameters

在由步態周期時間變異性指標與外側腓腸肌變異性指標構成的二維空間中，特征點變化呈現出一定的趨勢，如圖9、10所示。圖9、10中藍色直線為最小二乘法線性擬合結果，以表示二者變化趨勢，陰影部分為95%置信區間。由此可見，患者的特征點離原點越近，即步態周期時間變異系數與所提肌電變異系數越小，患者步態變異性越小，步態越穩定。所提指標可以為步態量化測評提供參考。

采用皮爾遜相關性系數分析所提出的神經激活度和肌肉激活度變異性指標與步態周期時間變異性之間的關系，相關性系數分別為0.602(p=0.000)和0.591(p=0.000)，具有顯著相關性。由此，說明了本文所提根據外側腓腸肌肌電信號提取神經激活度與肌肉激活度變異性指標的方法的可行性，確切證實了所提變異性指標的有效性。且本文所提的2種指標相較于傳統的步態周期時間變異性指標，對患者的狀態具有更好的分辨率。如圖9、10所示，步態周期時間變異性的變化范圍較小，而所提表面肌電信號變異性指標則具有更大的變化范圍，可以更好地對步態變異性進行分辨。

圖9 外側腓腸肌神經激活度變異性指標與步態周期時間變異性指標呈線性趨勢Fig.9 Linear regression between the CV of step time and LG neural activation

圖10 外側腓腸肌肌肉激活度變異性指標與步態周期時間變異性指標呈線性趨勢Fig.10 Linear regression between the CV of step time and LG muscle activation

在患者的連續穩定步態周期中，基于肌電信號提取的神經激活度、肌肉激活度在根據步相劃分后，對患者每個步態周期的量化，如圖11所示。圖11左側為患者在連續穩定步態過程中的足底壓力信號與對應的神經激活度、肌肉激活度曲線，右側為對應的計算變異性指標的步態周期時間、平均神經激活度和平均肌肉激活度柱狀圖，與對應提取的變異性指標。相較于步態周期時間長度，所提指標具有更好的區分度。

圖11 一名帕金森患者的足底壓力信號、神經激活度曲線、肌肉激活度曲線與提取出的變異性指標Fig.11 Plantar pressure signal, neural activation curve and muscle activation curve and the extracted variation index of a PD patient

針對帕金森患者神經調控步態量化測評中步態變異性的評估，所提的基于外側腓腸肌sEMG信號采集、無線傳輸與分析的步態量化測評方法，能為醫生提供更加有區分度的指標，為臨床診斷提供輔助。

4 結束語

本文針對帕金森患者腦深部電刺激神經調控治療步態量化測評需求，分析帕金森病理與步態特點，提出了一種基于外側腓腸肌sEMG數據的步態量化測評方法。建立一個肌電和足底壓力信號同步采集與無線傳輸系統，基于在臨床廣泛應用的步行實驗范式，開展臨床研究并進行有效性驗證。通過時間同步的足底壓力信號，驗證了所提測評方法和指標在帕金森病患者步態量化測評中的有效性。且所提神經激活度與肌肉激活度變異性指標，相較于傳統指標具有更好的區分度。由此，為帕金森神經調控治療提供了一種重要臨床手段。