基于Kinect的 帕金森病步態不對稱性識別方法①

張幼安，侯振杰，坎標，姚恩，張家瑋

1.常州大學華羅庚學院，江蘇常州市213164；2.常州大學信息數理學院，江蘇常州市213164；3.常州大學機械工程學院，江蘇常州市213164

帕金森病(Parkinson's disease,PD)是一種特征性運動障礙，臨床存在多種表現，如慌張步態、凍結步態、靜止性震顫、肌僵直等。作為老年人多發病，國內外學者做了很多努力預防或治療此病及其并發癥[1-9]。Gray等[10]研究顯示，帕金森病患者軀體運動功能下降，尤其是步態障礙與高死亡風險相關。

臨床對帕金森病患者病情的評估主要采用統一帕金森病分級量表(Unified Parkinson's Disease Rating Scale,UPDRS)，但該方法參數較多，較為繁瑣[11-13]。實驗室一般通過穿戴設備或測量足底壓力，需要鋪設專業行道或在腳底墊測力墊，設備價格昂貴，也影響行走步態[14]。有很多學者通過視頻圖像提取數據，但建模過程一般較為繁瑣，無法推廣[15-17]。

鑒于步態障礙的主要特征為步態不對稱性，本研究提出一種基于Kinect骨骼數據的帕金森病步態不對稱性識別方法；針對帕金森病“N”字形發病規律，結合Kinect采集患者骨骼點數據，通過聚類算法得到行走加速度分布，使用字符串相似度的矩陣算法[18]，比較左右腳加速度分布，量化步態對稱性指標，最后建立對稱步態隱馬爾科夫模型(Hidden Markov Model,HMM)完成常規識別工作。

1 對象與方法

1.1 研究對象

常州市武進區南夏墅敬老院及湖塘鎮南夏墅村志愿者18例，采集時間為2016年7月至8月，調研過程為常州大學重點暑期社會實踐項目。

健康組納入標準：①年齡≥40周歲；②無神經內科疾病史；③無腰痛、關節痛等導致下肢活動異常的癥狀；④無明顯殘疾或步態障礙。

患者組納入標準：①年齡≥40周歲；②經確診有帕金森病或有明顯帕金森綜合征癥狀，符合中國帕金森病的診斷標準；③可順利步行。

健康組10例，其中男性7例，女性3例；年齡40～65歲；身高158～176 cm。患者組8例，均為男性；年齡55～70歲；身高168～177 cm。

所有研究對象均知情同意。

1.2 方法

采用微軟公司Kinect V2.0，包含RGB彩色攝像頭、紅外線發射器和紅外線CMOS攝像頭；每秒最多獲取30幀圖像；可識別距離0.5～4.5 m。

Kinect骨骼跟蹤技術通過處理景深數據，建立人體各個關節的坐標。骨骼追蹤能確定人體各個部分，如頭、腿等。骨骼跟蹤產生X、Y、Z軸坐標，確定這些骨骼點的位置。因此，Kinect提供的骨骼數據為行走過程中身體各關節點的空間坐標序列。使用時通過調用自帶SDK進行實驗對象行走過程中的全身骨骼點采集。

設立步道長7 m、寬2 m，排除外界干擾。將Kinect平放在步道一端。受試者以正常舒服姿勢站立于Kinect正前方6 m處，操作者對Kinect進行距離和角度等參數校準。當Kinect指示燈亮起，受試者以最舒適狀態面向Kinect走去，直線行走6 m或至無法識別人體關節點坐標后停止。每例受試者測試6次。

采集數據為全身25個關節點的三維坐標序列，通過該骨骼點序列提取步態周期，繼而提取步態參數，進行對比實驗，驗證設備可靠性。計算行走過程中位移和加速度序列，聚類分析患者左右腳步態對稱性，最后完成常規識別。

1.3 觀察指標

采用MATLAB軟件編程提取步態特征。提取步態參數包括步長、步幅、步寬和聚類后的加速度序列。

1.3.1 步態特征提取

先由空間坐標序列的周期性變化計算步態周期，再通過步態周期計算各項步態參數。

1.3.1.1 步態周期計算

使用Kinect采集的人體骨骼點空間變化坐標序列。由于Kinect攝像頭擺放無法使X軸與地面完全平行，也無法以Y軸的0值視為落地點，所以不能根據返回坐標(如Y軸的0值)直接確定雙腳落地點來確定周期。

行走過程中，人體中心點呈現左右、上下擺動的周期性變化。本文將每次擺動的最低點作為特殊點，出現特殊點的幀作為特殊幀。通過求特殊幀位置即可計算步態周期。步態周期計算流程圖如圖1。

1.3.1.1.1選取中心點

行走過程中，上半身產生的彈性變形接近0，彎腰程度接近0，上半身可近似剛體；而下半身由于運動過于復雜，無法準確描繪中心點周期性變化；結合Kinect V2.0的采集點，取脊柱底部點作為中心點(圖2)。

1.3.1.1.2數據低通濾波

采集取中心點Y軸坐標序列。由于Kinect采樣頻率默認為30幀/s，因此時間域和頻域一致。做該點Y軸坐標頻域變化曲線圖，手動提取步態穩定區間，一般為前80%點。

第1次濾波：做相鄰數據差值，若該值大于前n組差值和，則將該值置為前n組差值的均值。由于Kinect采集性能較穩定，n置為5。

第2次濾波：對坐標序列行RLOWESS濾波，得到新的Y坐標變化序列。

1.3.1.1.3求特殊幀

根據濾波后Y坐標序列，做變化曲線圖。計算曲線上每一點的瞬時變化率，得到瞬時變化序列。在13～17之間依次取數作為周期基數，套取每個周期內的瞬時變化率最小值。根據相鄰套取結果的差值判斷：若與其他相鄰套取結果差之間存在倍數差關系，則該次套取失敗，重新選循環基數，直致篩選出正確特殊點，并求得特殊點序列。特殊點所在幀數即特殊幀位置。

人正常行走速度 V=4～4.7 km/h，即 1.1～1.3 m/s，以S=0.6 m作為正常跨步長，計算步頻為1.85～2.20/s。

Kinect采集頻率30幀/s，所以n=13.64～16.22幀/周期，因此Kinect采集的步態周期數在13～17之間，所以將13～17作為周期基數。

1.3.1.1.4計算步態周期

步態周期(T)為行走過程中一側足跟著地至該側足跟再次著地時所經過的時間，即連續3幀特殊幀出現所經歷的時間。

其中，tk為第k個特殊幀出現的時間，tk+2為第k+2個特殊幀出現的時間。

圖1 步態周期計算流程圖

圖2 Kinect V2.0的采集點

1.3.1.2 步態參數

步態參數定義見圖3，步長、步幅、步寬用來驗證Kinect采集步態的可靠性，加速度序列用來進行對稱性檢驗和步態識別。

圖3 步態參數定義

1.3.1.2.1步長

行走時一側足跟著地到緊接著的對側足跟著地的距離。使用空間點坐標距離公式計算：

其中，spi為第i次跨步步長，Xk、Yk、Zk為特殊點序號為k的空間坐標，Xk+1、Yk+1、Zk+1為特殊點序號為k+1的空間坐標。

1.3.1.2.2步幅

行走時由一側足跟著地到該側足跟二次著地所行進的距離。

1.3.1.2.3步寬

先以一側足(本文選右足)的2個落地點和1個懸空點，通過三點式方程構造平面空間；然后計算對側足(本文選左足)落地點至該平面的距離。

兩個右腳落地點R1(X1,Y1,Z1)和R2(X2,Y2,Z2)，懸空點R3(X3,Y3,Z3)建立的空間平面方程：

對側足落地點到平面的距離，即步寬：

其中，(X0,Y0,Z0)為對側落地點坐標。

1.3.1.2.4位移

位移的計算使用下一時刻點坐標至該點坐標的直線距離。

1.3.1.2.5速度

每次采集捕獲骨骼點的瞬時速度。由于采集頻率為30幀/s，所以相鄰幀時間間隔為0.033 s，用相鄰3個點的平均速度代替中點的瞬時速度，中點的瞬時速度：

其中，fc為 Kinect采集頻率，δt-1為點 Bt-1到點 Bt的距離，δt為點 Bt到Bt+1的距離。

1.3.1.2.6加速度

相鄰幀速度的變化率：

根據采集到的m+1個點的坐標，可以確定時間參數序列，即位移序列X1,X2,X3,……Xm；速度序列V1,V2,V3,……Vm-1；加速度序列a1,a2,a3,……am-3。

1.4 對比實驗

在地面鋪設白紙，使用米尺測距、秒表測時，采集一次步態周期實際步長、步幅、步寬，并計算速度。

1.5 對稱性檢驗

臨床研究表明[19]，出現凍結步態的帕金森病患者姿勢控制更差。帕金森病患病呈“N”字型進展：首先單側腿出現病癥，隨后進展至同側上肢，再發展至對側腿部，再至對側上肢。推測初期單腿患病時，肌肉控制力下降，患者會多使用健康腿，導致步行周期左右腿加速度不同；對比兩腿加速度序列便可以比較兩腿步態對稱性。

本研究通過比較左右腳加速度分布字符串的相似度，確定步態的對稱性，對比參數為左右膝的加速度序列。在確保加速度序列誤差足夠小的基礎上，分別對左、右序列進行K-means聚類，得到左、右一次聚類字符串；取一側聚類結果不變，用另一側所有可能的聚類結果字符串與前一側聚類字符串用字符串相似度矩陣算法，求相似度，取相似度最大值作為實際對稱性指標。

1.6 常規識別

建立對稱步態HMM模型進行常規識別。根據中心點波形分析，左右腳加速度分布會表現在中心點的加速度分布上。

收集健康組中心點加速度序列，使用Baum-Welch算法訓練HMM模型參數；通過前向算法統計不同程度帕金森病患者步態序列的發生概率，提取最小概率作為閾值。計算受試者步態序列相對于健康人HMM模型的發生概率，若在患者閾值范圍內，則記為患者。

在應用HMM模型預測步態加速度分布時，我們將已知加速度序列看作是可見狀態鏈，將人體腿部各動作激勵的肌肉或神經信號作為隱含狀態鏈。

2 結果

2.1 可靠性檢驗

使用Kinect采集的人體步行1次骨骼數據，其中心點Y軸坐標頻域如圖4。

圖4 中心點Y軸坐標頻域變化

取前80%點作為穩態序列，對穩態序列兩次濾波，作出穩態序列圖像，并在同一圖上作出濾波后的曲線，如圖5。

對曲線中每一點求斜率，并作斜率圖，如圖6。

圖5 穩態數據及濾波圖

圖6 斜率圖

一名測試者Kinect采集的步態參數，與手工測量步態參數比較見表1。Kinect采集的三種步態參數誤差均在5%以內，提取步態參數的方法可以信任。見表1。

表1 兩種測量方式比較(cm)

2.2 步態對稱性檢驗

展示1例健康人數據和3例最有區別度的患者數據。由于每次實驗采集骨骼點超過100條序列，僅顯示前7次數據，如表2。其中患者1為拖沓步態，患者2為輕度跛行，患者3為單腿萎縮。

表2 膝關節加速度原始數據(m/s2)

加速度狀態可分為5類：加速度=0、加速度減小的加速運動、加速度增大的加速運動、加速度減小的減速運動、加速度增加的減速運動。因此，K-means聚類時設K=5，結果見表3。

表3 膝關節加速度聚類結果

左右腳狀態相同的速度聚為同一類，將一側聚類后的結果進行循環種類號替換，聚類結果不同的可以任意替換(如1234可任意替換為1324、4123等)，聚類結果相同的只能替換成相同的聚類效果(如1233只能替換為2144、4311等)。一側加速度序列有120種相似度結果，提取相似度最大值作為該受試者的實際相似度。結果見表4。

表4 健康人與帕金森病患者左右腳加速度分布相似度比較

健康人左右相似度最大，患者1、患者2、患者3病情依次加重，左右相似度依次減小。

2.3 常規性識別

采用Baum-Welch算法學習HMM模型，加速度序列通過DBSCAN聚類，結果有6類，即HMM模型中Markov鏈每個狀態對應可能的觀察值數目為6；預設Markov鏈狀態數目N，初步訓練得到模型參數，再運用前向算法計算由該初步模型產生加速度序列的概率P；不斷重設N值，重復該步驟至概率收斂，得最佳N值為5。

結果：M=6，N=5

其中，M為每個狀態對應的可能的觀察值數；N為模型中Markov鏈狀態數；A為狀態轉移概率矩陣；B為觀察值概率矩陣；Pi為初始狀態概率矢量。

采集計算多組帕金森病癥患者的加速度序列，聚類標準化，前向算法多次計算后

某受試者經加速度序列聚類標準化后，前向算法計算得logP=-4.2×10-2，高于閾值，則判斷該受試者有可能患帕金森病。

對8例患者重新采集數據，共8組26條數據(一些患者由于身體原因無法重復多次采集)，并提取出能代表其特征的7組14條數據，成功識別出12條，成功率85.71%。提取了10組60條健康人數據，在有效的46條數據中，成功區別出35條，區別度76.09%。

3 討論

根據帕金森病不對稱步態進行患病識別一直受到關注。王人成等[20]采用基于FITTS定律提出的二維極坐標下的步態對稱性指標；Shorter等[21]采用偏差閾值分析，比較不對稱行走的屈伸行為與未裝支具組的標準，使用該方法可量化雙側關節的對稱性和不對稱性；鄒亞等[22]根據可穿戴慣性傳感器觀察參數提出對稱性計算式。

多數學者主要以左右腳參數比值做主要因式[23-24]，但對于不對稱步態的識別依據一直缺乏獨到且可執行的解釋，尤其是帕金森病患者不對稱性步態。廖福元等[25]將患者左、右腿支撐時間和跨步時間排成序列對，以多個尺度逼近序列對的相似度來量化步態的對稱性，但沒有結合病理提出識別原因。

本研究基于帕金森病“N”字形發病規律，將使患者患病初期產生較為明顯的不對稱步態，通過Ki-nect體感設備建立人體連桿模型，采集關節點，使設備對正常行走步態的影響降到最低。

本研究通過骨骼點三維坐標序列類正弦曲線的周期性特征，結合曲線斜率大小求波峰、波谷的方法提取周期，然后提取步態特征參數。對比足跡測量，證明該提取方法可行。周期提取還可以通過計算斜率方差的方式提取波峰、波谷，再提取周期。在提取周期后，除可計算本文所用的特征參數之外，還可以通過空間坐標的時間變化曲線提取其他步態參數，如步頻、步速、支撐相時間等。

對于帕金森病的步態，本文基于帕金森病患者左右腳加速度分布不同的假設，通過字符串相似度的矩陣算法識別不對稱步態。由于K-means聚類方法對噪聲不敏感，誤差可能較大；也可選擇使用DBSCAN聚類，但需人工多次干預，直到聚類到所需的聚類書目。在使用字符串相似度的矩陣算法時，由于需要對一側聚類結果進行窮舉，再與另一側進行匹配算法，較為復雜，有待改進。最后，我們使用HMM預測模型進行患者和非患者的識別，證明有一定效果。但HMM模型的訓練常需要較大數據，因此需要插值算法，或需要將幾次采集結果組合，會產生一定誤差，還需要進一步研究。