CN110151189A - 用于帕金森步态风险评估的非线性步态动力学判别方法 - Google Patents

Abstract

本发明公开了一种用于帕金森步态风险评估的非线性步态动力学判别方法。本发明方法提出了基于非线性动力学特征的多维特征判别图以用于帕金森患者的步态风险评估，联合了支持向量机、非线性动力学特征和图形化显示方法，将5维的非线性动力学特征在二维平面中以图形的方式展示。本发明方法适用于帕金森患者的步态风险评估。在正常的人群中，如果非线性步态动力学判别图的多个指标都超出判别图的分界点，则说明其患病的概率或者病症的严重程度将大大的提升，需要引起个体自身的密切关注以及医生的进一步检查。

Description

用于帕金森步态风险评估的非线性步态动力学判别方法

技术领域

本发明属于模式识别技术领域，具体涉及一种用于帕金森步态风险评估的非线性步态动力学判别图。

背景技术

帕金森病是一种多发于中老人群、以运动障碍为主要临床症状的神经退行性疾病，其临床上以静止性震颤、运动迟缓、肌僵直和姿势步态障碍为四大主征。目前帕金森病的治疗方式主要包括药物治疗和手术治疗，在患病早期，通常在临床中可以通过药物很好地控制由于患有帕金森病而导致的物理行为，如静止性震颤等临床症状，但是随着帕金森病患者病情的发展，患者的步态异变、震颤等症状会给患者的生活质量带来严重的影响，需要不同程度的手术介入治疗。如果能够正确评估患病风险及严重程度，依此及时、合理地指导治疗，这样对帕金森患者及可疑患者的病症控制、生活质量提高等方面都具有显著性的意义。

在基于步态信号的帕金森病诊断和判别领域中，已有方法大多是提取与病症有关联的时频参数，输入到各种不同的分类器，以实现帕金森诊断和风险判别的目的。而人体步态信号本质上可以看作是由复杂动力学系统产生的复杂非平稳信号，仅仅采用时频参数显然不足以刻画帕金森步态中蕴含的动态本质特性。

发明内容

针对上述背景和存在的问题，本发明提出一种用于帕金森步态风险评估的非线性步态动力学判别方法，提供一种更为形象的直观的图形，将帕金森患者的多维非线性步态动力学特征参数映射到二维空间中以图形方式显示，可在一定程度上提醒临床医生帕金森病的高危人群及其危险程度。

本发明的具体技术方案通过如下步骤实现：

步骤一、步态时间参数数据获取：

利用放置在试验者鞋底的力敏薄膜开关来获取人体步态时间参数数据，包括左脚跨步间隔、右脚跨步间隔、左脚站立间隔、右脚站立间隔、双脚支撑间隔这5列时间序列步态参数；具体包括左脚跨步时长L_stride,右脚跨步时长R_stride,左脚站立时长L_stance，右脚站立时长R_stance,双脚支撑时长D_support，构成一组特征矢量x：x＝[x₁,x₂,x₃,x₄,x₅]^T＝[L_stride.R_stride.L_stance,R_stance,D_support]^T(T为步态时间)。

步骤二、数据预处理：

将步骤一中采集到的步态时间参数数据进行中值滤波，即将所选择的参数变量用y＝medfilt1(x,n)将n阶一维中值滤波器应用到特征矢量x上，获取滤波后的步态时间参数序列。

步骤三、非线性步态动力学特征提取：

对步态时间参数序列进行非线性步态动力学特征提取，提取每个序列对应的模糊熵值、lz复杂度值以及C0复杂度值。

其中，模糊熵值衡量的是步态时间参数序列新模式产生的概率大小，测度值越大，新模式产生的概率越大，即序列复杂度越大。该特征的提取方法如下：

对于给定的N维步态时间参数序列[u(1),u(2),…,u(N)]，

定义相空间维度m和相似容度r，其中m≤N-2；重构相空间：

x(i)＝[u(i),u(i+1),…,u(i+m-1)]-u₀(i),i＝1,2,…,N-m+1,

其中，参数N＝1,2,3,4…N，

引入模糊隶属函数

其中，对于i＝1,2，…，N-m+1,计算：

为窗口向量X(i)和X(j)之间的最大绝对距离。

针对每个i，求其平均值，得到:

则原步态时间参数序列的模糊熵为：

其中

因此针对有限数据集，模糊熵估计为：

fuzzyEn(m,r,N)＝lnΦ^m(r)-lnΦ^m+1(r)。

Lz复杂度是一种定义在非线性领域度量符号序列复杂性的简单算法，对于步态时间参数序列来说，复杂度代表了随时间增长新模式的产生速率。该特征的提取方法如下：

对于给定的步态时间参数序列x(i)＝{x_i，i＝1，…n}，需要转换成0-1二值符号序列u(i)(i＝1,…,n),其中：

对于重构形成的(0,1)字符串，需要对其进行子串的划分与界定，对于此算法随着步态时间参数序列从前到后计算共有两个基本操作：插入和拷贝。一个长度为n的步态时间参数序列的Lz复杂度c(n)是使用插入和拷贝这两个基本操作来使的原始步态时间参数序列重复构造的最短的字串组合而得来；当序列的长度趋于无限时，几乎所有的符号序列复杂度都趋于一个确定的值，

即：则复杂度c(n)可以通过b(n)进行归一化；即：归一化后的C(n)可以用来对步态时间参数序列的复杂度进行测度。

C0复杂度算法对数据的要求不高，避免了LZ复杂度计算过程中令人担心的粗粒化过程。C0复杂度的算法思想主要是将要分析的步态时间参数序列进行分解，分为随机序列的部分以及规则序列部分，对序列的复杂度预测就被定义为做随机运动部分的步态时间参数序列长度和全部原始步态时间参数序列长度分别与时间轴之间围成的面积之比。该特征的提取方法如下：对于给定的一个长度为N的步态时间参数序列{x(t),t＝0,1,2,…,N-1}进行离散傅里叶变换：

式中i为虚数单位，求出f(k)的均方值：

定义一个大于1的正常数δ，对f(k)序列中超过均方值δ倍的频谱予以保留，反之则将其置零处理，得到转换后的序列：

对上式进行傅里叶逆变换得到：

其中根据C0复杂度算法定义，便可得到C0复杂度计算值：

步骤四、分类识别：

将步骤三提取到的15维特征左脚跨步时长L_stride提取到的模糊熵值、lz复杂度、C0复杂度，右脚跨步时长R_stride提取到的模糊熵值、lz复杂度、C0复杂度,左脚站立时长L_stance提取到的模糊熵值、lz复杂度、C0复杂度，右脚站立时长R_stance提取到的模糊熵值、lz复杂度、C0复杂度,双脚支撑时长D_support提取到的模糊熵值、lz复杂度、C0复杂度，从其中选择特征组合，每个组合均采用支持向量机(svm)进行训练模型，预测测试数据集中的帕金森患者和正常对照者。本发明采用的是线性核函数：K(x,x_i)＝(x^Tx_i)。最终得到最佳识别率对应的特征组合为:左脚跨步间隔模糊熵、右脚站立间隔模糊熵、双脚支撑间隔C0复杂度、左脚站立间隔lz复杂度、双脚支撑间隔lz复杂度这5列步态非线性动力学特征。

步骤五、绘制非线性步态动力学判别图：

通过支持向量机的训练模型可得到支持向量、偏制项b、和拉格朗日乘子α，再通过支持向量和和拉格朗日乘子α求得权值w，因此得到支持向量机的最优分界面wx+b＝0，通过求得的分界线可求得分界点x＝-b/w。绘制一个圆，并把圆周5等分；连接圆心和各等分点，并标以适当的刻度，将这5条半径定义为左脚跨步间隔模糊熵、右脚站立间隔模糊熵、双脚支撑间隔C0复杂度、左脚站立间隔lz复杂度、双脚支撑间隔lz复杂度这5个步态非线性动力学特征的坐标轴；将所求特征分界点在坐标轴上标识并连接成五边形，形成非线性步态动力学特征分界面。对于未知受试者的测量数据，将上述5个特征值分别取在对应的坐标轴上，连接相邻特征值得到的图形与特征分界面对比，即可以直观简易地对帕金森步态进行判别。

本发明与现有技术相比，具有如下优点和有益效果：

1、与传统的基于影像信号的帕金森检测方法相比，本发明通过鞋底压力传感器获取步态时间参数特征，利用特征提取及智能分类识别方法可以方便简单、非侵入的实现帕金森疾病风险评估。

2、与传统的基于步态信号的帕金森分类识别方法相比，本发明方法对帕金森患者的非线性步态动力学特征进行提取，利用帕金森患者和健康正常个体的步态模式之间在非线性动力学上的差异实现帕金森风险评估。

3、本发明方法以特征分布判别图的形式直观呈现了非线性步态动力学特征，正常对照组和帕金森患者不同类别数据对应的特征图有着不同的图形大小，临床医生可以通过判别图直观地观察特征图形的具体大小，大致估计患者所处的帕金森病症严重程度。如果非线性步态动力学判别图的多个指标都超出判别图给出的临界范围，则说明在该受试者中，其患病的风险或者病变严重的程度较高，需要引起个体自身的密切关注以及医生的进一步检查。

附图说明

图1是本发明实施例的帕金森步态风险评估的非线性步态动力学判别图方法的流程图。

图2a实施例中所用PD患者左脚跨步时长特征示意图。

图2b实施例中所用PD患者右脚跨步时长特征示意图。

图2c实施例中所用PD患者左脚站立时长特征示意图。

图2d实施例中所用PD患者右脚站立时长特征示意图。

图2e实施例中所用PD患者双脚支撑时长特征示意图。

图3是本发明实施例中帕金森患者非线性动力学特征数据的分界点。

图4是本发明实施例中帕金森步态风险评估的非线性动力学判别图。

图5是本发明实施例中某一帕金森患者的非线性动力学特征和判别图的综合显示图。

图6是本发明实施例中某一健康个体的非线性动力学特征和判别图的综合显示图。

具体实施方式

下面结合实施例及附图对本发明的具体实施方式作进一步详细的描述，但本发明的实施和保护范围不限于此。

实施例

如图1所示，为本发明实施例的用于帕金森步态风险评估判别的流程图，包括以下步骤：

步骤1获取帕金森患者和对照组的步态步态时间参数序列数据：本发明为使研究结果真实可靠，验证其较好的效果，数据集是由(http://www.physionet.org/physiobank/database/gaitndd)提供的步态动力学数据库。其步态信号是通过放置在试验者鞋底的压力传感器来获取的。目前该数据库提供如下的步态参数：经过的时间(s)、左脚跨步时长(s)、右脚跨步时长(s)、左脚摆动时长(s)、右脚摆动时长(s)、左脚摆动时长(步幅百分比)、右脚摆动时长(步幅百分比)、左脚站立时长(s)、右脚站立时长(s)、左脚站立时长(步幅百分比)、右脚站立时长(步幅百分比)、双脚支撑时长(s)、双脚支撑时长(步幅百分比)。本发明选择的是左脚跨步时长(如图2a)、右脚跨步时长(如图2b)、左脚站立时长(如图2c)、右脚站立时长(如图2d)、双脚支撑时长(如图2e)这5列步态时间参数序列步态参数。

步骤2数据预处理：实验设计是受试者要在77米的跑道上走5分钟，所以必然要进行来回走，期间会有多次转身。为了消除在测试过程中由于受试者由静止开始起步而产生的不稳定数据，要剔除掉起步后前20秒的数据，另外对得到的全过程的数据进行中值滤波，以剔除一些不合理的极大干扰值。其计算方法如下：y＝medfilt1(x,n)将n阶一维中值滤波器应用到矢量x上。该函数将终点以外的信号考虑为0。输出y与x有相同的长度。对n为奇数的情况，y(k)是x(k-(n-1)/2:k+(n-1)/2)的中值；对n为偶数的情况，y(k)是x(k-n/2),x(k-(n/2)+1),...,x(k+(n/2)-1)的中值。

步骤3特征提取：将预处理的数据运用模糊熵、lz复杂度、C0复杂度等算法求出已选择步态时间参数序列的特征。对于给给的N维步态步态时间参数序列[u(1),u(2),…,u(N)]，求模糊熵时，定义相空间维度m为2和相似容度r为0.2倍步态时间参数序列标注差，重构相空间：x(i)＝[u(i),u(i+1),…,u(i+m-1)]-u₀(i),i＝1,2,…,N-m+1,其中，

引入模糊隶属函数

其中，r为相似容限度。对于i＝1,2，…，N-m+1,计算

j＝1,2,…,n＝N-m+1,and j≠i，

为窗口向量X(i)和X(j)之间的最大绝对距离。针对每个i，求其平均值，得到定义原步态时间参数序列的模糊熵为针对有限数据集，模糊熵估计为fuzzyEn(m,r,N)＝lnΦ^m(r)-lnΦ^m+1(r)。Lz复杂度特征的提取方法如下：对于给定的步态步态时间参数序列x(i)＝{x_i，i＝1，…n}需要转换成0-1二值符号序列u(i)(i＝1,…,n),其中其中对于重构的形成的(0,1)字符串之后，需要对其进行子串的划分与界定，对于此算法随着步态时间参数序列从前到后计算共有两个基本操作：插入和拷贝。一个长度为n的步态时间参数序列的Lz复杂度c(n)是使用插入和拷贝这两个基本操作来使的原始序列重复构造的最短的字串组合而得来；当序列的长度趋于无限时，几乎所有的符号序列复杂度都趋于一个确定的值，即：则复杂度c(n)可以通过b(n)进行归一化，即：归一化后的C(n)可以用来对步态时间参数序列的复杂度进行测度。将步态时间参数序列输入到lz复杂度算法可求得步态步态时间参数序列的lz复杂度特征。C0复杂度特征的提取方法如下：对于给定的一个长度为N的步态时间参数序列{x(t),t＝0,1,2,…,N-1}进行离散傅里叶变换：式中i为虚数单位，求出f(k)的均方值：定义一个大于1的正常数δ，对f(k)序列中超过均方值δ倍的频谱予以保留，反之则将其置零处理，得到转换后的序列：对上式进行傅里叶逆变换得到：其中根据C0复杂度算法定义，便可得到C0复杂度计算值：将步态时间参数序列输入到C0复杂度算法可求得步态时间参数序列的C0复杂度特征。

步骤4分类识别：将得到的15个特征进行自由组合，并采用具有小样本、非线性等优点的SVM分类算法进行分类实验，SVM算法采用的是线性核函数K(x,x_i)＝(x^Tx_i)。通过特征组合比较最终左脚跨步间隔模糊熵、右脚站立间隔模糊熵、双脚支撑间隔C0复杂度、左脚站立间隔lz复杂度、双脚支撑间隔lz复杂度融合效果最佳。特征融合前须归一化处理，所采用的方法是min-max方法。其计算方法如下：假设对于序列x₁,x₂,···,x_n进行变换：则新序列为y₁,y₂,···,y_n∈[0,1]，融合后特征的识别率达到了92.8571％。

步骤5绘制非线性步态动力学特征判别图：将左脚跨步间隔模糊熵、右脚站立间隔模糊熵、双脚支撑间隔C0复杂度、左脚站立间隔lz复杂度、双脚支撑间隔lz复杂度分别通过支持向量机进行识别，通过得到的的训练模型可得到权值w、拉格朗日乘子α和偏置项b，因此可训练出分界线wx+b＝0，由此可得每个单特征的分界点为x＝-b/w。对单个特征进行SVM分类计算，将获得的分界点映射到二维图形的各属性坐标轴。如图3所示，这里将31例数据(17例数据训练，14例数据测试)对每个特征通过SVM训练出分界线w.x+b＝0并计算相应的分界点。绘制一个圆，并把圆周5等分；连接圆心和各等分点，并标以适当的刻度，将这5条半径定义为各属性的坐标轴；对测量数据的5个属性值分别取在对应的坐标轴上，然后连接相邻的点形成一个5边形。对于步态动力学的5个非线性步态动力学特征，连接各个分界点组成一个多边形分界面，以此分界面和帕金森患者构成的多边形、正常人构成的多边形对比为医生提供了直观的判别方式。如图4所示，我们认为在分界面范围之内的可以很直观的判定为帕金森病低危人群，完全在分界面范围之外被判定为帕金森病高危人群。若受试者有某一项或某几项非线性动力学指标超出了分界点，则被评估为高危须引导自身个体的密切注意以及进一步的检查。图5是选取的帕金森患者的对应的非线性动力学特征图和判别图的综合显示，图6是选取的正常对照组对应的非线性动力学特征图和判别图的综合显示，可以明显发现，帕金森患者的面积远远大于正常对照组的面积。

尽管已经示出和描述了本发明的实施例，但本发明的实施方式并不受上述实施例的限制，其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化，均应为等效的置换方式，都包含在本发明的保护范围之内。

Claims (6)

1.用于帕金森步态风险评估的非线性步态动力学判别方法，其特征在于，包括如下步骤：

步骤一、步态时间参数获取：提取五个不同的步态时间参数，构成一组步态参数变量；

利用放置在试验者鞋底的力敏薄膜开关来获取人体步态时间参数序列数据，具体包括左脚跨步时长L_stride,右脚跨步时长R_stride,左脚站立时长L_stance，右脚站立时长R_stance,双脚支撑时长D_support，构成一组特征矢量x：x＝[x₁,x₂,x₃,x₄,x₅]^T＝[L_stride.R_stride.L_stance,R_stance,D_support]^T，其中T为步态时间；

步骤二、数据预处理：对步骤一提取出来的步态时间参数序列据进行中值滤波处理以去除极大值干扰；即将所选择的参数变量用y＝medfilt1(x,n)将n阶一维中值滤波器应用到特征矢量x上，获取滤波后的步态时间参数序列；

步骤三、非线性步态动力学特征提取：

对每个步态时间参数序列进行非线性步态动力学特征提取，提取每个序列对应的模糊熵值、lz复杂度值以及C0复杂度值；

步骤四、分类识别：将步骤三提取到的15个非线性步态动力学特征进行特征组合，对正常健康对照组和帕金森病患者采用支持向量机(SVM)训练分类；根据训练模式与测试模式之间关于非线性动力学指标间的差异，把测试模式准确分类识别出来，实现异常步态检测；通过特征组合比较，最终确定分类识别效果最佳的一组特征组合，用于后续非线性动力学判别图的绘制；

步骤五、绘制非线性步态动力学判别图：由支持向量机训练集得到的支持向量机模型求出最佳特征组合中每个非线性步态动力学指标最优超平面wx+b＝0；结合最优超平面的分界点和非线性步态动力学特征参数绘制非线性步态动力学判别图。

2.根据权利要求1所述的用于帕金森步态风险评估的非线性步态动力学判别方法，其特征在于步骤三所述的模糊熵值提取方法，具体：

对于给定的N维步态时间参数序列[u(1),u(2),…,u(N)]，

定义相空间维度m和相似容度r，其中m≤N-2；重构相空间：

x(i)＝[u(i),u(i+1),…,u(i+m-1)]-u₀(i),i＝1,2,…,N-m+1,

其中，参数N＝1,2,3,4…N，

引入模糊隶属函数

其中，对于i＝1,2，…，N-m+1,计算：

为窗口向量X(i)和X(j)之间的最大绝对距离；

针对每个i，求其平均值，得到:

则原步态时间参数序列的模糊熵为：

其中

因此针对有限数据集，模糊熵估计为：

fuzzyEn(m,r,N)＝lnΦ^m(r)-lnΦ^m+1(r)。

3.根据权利要求2所述的用于帕金森步态风险评估的非线性步态动力学判别方法，其特征在于步骤三所述lz复杂度值提取具体如下：

对于给定的步态时间参数序列x(i)＝{x_i，i＝1，…n}，需要转换成0-1二值符号序列u(i)(i＝1,…,n),其中：

对于重构形成的(0,1)字符串，需要对其进行子串的划分与界定，对于此算法随着步态时间参数序列从前到后计算共有两个基本操作：插入和拷贝；一个长度为n的步态时间参数序列的Lz复杂度c(n)是使用插入和拷贝这两个基本操作来使的原始步态时间参数序列重复构造的最短的字串组合而得来；当序列的长度趋于无限时，几乎所有的符号序列复杂度都趋于一个确定的值，

即：则复杂度c(n)通过b(n)进行归一化；即：归一化后的C(n)可以用来对步态时间参数序列的复杂度进行测度。

4.根据权利要求3所述的用于帕金森步态风险评估的非线性步态动力学判别方法，其特征在于步骤三所述C0复杂度提取具体实现如下：

对于给定的一个长度为N的步态时间参数序列{x(t),t＝0,1,2,…,N-1}进行离散傅里叶变换：

式中i为虚数单位，求出f(k)的均方值：

定义一个大于1的正常数δ，对f(k)序列中超过均方值δ倍的频谱予以保留，反之则将其置零处理，得到转换后的序列：

对上式进行傅里叶逆变换得到：

其中根据C0复杂度算法定义，便可得到C0复杂度计算值：

5.根据权利要求4所述的用于帕金森步态风险评估的非线性步态动力学判别方法，其特征在于步骤四所述的支持向量机分类识别方法是将得到的15个特征进行自由组合，并采用SVM分类算法进行分类实验，SVM算法采用的是线性核函数K(x,x_i)＝(x^Tx_i)；最终选择结果为左脚跨步间隔模糊熵、右脚站立间隔模糊熵、双脚支撑间隔C0复杂度、左脚站立间隔lz复杂度、双脚支撑间隔lz复杂度这5列步态非线性动力学特征。

6.根据权利要求5所述的用于帕金森步态风险评估的非线性步态动力学判别方法，其特征在于步骤五具体如下：

通过支持向量机的训练模型得到支持向量、偏制项b、和拉格朗日乘子α，再通过支持向量和和拉格朗日乘子α求得权值w，因此得到支持向量机的最优分界面wx+b＝0，通过求得的分界线可求得分界点x＝-b/w；绘制一个圆，并把圆周5等分；连接圆心和各等分点，并标以设定的刻度，将这5条半径定义为左脚跨步间隔模糊熵、右脚站立间隔模糊熵、双脚支撑间隔C0复杂度、左脚站立间隔lz复杂度、双脚支撑间隔lz复杂度这5个步态非线性动力学特征的坐标轴；将所求特征分界点在坐标轴上标识并连接成五边形，形成非线性步态动力学特征分界面；对于未知受试者的测量数据，将上述5个特征值分别取在对应的坐标轴上，连接相邻特征值得到的图形与特征分界面对比，从而直观简易地对帕金森步态进行判别。