CN101816822A - 功能性电刺激pid参数双源特征融合微粒子群整定方法 - Google Patents

Abstract

本发明涉及以电脉冲刺激进行肢体康复的器械领域。提供一种功能性电刺激中PID参数的双源特征融合混沌微粒子群整定方法，能够准确稳定实时地控制FES系统的电流强度，有效地提高FES系统准确性和稳定性。本发明采用的技术方案是，首先，利用助行过程的柄反作用矢量HRV预测膝关节角度；其次，利用混沌微粒群算法整定比例微积分PID参数，实时调控FES电流水平强度，最终实现比例微积分PID控制参数的自适应在线整定，并用于功能性电刺激FES系统。本发明主要应用于整定功能性电刺激中PID参数。

Description

功能性电刺激PID参数双源特征融合微粒子群整定方法

技术领域

本发明涉及以电脉冲刺激进行肢体康复的器械领域，尤其是功能性电刺激中PID参数的双源特征融合混沌微粒子群整定方法。

背景技术

功能性电刺激(Functional Electrical Stimulation，FES)是通过电流脉冲序列来刺激肢体运动肌群及其外周神经，有效地恢复或重建截瘫患者的部分运动功能的技术。据统计，由于脊髓再生能力微弱，针对脊髓损伤瘫痪患者，目前尚未有可直接修复损伤的有效医治方法，实施功能康复训练是一有效的措施。脊髓损伤瘫痪患者人数逐年增多，功能康复训练是亟待需求的技术。20世纪60年代，Liberson首次成功地利用电刺激腓神经矫正了偏瘫患者足下垂的步态，开创了功能性电刺激用于运动和感觉功能康复治疗的新途径。目前，FES已经成为了恢复或重建截瘫患者的部分运动功能，是重要的康复治疗手段。然而，如何精密控制FES的触发时序和脉冲电流强度以保证电刺激作用效果能准确完成预定的功能动作仍是FES的技术关键。据统计，目前FES的触发控制的方式研究尚少，而且根据作用效果与预定动作偏差，用闭环控制来自动调整FES刺激强度和时序参数，从而大大提高了FES系统的准确性和稳定性，但是现在有效的控制方法仍然在探索之中。

柄反作用矢量(handle reactions vector，HRV)是根据步行器帮助下的站立及行走的过程中，步行器提供给患者的效用实际上可以分为明确独立的3个部分：前后向的力推进，左右向的力平衡和上下向的力支持，这其实也可理解为患者为维持自身正常站立行走对外界所需的附加力学诉求提出的新概念，即是患者在站立行走过程中对步行器的作用合成简化为集中载荷，分别用手柄中点横截面形心处的两个力学矢量来表示，如图1所示，矢量在x，y，z轴上的方向分量合力大小可以分别表征患者借助步行器所获得的力推进，力平衡和力支持水平。其中，定义坐标系所设定的x轴正向为患者的右向，y轴正向为患者的前向，z轴正向为患者的上向。这样，HRV的定义公式也可以写为：

[HRV]＝[HRV₁，HRV_r]^T＝[F_lx，F_ly，F_lz，F_rx，F_ry，F_rz]^T              (1)

目前，HRV被广泛地应用在监视在电刺激过程中病人行走时的状况，继而防止病人摔倒，造成二次伤害。本专利提出利用此参数预测膝关节角度，继而精密控制FES系统的电流水平强度，保证电刺激作用效果能准确完成预定的功能动作，并且防止肌疲劳。

比例微积分(proportional-integral-differential，PID)是一种非常实用的反馈调节算法，它根据系统检测或操作偏差，利用比例、积分、微分运算获得所需调节量以对系统进行反馈控制，因其操作方便而广泛用于工程实践。尤其当被控系统特性参数不明确或难以及时在线测定时，稳妥的闭环控制即可采用PID整定算法。面对肌肉的复杂性和时变性操作环境，由于PID的稳定性好、工作可靠，目前仍在功能性电刺激领域得到了广泛的应用。PID核心技术是精密确定其中比例、积分、微分系数，尤其在FES领域，对系统稳定性要求极为严格，所以对PID参数选择尤为重要。PID控制要取得较好的控制效果，必须调整好比例、积分和微分三种控制作用，形成控制量中既相互配合又相互制约的关系。

发明内容

为克服现有技术的不足，提供一种功能性电刺激中PID参数的双源特征融合混沌微粒子群整定方法，能够准确稳定实时地控制FES系统地电流强度，有效地提高FES系统准确性和稳定性，并获得可观的社会效益和经济效益。为达到上述目的，本发明采用的技术方案是：功能性电刺激中PID参数的双源特征融合混沌微粒子群整定方法，包括：

首先，利用助行过程的柄反作用矢量HRV预测膝关节角度；

其次，利用混沌微粒群算法整定比例微积分PID参数，实时调控FES电流水平强度，整定流程为：首先根据比例微积分PID的三个决策变量K_p、K_i和K_d取值范围的上下界，确定粒子群群体规模、搜索空间维数等参数，并初始化粒子群体的速度和位置，然后利用通过实际关节角度与肌肉模型输出关节角度的相应关系作为适度评价函数计算粒子群中每一个粒子的适应度值，并将其适应度与本身的最佳位置适应度值作比较，并将其作为粒子代表值，然后在调整粒子的速度及其他参数，改变粒子的最佳位置，直到稳定为止，计算最终最佳的位置即得比例微积分PID的K_p、K_i和K_d三个系数，在新的比例微积分PID系数下计算系统输出yout及其与肌肉模型输出关节角度的偏差后再进入下一步神经网络的自学习与加权系数自调整，反复此过程，最终实现比例微积分PID控制参数的自适应在线整定，并用于功能性电刺激FES系统。

所述肌肉模型输出关节角度是采用偏最小二乘回归的方法，即：

设有m个HRV变量HRV1，…，HRVm，p个M变量，M1，…，Mp，共i(i＝1，…，n)个观测值的数据集，T、U分别为从HRV变量与M变量中提取的成分，称为偏最小二乘因子，

从原始变量集中提取第一对成分T1、U1的线性组合为：

T₁＝ω₁₁HRV₁+…+ω_1mHRV_m＝ω₁′HRV                   (4)

U₁＝v₁₁M₁+…+v_1pM_p＝v₁′M                            (5)

其中ω₁＝(ω₁₁，…，ω_1m )′为模型效应权重，v₁＝(v₁₁，…，v_1p)′为M变量权重，将上述提取第一成分的要求转化为求条件极值问题：

其中t1、u1为由样本求得的第一对成分的得分向量，HRV0、M0为初始变量，利用拉格朗日乘子法，上述问题转化为求单位向量ω₁和v₁，使θ₁＝ω₁′ HRV₀′M₀V₁最大，即求矩阵HRV₀′M₀M₀′HRV₀的特征值和特征向量，其最大特征值为θ₁ ²，相应的单位特征向量就是所求的解ω₁，而v₁由公式

得到；

其次建立初始变量对T1的方程

$\left\{\begin{array}{c}{\mathrm{HRV}}_0=t_1{\mathrm{\α}}_1^{\′}+E_1\\ M_0=t_1{\mathrm{\β}}_1^{\′}+F_1\end{array}\right.---\left(7\right)$

其中t₁意义同前，α₁′＝(α₁₁，…，α_1m)，β₁′＝(β₁₁，…，β_1p)为仅一个M量t1时的参数向量，E1、F1分别为n×m和n×p残差阵，按照普通最小二乘法可求得系数向量α₁和β₁，其中α₁成为模型效应载荷量；

如提取的第一成分不能达到回归模型的精度，运用残差阵E1、F1代替X0、Y0，重复提取成分，依次类推，假设最终提取了r个成分，HRV0、M0对r个成分的回归方程为：

$\left\{\begin{array}{c}{\mathrm{HRV}}_0=t_1{\mathrm{\α}}_1^{\′}+...+t_r{\mathrm{\α}}_r^{\′}+E_r\\ M_0=t_1{\mathrm{\β}}_1^{\′}+...+t_r{{\mathrm{\β}}_r}^{\′}+F_r\end{array}\right.---\left(8\right)$

把第一步分析所得HRV量中提取成分Tk(k＝1，…，r)线性组合带入M量对r个成分建立的回归方程，即t_r＝ω_k1HRV₁+…+ω_kmHRV_m代入M_j＝t₁β_1j+…+t_rβ_rj(j＝1，…，p)，即得标准化变量的回归方程M_j＝α_j1HRV₁+…+α_jmHRV_m；

最后根据公式L＝M×HRV^-1，即可求出L，M表示膝关节角度，HRV表示使用者施加在步行器上力的柄反作用矢量，L表示HRV与M之间的关系。

所述利用混沌微粒群算法整定PID参数，进一步细化为：

比例微积分PID采用比例单元P、积分单元I和微分单元D三部分组成，根据系统的误差，通过设定的K_p、K_i和K_d三个参数对系统进行控制：

$\mathrm{yout}\left(t\right)=K_p\mathrm{error}\left(t\right)+K_i\underset{j=0}{\overset{t}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}\left(t\right)-\mathrm{error}(t-1)]---\left(9\right)$

其中K_p是比例系数，K_i是积分系数，K_d是微分系数，error为预设输出与实际输出的偏差，u(t)为PID的输出，同时又是受控系统的输入，由PID输出公式 $\mathrm{yout}\left(t\right)=K_p\mathrm{error}\left(t\right)+K_i\underset{j=0}{\overset{t}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}\left(t\right)-\mathrm{error}(t-1)]$ 可以得到

$u(t-1)=K_p\mathrm{error}(t-1)+K_i\underset{j=0}{\overset{t-1}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}(t-1)-\mathrm{error}(t-2)]---\left(10\right)$

根据：

Δu(t)＝u(t)-u(t-1)

＝K_p(error(t)-error(t-1))+K_ierror(t)+K_d(error(t)-2error(t-1)+error(t-2))

……………………………………………………………                        (11)

有：

u(t)＝Δu(t)+u(t-1)＝

u(t-1)+K_p(error(t)-error(t-1))+K_ierror(t)+K_d(error(t)-2error(t-1)+error(t-2))

………………(12)

采用混沌微粒群算法进行比例微积分PID控制参数的自适应优化，选择收索空间为3维，即分别为PID控制器的三参数，选取群体规模m＝20，群体的初始速度和位置在一定的空间范围内随机产生，分别表示为：v_i＝(v_i1，v_i2，v_i3)，x_i＝(x_i1，x_i2，x_i3)，记第i个粒子迄今搜索到最优位置为p_i＝(p_i1，p_i2，p_i3)，整个粒子群迄今为止搜索到得最优位置为p_gi＝(p_gi1，p_gi2，p_gi3)，其中，i＝1，2，…，20，粒子群优化算法采用以下公式对粒子群操作，

v_id←v_id+c₁r₁(p_id-x_id)+c₂r₂(p_gid-x_id)+c₃r₃(q_id-x_id)              (13)

x_id←x_id+v_id                                                     (14)

其中，i＝1，2，…，20；学习因子c₁、c₂和c₃是非负数，一般取值为0.5；r₁、r₂和r₃是介于[0，1]之间的随机数，q_id是随机选取粒子的位置；

实现的具体步骤为：

1、确定参数：学习因子c₁、c₂和c₃，和群体的规模N，进化次数以及混沌寻优次数；

2、随机产生N个粒子进行操作；

3、按公式(13)和(14)对粒子进行操作；

4、对最优位置p_gi＝(p_gi1，p_gi2，p_gi3)进行混沌优化，将p_gid(i＝1，2，…，20)，映射到

Logistic方程z_i+1＝μz_i(1-z_i)i＝0，1，2，…的定义域[0，1]；

i＝1，2，…，20，然后，用Logistic方程进行迭代产生混沌变量

m＝1，2，…再把产生的混沌变量序列(m＝1，2，…)通过逆映射

(m＝1，2，…)返回到原解空间，得

(m＝1，2，…)

在原解空间对混沌变量经历的每一个可行解

(m＝1，2，…)计算其适应值，保留性能最好的可行解p^*；

5、随机从当前群体中选出的一个粒子用p^*取代；

6、若达到最大代数或者得到满意解，则优化过程结束，否则返回步骤3。

本发明的特点在于：利用助行器的HRV变化预测膝关节角度变化，然后通过混沌粒子群算法优化PID的比例系数、微分系数以及积分系数，继而控制FES系统的电流脉冲强度，有效地提高了FES系统准确性和稳定性。

附图说明

图1柄反作用矢量(HRV)定义示意图。

图2基于HRV的FES系统结构框图。

图3混沌粒子群算法整定PID参数控制方法的结构框图。

图4助行功能性电刺激中的人体模型。

图5实验场景。

图6混沌粒子群算法整定的PID控制追踪结果。

图7混沌粒子群算法整定PID参数控制下预设输入关角度与实际输出的相对误差。

具体实施方式

基于HRV的功能性电刺激助行中的精密控制新技术的应用的结构如图2所示，其工作流程为：首先，利用助行过程的HRV预测膝关节角度，其次，利用混沌微粒群算法整定PID参数，实时调控FES电流水平强度。其整定结构示意图如图3所示，为：首先根据PID的三个决策变量K_p、K_i和K_d取值范围的上下界，确定粒子群群体规模、搜索空间维数等参数，并初始化粒子群体的速度和位置，然后利用通过实际关节角度与肌肉模型输出关节角度的相应关系作为适度评价函数计算粒子群中每一个粒子的适应度值，并将其适应度与本身的最佳位置适应度值作比较，并将其作为粒子代表值，然后在调整粒子的速度等其他参数，改变粒子的最佳位置，直到稳定为止，计算最终最佳的位置即得PID的K_p、K_i和K_d三个系数。在新的PID系数下计算系统输出yout及其与肌肉模型的偏差后再进入下一步神经网络的自学习与加权系数自调整。反复此过程，最终实现PID控制参数的自适应在线整定，并用于FES系统。

一、HRV预测膝关节角度模型

助行过程中，当使用者在功能性电刺激作用下，抬腿迈步时，为了支持身体稳定，使用者在步行器上所施加的力则有所不同，因为关节的大小不同会使人体重心处于不同位置，则克服重力所施加的力也不同，同时人体所处的平面位置也有所改变，为了位置倾翻则在平面上所施加的力也有所变化，因此，关节角度和使用者对步行器所施加的力有一定的关系，如图4所示。

M＝L·HRV+wPW                        (1)

其中，M表示膝关节角度，HRV表示使用者施加在步行器上力的柄反作用矢量，L表示HRV与M之间的关系，w表示系数，W表示上臂、躯干和下肢的重心，P表示三重心与M之间的关系。

实际中，由于步行器的作用，人体重心移动较小，膝关节角度则可表示成

M＝L·HRV                            (2)

其中，M表示膝关节角度，HRV表示使用者施加在步行器上力的柄反作用矢量，L表示HRV与M之间的关系。根据式2所示，确定L就可以利用HRV取出相应时刻的膝关节角度。

L＝M□HRV^-1                           (3)

本专利求解L时，采用了偏最小二乘回归的方法。

设有m个HRV变量HRV1，…，HRVm，p个M变量，M1，…，Mp，共i(i＝1，…，n)个观测值的数据集。T、U分别为从HRV变量与M变量中提取的成分，这里提取的成分通常称为偏最小二乘因子。

从原始变量集中提取第一对成分T1、U1的线性组合为：

T₁＝ω₁₁HRV₁+…+ω_1mHRV_m＝ω₁′ HRV    (4)

U₁＝v₁₁M₁+…+v_1pM_p＝v₁′M              (5)

其中ω₁＝(ω₁₁，…，ω_1m)′为模型效应权重，v₁＝(v₁₁，…，v_1p)为M变量权重。为保证T1、U1各自尽可能多地提取所在变量组的变异信息，同时保证两者之间的相关程度达到最大，据成分的协方差可由相应成分的得分向量的内积来计算的性质，上述提取第一成分的要求转化为求条件极值问颗。

其中t1、u1为由样本求得的第一对成分的得分向量，HRV0、M0为初始变量。利用拉格朗日乘子法，上述问题转化为求单位向量ω₁和v₁，使θ₁＝ω₁ HRV₀′M₀v₁最大，即求矩阵HRV₀′M₀M₀′HRV₀的特征值和特征向量，其最大特征值为θ₁ ²，相应的单位特征向量就是所求的解ω₁，而v₁由公式

得到。

其次建立初始变量对T1的方程

$\left\{\begin{array}{c}{\mathrm{HRV}}_0=t_1{\mathrm{\α}}_1^{\′}+E_1\\ M_0=t_1{\mathrm{\β}}_1^{\′}+F_1\end{array}\right.---\left(7\right)$

其中t1意义同前，α₁′＝(α₁₁，…，α_1m)，β₁′＝(β₁₁，…，β_1p)为仅一个M量t1时的参数向量，E1、F1分别为n×m和n×p残差阵。按照普通最小二乘法可求得系数向量α₁和β₁，其中α₁成为模型效应载荷量。

如提取的第一成分不能达到回归模型的精度，运用残差阵E1、F1代替X0、Y0，重复提取成分，依次类推。假设最终提取了r个成分，HRV0、M0对r个成分的回归方程为：

$\left\{\begin{array}{c}{\mathrm{HRV}}_0=t_1{\mathrm{\α}}_1^{\′}+...+t_r{\mathrm{\α}}_r^{\′}+E_r\\ M_0=t_1{\mathrm{\β}}_1^{\′}+...+t_r{{\mathrm{\β}}_r}^{\′}+F_r\end{array}\right.---\left(8\right)$

把第一步分析所得HRV量中提取成分Tk(k＝1，…，r)线性组合带入M量对r个成分建立的回归方程，即t_r＝ω_k1HRV₁+…+ω_kmHRV_m代入M_j＝t₁β_1j+…+t_rβ_rj(j＝1，…，p)，即得标准化变量的回归方程M_j＝α_j1HRV₁+…+α_jmHRV_m。

最后根据式3，即可求出L。

二、混沌粒子群算法整定PID参数的控制

PID由比例单元P、积分单元I和微分单元D三部分组成，根据系统的误差，通过设定的K_p、K_i和K_d三个参数对系统进行控制。

$\mathrm{yout}\left(t\right)=K_p\mathrm{error}\left(t\right)+K_i\underset{j=0}{\overset{t}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}\left(t\right)-\mathrm{error}(t-1)]---\left(9\right)$

其中K_p是比例系数，K_i是积分系数，K_d是微分系数，error为预设输出与实际输出的偏差，u(t)为PID的输出，同时又是受控系统的输入。

由PID输出公式(1)可以得到

$u(t-1)=K_p\mathrm{error}(t-1)+K_i\underset{j=0}{\overset{t-1}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}(t-1)-\mathrm{error}(t-2)]---\left(10\right)$

根据：

Δu(t)＝u(t)-u(t-1)

＝K_p(error(t)-error(t-1))+K_ierror(t)+K_d(error(t)-2error(t-1)+error(t-2))

……………………………………………………………(11)

有：

u(t)＝Δu(t)+u(t-1)＝

u(t-1)+K_p(error(t)-error(t-1))+K_ierror(t)+K_d(error(t)-2error(t-1)+error(t-2))

………………(12)

本发明采用混沌微粒子群算法进行PID控制参数的自适应优化，选择收索空间为3维，即分别为PID控制器的三参数，选取群体规模m＝20，群体的初始速度和位置在一定的空间范围内随机产生，分别表示为：v_i＝(v_i1，v_i2，v_i3)，x_i＝(x_i1，x_i2，x_i3)，记第i个粒子迄今搜索到最优位置为p_i＝(p_i1，p_i2，p_i3)，整个粒子群迄今为止搜索到得最优位置为p_gi＝(p_gi1，p_gi2，p_gi3)。其中，i＝1，2，…，20，粒子群优化算法采用以下公式对粒子群操作。

v_id←v_id+c₁r₁(p_id-x_id)+c₂r₂(p_gid-x_id)+c₃r₃(q_id-x_id)                   (13)

x_id←x_id+v_id                                                          (14)

其中，i＝1，2，…，20；学习因子c₁、c₂和c₃是非负数，一般取值为0.5；r₁、r₂和r₃是介于[0，1]之间的随机数，q_id是随机选取粒子的位置。

其实现的具体步骤为：

1、确定参数：学习因子c₁、c₂和c₃，和群体的规模N，进化次数以及混沌寻优次数；

2、随机产生N个粒子进行操作；

3、按公式(13)和(14)对粒子进行操作；

4、对最优位置p_gi＝(p_gi1，p_gi2，p_gi3)进行混沌优化，将p_gid(i＝1，2，…，20)，映射到Logistic方程z_i+1＝μz_i(1-z_i)i＝0，1，2，…的定义域[0，1]；

i＝1，2，…，20，然后，用Logistic方程进行迭代产生混沌变量

m＝1，2，…，再把产生的混沌变量序列

(m＝1，2，…)通过逆映射

(m＝1，2，…)返回到原解空间，得

(m＝1，2，…)

在原解空间对混沌变量经历的每一个可行解

(m＝1，2，…)计算其适应值，保留性能最好的可行解p^*。

5、随机从当前群体中选出的一个粒子用p^*取代。

6、若达到最大代数或者得到满意解，则优化过程结束，否则返回步骤3.

三、实验方案

实验装置采用无线传输的步行器系统和美国SIGMEDICS公司生产的Parastep功能性电刺激系统，该系统包含微处理器和刺激脉冲发生电路，含六条刺激通道，电池供电。实验内容为：利用FES系统对下肢相关肌群进行刺激，使受试者按照预定的动作，同时记录施加在步行器上的HRV首先通过安装在步行器上的12导联应变片(BX3506AA)电桥网络转化成的电压信号以及膝关节角度运动轨迹。要求受试者身体健康，无下肢肌肉、骨骼疾患，无神经疾患及严重心肺疾患。实验时受试者安坐于步行器前，将刺激电极固定于相应的部位，未施加电刺激时，受试者保持轻松。FES实验场景如图5所示。电刺激脉冲序列采用经典的Lilly波形，脉冲频率为25Hz、脉宽150μs，脉冲电流在0～120m范围内可调。实验中，实时地记录HRV以及可通过改变脉冲电流大小来调整刺激强度以改变由刺激产生的膝关节角度。实验前，设定期望的膝关节角度运动轨迹，实验中利用角度测量计实时检测膝关节张角变化。实验数据采样率为128Hz，数据记录时长为60s。

有益效果

混沌微粒子群算法整定PID参数的新算法对FES脉冲电流幅值进行测算和调整，使FES作用所产生的膝关节角度运动贴近预期的运动轨迹。图6为混沌微粒子群算法整定的PID控制追踪结果。图中红线表示预期运动轨迹、蓝线为实际输出关节角度。X轴为时间，Y轴为膝关节运动角度。为更清楚地观察混沌粒子群算法整定PID的控制误差，如图7混沌微粒子群整定PID控制下预设输入膝关节角度与实际膝关节角度的相对误差所示，则可以看出误差均在5％之内，可以达到精确的控制。

本发明的主旨是提出一种新的FES的精密控制方法，利用步行器的HRV参数预测的膝关节角度与实际的膝关节角度预测的关节角度的误差，通过混沌微粒群算法优化PID的比例系数、积分系数以及微分系数，继而准确稳定实时地控制FES系统地电流强度。该项发明可有效地提高FES系统准确性和稳定性，并获得可观的社会效益和经济效益。最佳实施方案拟采用专利转让、技术合作或产品开发。

Claims (3)

1.一种功能性电刺激PID参数双源特征融合微粒子群整定方法，其特征是，包括下列步骤：

首先，利用助行过程的柄反作用矢量HRV预测膝关节角度；

其次，利用混沌微粒群算法整定比例微积分PID参数，实时调控FES电流水平强度，整定流程为：首先根据比例微积分PID的三个决策变量K_p、K_i和K_d取值范围的上下界，确定粒子群群体规模、搜索空间维数等参数，并初始化粒子群体的速度和位置，然后利用通过实际关节角度与肌肉模型输出关节角度的相应关系作为适度评价函数计算粒子群中每一个粒子的适应度值，并将其适应度与本身的最佳位置适应度值作比较，并将其作为粒子代表值，然后在调整粒子的速度及其他参数，改变粒子的最佳位置，直到稳定为止，计算最终最佳的位置即得比例微积分PID的K_p、K_i和K_d三个系数，在新的比例微积分PID系数下计算系统输出yout及其与肌肉模型输出关节角度的偏差后再进入下一步神经网络的自学习与加权系数自调整，反复此过程，最终实现比例微积分PID控制参数的自适应在线整定，并用于功能性电刺激FES系统。

2.根据权利要求1所述的一种功能性电刺激PID参数双源特征融合微粒子群整定方法，其特征是，肌肉模型输出关节角度是采用偏最小二乘回归的方法，即：

设有m个HRV变量HRV1，…，HRVm，p个M变量，M1，…，Mp，共i(i＝1，…，n)个观测值的数据集，T、U分别为从HRV变量与M变量中提取的成分，称为偏最小二乘因子，

从原始变量集中提取第一对成分T1、U1的线性组合为：

T₁＝ω₁₁HRV₁+…+ω_1mHRV_m＝ω′₁HRV    (4)

U₁＝v₁₁M₁+…+v_1pM_p＝v′₁M             (5)

其中ω₁＝(ω₁₁，…，ω_1m)′为模型效应权重，v₁＝(v₁₁，…，v_1p)′为M变量权重，将上述提取第一成分的要求转化为求条件极值问题：

其中t1、u1为由样本求得的第一对成分的得分向量，HRVO、MO为初始变量，利用拉格朗日乘子法，上述问题转化为求单位向量ω₁和v₁，使θ₁＝ω′₁HRV′₀M₀v₁最大，即求矩阵HRV′₀M₀M′₀HRV₀的特征值和特征向量，其最大特征值为θ₁ ²，相应的单位特征向量就是所求的解ω₁，而v₁由公式

得到；

其次建立初始变量对T1的方程

$\left\{\begin{array}{c}{\mathrm{HRV}}_0=t_1{\mathrm{\α}}_1^{\′}+E_1\\ M_0=t_1{\mathrm{\β}}_1^{\′}+F_1\end{array}\right.---\left(7\right)$

其中t1意义同前，α′₁＝(α₁₁，…，α_1m)，β′₁＝(β₁₁，…，β_1p)为仅一个M量t1时的参数向量，E1、F1分别为n×m和n×p残差阵，按照普通最小二乘法可求得系数向量α₁和β₁，其中α₁成为模型效应载荷量；

如提取的第一成分不能达到回归模型的精度，运用残差阵E1、F1代替XO、YO，重复提取成分，依次类推，假设最终提取了r个成分，HRVO、M0对r个成分的回归方程为：

$\left\{\begin{array}{c}{\mathrm{HRV}}_0=t_1{\mathrm{\α}}_1^{\′}+...+t_r{\mathrm{\α}}_r^{\′}+E_r\\ M_0=t_1{\mathrm{\β}}_1^{\′}+...+t_r{{\mathrm{\β}}_r}^{\′}+F_r\end{array}\right.---\left(8\right)$

把第一步分析所得HRV量中提取成分Tk(k＝1，…，r)线性组合带入M量对r个成分建立的回归方程，即t_r＝ω_k1HRV₁+…+ω_kmHRV_m代入M_j＝t₁β_1j+…+t_rβ_rj(j＝1，…，p)，即得标准化变量的回归方程M_j＝α_j1HRV₁+…+α_jmHRV_m；

最后根据公式L＝M×HRV^-1，即可求出L，M表示膝关节角度，HRV表示使用者施加在步行器上力的柄反作用矢量，L表示HRV与M之间的关系。

3.根据权利要求1所述的一种功能性电刺激PID参数双源特征融合微粒子群整定方法，其特征是，利用混沌微粒群算法整定PID参数，进一步细化为：

比例微积分PID采用比例单元P、积分单元I和微分单元D三部分组成，根据系统的误差，通过设定的K_p、K_i和K_d三个参数对系统进行控制：

$\mathrm{yout}\left(t\right)=K_p\mathrm{error}\left(t\right)+K_i\underset{j=0}{\overset{t}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}\left(t\right)-\mathrm{error}(t-1)]---\left(9\right)$

其中K_p是比例系数，K_i是积分系数，K_d是微分系数，error为预设输出与实际输出的偏差，u(t)为PID的输出，同时又是受控系统的输入，由PID输出公式

可以得到

$u(t-1)=K_p\mathrm{error}(t-1)+K_i\underset{j=0}{\overset{t-1}{\mathrm{\Σ}}}\mathrm{error}\left(j\right)+K_d[\mathrm{error}(t-1)-\mathrm{error}(t-2)]---\left(10\right)$

根据：

Δu(t)＝u(t)-u(t-1)

＝K_p(error(t)-error(t-1))+K_ierror(t)+K_d(error(t)-2error(t-1)+error(t-2))

……………………………………………………………(11)

有：

u(t)＝Δu(t)+u(t-1)＝

u(t-1)+K_p(error(t)-error(t-1))+K_ierror(t)+K_d(error(t)-2error(t-1)+error(t-2))

                                  ………………(12)

采用混沌微粒群算法进行比例微积分PID控制参数的自适应优化，选择收索空间为3维，即分别为PID控制器的三参数，选取群体规模m＝20，群体的初始速度和位置在一定的空间范围内随机产生，分别表示为：v_i＝(v_i1，v_i2，v_i3)，x_i＝(x_i1，x_i2，x_i3)，记第i个粒子迄今搜索到最优位置为p_i＝(p_i1，p_i2，p_i3)，整个粒子群迄今为止搜索到得最优位置为p_gi＝(p_gi1，p_gi2，p_gi3)，其中，i＝1，2，…，20，粒子群优化算法采用以下公式对粒子群操作，

v_id←v_id+c₁r₁(p_id-x_id)+c₂r₂(p_gid-x_id)+c₃r₃(q_id-x_id)    (13)

x_id←x_id+v_id                                           (14)

其中，i＝1，2，…，20；学习因子c₁、c₂和c₃是非负数，一般取值为0.5；r₁、r₂和r₃是介于[0，1]之间的随机数，q_id是随机选取粒子的位置；

实现的具体步骤为：

1、确定参数：学习因子c₁、c₂和c₃，和群体的规模N，进化次数以及混沌寻优次数；

2、随机产生N个粒子进行操作；

3、按公式(13)和(14)对粒子进行操作；

4、对最优位置p_gi＝(p_gi1，p_gi2，p_gi3)进行混沌优化，将p_gid(i＝1，2，…，20)，映射到Logistic方程z_i+1＝μz_i(1-z_i)i＝0，1，2，…的定义域[0，1]；i＝1，2，…，20，然后，用Logistic方程进行迭代产生混沌变量

m＝1.2.…，再把产生的混沌变量序列(m＝1，2，…)通过逆映射(m＝1，2，…)返回到原解空间，得

(m＝1，2，…)

在原解空间对混沌变量经历的每一个可行解(m＝1，2，…)计算其适应值，保留性能最好的可行解p^*；

5、随机从当前群体中选出的一个粒子用p^*取代；

6、若达到最大代数或者得到满意解，则优化过程结束，否则返回步骤3。

Images