自适应式跑步机是实现虚拟现实中自然行走交互的重要设备之一[1-3]。相比传统的定速跑步机,自适应式跑步机可以根据用户的行走意图变化而主动地改变跑步机速度,使得用户的步态更接近于地面行走步态。一般来说,用户行走意图主要包括开始行走、停止行走及行走速度3方面。其中,行走速度精确估计是保证自适应式跑步机速度跟随功能和性能的前提条件。
光学运动捕捉系统是精确测量用户行走速度的黄金标准,通过在脚部标记反射光点的方式,光学运动捕捉系统可以准确地测量用户的步幅和步频,进而计算出行走速度。但该系统的测量布置过程较为繁琐,而且价格高昂。而廉价的RGB-D相机如微软公司的Kinect尚存在着测量精度和可靠性的争议[4-6]。此外,利用测力板测量用户行走时前后方向的地面反力,结合积分环节和Kalman滤波的方法也可以有效地估计用户行走速度,不足之处是难以长时间连续地测量行走速度,这是因为布置测力板的空间有限。尽管通过将测力板嵌入到跑步机中能够克服该缺陷,但会导致此类跑步机的造价提高[7]。相比之下,惯性测量单元(inertial measurement unit, IMU)和足底压力传感器由于高性价比及良好的便携性能而受到了相关领域研究者的广泛关注。足底压力中心(center of pressure, CoP)的轨迹变化反映了站立相的过程中用户质心的动态变化,因此有学者提出了基于CoP移动速度的行走速度估计方法,并发现CoP移动速度与行走速度之间呈线性正相关[8-10]。但是,对于不同的受试者而言,CoP速度与行走速度之间的关系表达式的系数存在一定的差别,意味着此方法的普适性一般。此外,对于不同年龄段的人群而言,CoP移动速度与行走速度之间的相关性强弱也会发生变化,从而影响速度估计的精度。近几年,神经网络和深度学习的方法也被用于估计行走速度[11-12],尽管该类方法对于训练集内的数据有很高的拟合程度,但对不在训练集内的步态数据的估计精度较差。
基于IMU加速度的双重积分行走速度估计方法则克服了上述问题,因而得到了广泛应用[13-17]。该方法先将IMU固定在用户小腿、脚踝或者脚跟位置,以测量用户行走时的加速度;之后利用得到的角速度来识别脚掌贴合地面(foot-flat)事件,并划分步态周期;结合零速度更新法在一个步态周期内对加速度进行二次积分得到步幅;最终的行走速度由步幅除以步态周期得到。然而,由于IMU安装误差及积分漂移等问题[18-21]的存在,使得这种方法在慢速行走时会出现欠估计,而在快速行走时则会出现过估计。这意味着现有的零速度更新法在解决积分漂移问题时仍存在不足之处,尽管学者Li等[14]提出的增加积分补偿项和步幅补偿项的方法能够在一定程度上解决积分漂移问题,但是补偿项的参数大小还是需要根据速度估计结果人为进行调整,这无疑增加了速度估计的复杂度。加速度二次积分的另外一个问题是噪声累积误差,学者Krzysztof[20]通过引入信息熵和Hilbert变换的方法,对积分过程产生的噪声累积误差进行估计,并将其从预测的行走速度当中剔除,得到了更小的均方根误差和更高的信噪比。但该方法需要复杂的数学计算过程,而且其结果的好坏受步态周期划分准确度的影响较大。Zhang等[21]则提出了一种自适应零速度更新法以克服累积误差问题,但该方法需要一个额外的胸部加速度计。
为克服上述不足,本文提出了一种基于经验模态分解法和摆动脚参数的行走速度估计方法,以解决双重积分估计行走速度方法(文后统称为传统方法)中存在的积分漂移和累积误差问题。首先利用经验模态分解(empirical mode decomposition, EMD)方法去除加速度一次积分后的速度信号的漂移分量,并选择与步频相近的信号分量重构速度信号;其次,针对重构的速度信号积分过程产生的累积误差而导致的EMD模态混叠问题,采用了集成经验模态分解(ensemble EMD, EEMD)方法对二次积分得到的位移信号进行去噪,以获取准确的摆动脚位移。最后,依据该信号来划分脚尖离地和脚跟着地的步态事件,同时结合摆动相的位移计算出摆动脚的平均速度,通过摆动脚的速度估计进行拟合,获得最终的速度估计结果。
1 理论分析 1.1 基于摆动脚参数的行走速度估计原理根据Yoon等[22]提出的理论,利用摆动脚的特征来估计行走速度需要作如下近似:
1) 脚部摆动的速度可以用正弦曲线拟合;
2) 在步态摆动相,脚部移动的距离是用户腰部位移的2倍。
由近似1可以得到摆动相内脚部摆动速度的表达式:
$ V_{\mathrm{swg}}=V_{\mathrm{mag}} \sin \left(\frac{\rm{\mathtt{π} }}{T_{\mathrm{sw}}} t\right), 0 \leqslant t \leqslant T_{\mathrm{sw}} . $ | (1) |
其中:Vswg是脚部摆动速度,Vmag是脚部摆动速度的最大值,Tsw是摆动相时间。
在摆动相内对式(1)进行定积分,可以得到用户步长(step length)的表达式:
$ S_{\mathrm{sw}}=\int_{0}^{T_{\mathrm{sw}}} V_{\mathrm{swg}} \mathrm{d} t=2 V_{\mathrm{mag}} \frac{T_{\mathrm{sw}}}{\mathtt{π}}. $ | (2) |
进一步地,摆动相内的平均速度可以由步长除以摆动相的时间得到:
$ V_{\mathrm{sw}_{-} \mathrm{avg}}=\frac{S_{\mathrm{sw}}}{T_{\mathrm{sw}}}=\frac{2 V_{\mathrm{mag}}}{\mathtt{π}} . $ | (3) |
此外,由近似2(见图 1)可得用户腰部位置的平均速度为
$ V_{\text {pelv_avg }}=\frac{V_{\text {sw_avg }}}{2}=\frac{V_{\text {mag }}}{\mathtt{π}} \text {. } $ | (4) |
![]() |
图 1 行走时用户站立脚、摆动脚、腰部骨盆之间的相对位置 |
式(4)表明在一个步长范围内,用户腰部位置的平均速度可以由摆动脚的平均速度或其最大速度来计算得到。当用户在跑步机上匀速行走时,其相对于地面的位置变化很小,因而可以认为Vpelv_avg与跑步机速度Vtreadmill大小相等但方向相反,则有
$ V_{\text {treadmill }}=-V_{\text {pelv_avg }} . $ | (5) |
另外,相比在地面行走,在跑步机上行走过程中测量得到的用户摆动脚的最大速度Vmag_t和平均速度Vsw_avg_t也受到Vtreadmill的影响,因而有
$ \left\{\begin{array}{l} V_{\text {mag_t }}=V_{\text {mag }}+V_{\text {treadmill }}, \\ V_{\text {sw_avg_t }}=V_{\text {sw_avg }}+V_{\text {treadmill }} . \end{array}\right. $ | (6) |
将式(4)、(5)代入(6)中可以得到
$ V_{\text {pelv_avg }}=\frac{V_{\operatorname{mag}_{-} \mathrm{t}}}{\mathtt{π}-1}=V_{\mathrm{sw_avg_{ }} \mathrm{t}} . $ | (7) |
由式(7)可以发现,当用户在跑步机上匀速行走时,可以根据摆动脚的最大或平均速度来估计用户的行走速度。由于用户在慢速和快速行走时其摆动相所占的步态周期比例会显著不同,这使得在不同速度下的估计值与实际值之间仍存在一定的偏差。为此,还需要对预测值和实际值进行曲线拟合,达到修正偏差的效果。需要注意的是,Yoon等[21]是利用摆动脚的最大速度来估计用户行走速度,而本文是采用摆动脚的平均速度与行走速度之间的关系来估计用户行走速度。此外,Yoon等采用的测量方法是基于光学定位系统,而本文的测量方案是基于IMU传感器,两者的测量原理存在着本质上的不同。
1.2 基于EMD方法消除摆动脚速度的积分漂移为获取摆动脚平均速度,直接的方法是将摆动脚的位移除以摆动相的时间。其中,摆动相的时间可以利用角速度的峰值特征划分步态事件得到;而摆动脚的位移则可以由加速度通过二次积分获得。然而,由于积分漂移的原因,上述方法获得的摆动脚位移会出现倾斜效应。为了解决该问题,本文针对二次积分过程引入了EMD方法[23]。EMD方法是一种自适应的数据处理和挖掘方法,尤其适用于非线性、非平稳时间序列的处理,而步态信号也存在固有的非线性特征,因而EMD方法在步态分析领域得到了较为广泛的应用。
EMD方法的基本思想是将一个频率不规则的波分解为多个单一频率的波和残波形式。其中,各个单一频率的波被称为本征模函数(intrinsic mode function, IMF),它们分别包含了原信号的在不同时间尺度的局部特征信号。针对一个时间序列f(t)的EMD筛选(sifting)过程如下:
首先,找到原时间序列f(t)的所有极大值和极小值点,并通过3次样条插值函数拟合,分别形成f(t)的上包络线和下包络线。由此可以计算出上下包络线的均值曲线m(t)。
之后,从原时间序列f(t)中剔除均值曲线m(t)得到
$ h_{1}(t)=f(t)-m(t) . $ | (8) |
验证h1(t)是否满足以下2个条件:
1) 极值点的数目和过零点的数目必须相等或者最多相差1个;
2) 在任意点,上包络线的局部最大值和下包络线的局部最小值之间的平均值等于0。
如果不满足上述条件,则说明h1(t)不是一个IMF,需要将h1(t)作为新时间序列重复式(8)的计算过程,直到其满足上述2个条件;若h1(t)满足上述条件,则定义h1(t)为一个IMF并表示为
$ c_{1}(t)=h_{1}(t) . $ | (9) |
单次的筛选过程至此结束。
在筛选出1个IMF之后,将其从原始信号f(t)中剔除以得到残余信号:
$ r_{1}(t)=f(t)-c_{1}(t) . $ | (10) |
将r1(t)作为新的序列并重复上述单次筛选过程,以得到其他IMF,直至得到的残余信号r(t)是一个单调函数则筛选过程停止。最终,原始信号f(t)被分解为
$ f(t) = \sum\limits_{i = 1}^n {{c_n}} (t) + {r_n}(t) . $ | (11) |
借助EMD方法,原始信号可以被分解成不同频率的IMF和残余信号,其中高频信号认为是噪声信号,接近正常行走步频的信号认为是步态特征的信号,残余信号则反映了变化趋势。可以发现,EMD方法能够直接提取出信号中的残余单调分量,该分量是积分漂移产生的倾斜效应的主要来源。因此,将残余单调分量从原信号当中剔除即可解决积分漂移问题。
1.3 基于EEMD方法克服摆动脚位置的累积误差节1.2的EMD过程是依据信号本身特征自适应分解的,不需要人为判定带宽范围。考虑到用户在不同速度行走时步频会产生相应的变化,相比传统的滤波器设计方法,EMD方法有更强的适应性。然而,当传感器信号被非连续噪声污染后,采用EMD方法会出现模态混叠现象。这种现象会使得IMF分量失去其具有原来单一特征尺度的特性,形成尺度混杂的振荡,进而失去其原有的物理意义。
速度估计中的加速度双重积分无疑会将原始的噪声信号累积放大,造成EMD的模态混叠现象,使得摆动脚位置测量结果不准确。为了克服模态混叠现象,Huang等[24]提出了EEMD方法。EEMD方法的基本原理是将白噪声加入待分解信号,利用白噪声频谱的均匀分布,当信号叠加在遍布整个时频空间且分布一致的白噪声背景上时,不同时间尺度的信号会自动分布到合适的参考尺度上,并且由于零均值噪声的特性,经过多次平均,噪声将相互抵消,集成均值的结果就可作为最终结果。
综上,所提出的行走速度估计方法包括3个步骤:
1) 坐标转换。如图 2所示,将IMU固定在受试者的脚踝上方。此时,在矢状面上,IMU本身的局部坐标系ax-az,与世界坐标系ah-av之间存在一定的角度差φy。其中,ah轴方向与用户前进方向保持一致,而av轴方向与重力加速度方向相反。为了通过双重积分计算脚部位置,需要先将局部坐标系下的加速度信号转换到世界坐标系下。
![]() |
图 2 本文行走速度估计方法步骤说明 |
2) 积分重构。对坐标变换后的加速度ar积分得到速度值vr,再利用EMD方法将信号分解成不同频率的IMF,剔除其中的残余漂移分量,并选择与正常行走步频接近的信号分量重构得到新的速度信号vemd。之后,采用EEMD方法对vemd信号进行分解,剔除其中的高频噪声分量后得到重构的位置信号seemd。
3) 速度估计。根据seemd完成摆动周期划分和摆动脚位移计算,从而计算出摆动脚平均速度,结合节1.1中基于摆动脚的速度估计原理可得到最终的速度估计结果。
需要注意的是,所提出的方法并未完全采用EEMD方法来同时消除积分漂移和累积误差问题。这是由于在一次积分后凸显的主要问题是积分漂移,此时信号的信噪比仍然较高,无须采用EEMD方法。而二次积分后导致的噪声累积误差问题不可忽视,此时采用EEMD方法更加合理。此外,EEMD方法分解过程比EMD方法更复杂,采用“先EMD后EEMD”的方式既能够有针对性地解决上述问题,同时与完全采用EEMD方法相比过程更简便。
2 实验验证为了验证所提方法的有效性,对一名受试者在5种不同行走速度下进行了测试。其中,实际行走速度由跑步机速度得到,5种行走速度分别为2.6、3.2、3.8、4.4、5.0 km/h。IMU固定在受试者脚踝正上方,采样频率为50 Hz。受试者在每种行走速度下保持匀速行走1 min,取第10 s到第50 s之间的行走过程作为样本数据。一个典型的坐标转换后的加速度曲线如图 3a所示, 对其进行梯形积分后得到的速度信号有明显的倾斜效应(skew effects),无法用于分析用户行走过程(见图 3b)。
![]() |
图 3 加速度直接积分产生的漂移现象 |
对图 3b的信号进行EMD,得到如图 4的IMF和残余信号。其中,IMF1虽然保留了部分步态准周期特性,但所含能量低且存在不规则的噪声干扰;而IMF4和IMF5则反映的是信号低频的变化,与步态的准周期特性关系不大。因此,本文将采用IMF2+IMF3重构得到的信号作为新的速度量(见图 5)。
![]() |
图 4 分解后的IMF和残余信号 |
![]() |
图 5 重构的摆动脚速度信号 |
对重构的摆动脚速度信号再次积分可以得到摆动脚的位移,对其进行EMD,并取IMF1作为重构后的摆动脚位移,得到如图 6a所示的结果,可以发现该信号中出现了模态混叠现象。采用EEMD消除模态混叠现象后,得到准周期性的摆动脚位移信号如图 6b所示。其中,蓝色星号对应了摆动脚位移的最小值点,此时用户脚部摆动到了身体最后方,对应了脚尖离地时刻;红色圆圈则对应了摆动脚位移的最大值点,此时用户脚部摆动到了身体最前方,对应了脚跟离地的时刻。根据相邻的蓝色星号和红色圆圈之间的时刻差和幅值差即可计算出摆动脚的摆动周期和位移,进而得到摆动脚的平均速度。
![]() |
图 6 (网络版彩图)有模态混叠与无模态混叠信号对比 |
根据上述方法,可以计算出不同行走速度下的摆动相时间、位移和平均速度如表 1所示。随着行走速度的增加,摆动周期逐渐减小,而摆动脚位移则逐渐增加,这使得摆动脚的平均速度也在逐渐增加。进一步地,对摆动脚平均速度和实际行走速度之间进行二次曲线拟合,以减小部分较大偏离点对估计结果造成的影响(见图 7a),得到表达式如下:
$ V_{\text {est }}=a V_{\text {sw_avg }}^{2}+b V_{\text {sw_avg }}+c . $ | (12) |
实际速度/ (m·s-1) | 摆动相时间/s | 摆动脚位移/m | 摆动脚平均速度/(m·s-1) |
0.72 | 0.63±0.04 | 0.48±0.03 | 0.76±0.06 |
0.89 | 0.58±0.02 | 0.53±0.02 | 0.91±0.03 |
1.06 | 0.53±0.03 | 0.54±0.03 | 1.02±0.08 |
1.22 | 0.48±0.01 | 0.60±0.03 | 1.25±0.05 |
1.39 | 0.46±0.01 | 0.65±0.04 | 1.41±0.08 |
![]() |
图 7 (网络版彩图)曲线拟合后的速度估计结果及不同方法间的结果对比 |
式中:Vest表示根据摆动脚平均速度估计得到的行走速度;a、b、c为拟合系数,分别为-0.41、1.83、-0.41。根据上式,最终的速度估计结果如图 7b所示,所提出的速度估计方法的平均误差为0.000 4 m/s,均方根误差为0.06 m/s。
由图 7b可以看出,传统方法在增益调整前会出现慢速欠估计和快速过估计的问题,手动调整增益后,传统方法的速度估计平均误差为0.002 7 m/s,均方根误差为0.055 m/s。不难发现,本文方法整体的平均误差更小,而均方根误差则略大于传统方法。此外,本文方法的优势是无须人工调整补偿增益值,也不需要精准识别foot-flat事件,比传统方法实现起来更为简洁。
综上,本文方法能够达到较高的速度估计精度,并克服了传统方法的存在的不足。由于EMD方法具有自适应性,因此所提速度估计策略具有较好的泛化能力。
3 结论本文提出了一种基于经验模态分解法和摆动脚参数的行走速度估计方法。利用EMD和EEMD分别解决IMU双重积分过程中存在的积分漂移和累积误差问题,并基于提取的摆动脚特征完成行走速度估计。相比传统方法,所提出的方法无须人工调整补偿增益值即可实现更高的速度估计精度。其次,EMD和EEMD方法能够针对信号本身的特征进行分解,使得所提出的方法具有较好的适应性。因此,该方法在步态分析、康复训练、老年人的健康监督等相关领域也有重要的参考价值。下一步将解决所提出的方法的实时性问题,以保证其能有效地应用于自适应式跑步机的控制系统。
[1] |
SOUMAN J L, ROBUFFO P, GIORDANO, et al. Cyberwalk: Enabling unconstrained omnidirectional walking through virtual environments[J]. ACM Transactions on Applied Perception, 2011, 8(4): 25. |
[2] |
FEASEL J, WHITTON M C, KASSLER L, et al. The integrated virtual environment rehabilitation treadmill system[J]. IEEE Transactions Neural System & Rehabilitation Engineering, 2011, 19(3): 290-297. |
[3] |
CHRISTENSEN R R, HOLLERBACH J M, XU Y, et al. Inertial-force feedback for the tread-port locomotion interface[J]. Presence, 2014, 9(1): 1-14. |
[4] |
FOSTY B, BEN-SADOUN G, SACCO G, et al. Accuracy and reliability of the RGB-D camera for measuring walking speed on a treadmill[J]. Gait & Posture, 2016, 48: 113-119. |
[5] |
XU X, MCGORRY R W, CHOU L S, et al. Accuracy of the Microsoft KinectTM for measuring gait parameters during treadmill walking[J]. Gait & Posture, 2015, 34(2): 1-7. |
[6] |
URSKA P, BRITTANY H, JUDITH E D. Validity and reliability of the Kinect for assessment of standardized transitional movements and balance: Systematic review and translation into practice[J]. Physical Medicine & Rehabilitation Clinics of North America, 2019, 30: 399-422. |
[7] |
RAY N T, KNARR B A, HIGGINSON J S. Walking speed changes in response to novel user-driven treadmill control[J]. Journal of Biomechanics, 2018, 78: 143-149. DOI:10.1016/j.jbiomech.2018.07.035 |
[8] |
CHA M, HAN S, KIM H, et al. User-driven treadmill using walking speed estimated from plantar pressure sensor[J]. Electronics Letters, 2017, 53(8): 524-526. DOI:10.1049/el.2016.4171 |
[9] |
KEIJSERS N L W, STOLWIJK N M, RENZENBRINK G J, et al. Prediction of walking speed using single stance force or pressure measurements in healthy subjects[J]. Gait & Posture, 2016, 43: 93-95. |
[10] |
FUCHIOKA S, IWATA A, HIGUCHI Y, et al. The forward velocity of the center of pressure in the midfoot is a major predictor of gait speed in older adults[J]. International Journal of Gerontology, 2015, 9(2): 119-122. DOI:10.1016/j.ijge.2015.05.010 |
[11] |
JOO S B, OH S E, SIM T, et al. Prediction of gait speed from plantar pressure using artificial neural networks[J]. Expert Systems with Applications, 2014, 41(16): 7398-7405. DOI:10.1016/j.eswa.2014.06.002 |
[12] |
WANG K L, XU J. A speed regression using acceleration data in a deep convolutional neural network[J]. IEEE Access, 2019, 7: 9351-9356. DOI:10.1109/ACCESS.2019.2890967 |
[13] |
MANNINI A, SABATINI A M. Walking speed estimation using foot-mounted inertial sensors: Comparing machine learning and strap-down integration methods[J]. Medical Engineering & Physics, 2014, 36(10): 1312-1321. |
[14] |
LI Q, YOUNG M, NAING V, et al. Walking speed estimation using a shank-mounted inertial measurement unit[J]. Journal of Biomechanics, 2010, 43(8): 1640-1643. DOI:10.1016/j.jbiomech.2010.01.031 |
[15] |
WANG L, SUN Y, LI Q, et al. Estimation of step length and gait asymmetry using wearable inertial sensors[J]. IEEE Sensors Journal, 2018, 18(9): 3844-3851. DOI:10.1109/JSEN.2018.2815700 |
[16] |
吴哲明, 孙振国, 张文增, 等. 基于惯性测量单元旋转的陀螺漂移估计和补偿方法[J]. 清华大学学报(自然科学版), 2014, 54(9): 1143-1147. WU Z M, SUN Z G, ZHANG W Z, et al. Gyroscope bias estimation and compensation by rotation of the inertial measurement unit[J]. Journal of Tsinghua University (Science and Technology), 2014, 54(9): 1143-1147. (in Chinese) |
[17] |
班朝, 任国营, 王斌锐, 等. 基于IMU的机器人姿态自适应EKF测量算法研究[J]. 仪器仪表学报, 2020, 41(2): 33-39. BAN C, REN G Y, WANG B R, et al. Research on self-adaptive EKF algorithm for robot attitude measurement based on IMU[J]. Chinese Journal of Scientific Instrument, 2020, 41(2): 33-39. (in Chinese) |
[18] |
SESSA S, ZECCA M, LIN Z, et al. A methodology for the performance evaluation of inertial measurement units[J]. Journal of Intelligent and Robotic Systems, 2013, 71: 143-157. DOI:10.1007/s10846-012-9772-8 |
[19] |
ANWARY A R, YU H, VASSALLO M. Optimal foot location for placing wearable IMU sensors and automatic feature extraction for gait analysis[J]. IEEE Sensors Journal, 2018, 18(6): 2555-2567. DOI:10.1109/JSEN.2017.2786587 |
[20] |
BRZOSTOWSKI K. Toward the unaided estimation of human walking speed based on sparse modeling[J]. IEEE Transactions on Instrumentation and Measurement, 2018, 67(6): 1389-1398. DOI:10.1109/TIM.2018.2800198 |
[21] |
ZHANG R, YANG H, HFLINGER F, et al. Adaptive zero velocity update based on velocity classification for pedestrian tracking[J]. IEEE Sensors Journal, 2017, 17(7): 2137-2145. DOI:10.1109/JSEN.2017.2665678 |
[22] |
YOON J, PARK H S, DAMIANO D. A novel walking speed estimation scheme and its application to treadmill control for gait rehabilitation[J]. Journal of NeuroEngineering and Rehabilitation, 2012, 9: 62. DOI:10.1186/1743-0003-9-62 |
[23] |
HUANG N E, SHEN Z, LONG S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J]. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, 1998, 454: 903-995. |
[24] |
WU Z H, HUANG N E. Ensemble empirical mode decomposition: A noise-assisted data analysis method[J]. Advances in Adaptive Data Analysis, 2009, 1(1): 1-41. |