刘晨晖 张为民 薛 峰
(同济大学机械与能源工程学院,上海 201804)
近年来,随着生产模式的演进,自动化、智能化和精细化模式已经成为了工业生产中的发展方向。机械臂以其自由度高,工作空间大等特性,广泛用于特征异构的个性化产品装配和物流等生产过程中,为了生产效率的提高发挥着越来越重要的作用。
在不同产品导入试生产的调试阶段,除了应对机械臂运动轨迹进行适用性评估及可能的重新规划外,还需要对现有或重新规划完成的轨迹运动参数进行最优化的微调、配置或整定,使机械臂以更少调优成本和更优的性能去适应新的生产任务。因此,在不更改机械臂结构设计参数的条件下,对其运动参数进行调优,以快速完成生产场景切换,仍是解决生产优化问题的重点有效手段。
目前对机械臂轨迹运动参数优化研究主要有黑盒优化与白盒优化两个方面。现有白盒优化可以通过数学建模,构建出因素—指标的近似数理变化关系,从而解决优化的机理性问题;但每次场景切换都需要重新完成建模,优化实验成本和门槛较高。因而通常将机械臂运动参数优化问题视作黑盒优化模型,以多组因素—指标对为训练样本,得到数据驱动的因素—指标映射模型,进而通过极值化得到最优解。
贝叶斯优化算法(Bayesian optimization, BO)作为一种强大、灵活和求解快速的参数优化工具,广泛地用于黑盒函数的高效全局优化[1],不少学者已将这种方法用于工业场景的优化研究中。 Calandra R等[2]针对仿生双足步行机器人优化问题,通过仿真和实际实验,评估、验证了贝叶斯优化在步态优化中的效果。Roveda L等[3]以执行器末端配备有视觉系统的机器人为研究对象,基于贝叶斯优化方法,提出了一种在有限的实验次数下优化机器人末端执行器位姿的方法,以使其在工业任务中灵活工作。在其他领域,王子涵等[4]将贝叶斯模型引入机床热变形的神经网络模型中,通过快速求解超参数,提高了神经网络的预测性能。
上述研究中的贝叶斯优化均使用了高斯过程概率代理模型,尽管体现出了较好的灵活性和优化效果,但其只能处理数据类型为连续或数值型的待优化参数,不能适用于实际生产中广泛存在的离散型或其他难以数值化的参数优化问题。为解决这一问题,大多数研究通过巧妙地改造高斯过程的核函数或者将离散的变量映射到连续空间[5],从而使得该方法能够支持更多优化问题的解决。因此,关于代理模型替换的研究就不断涌现。以随机森林为代表的基于序列模型的算法配置优化算法(sequential model-based algorithm configuration, SMAC)[6]的提出,解决了贝叶斯模型在高斯回归过程代理下,待优化的参数类型不能为离散型的情况,并在神经网络模型非连续型结构参数优化中得到了成功应用;在工业场景中尚未见文献报道。
综上所述,本研究以ABB IRB 1410机械臂为实验对象,选取了难以连续化、宜给定离散值的机械臂路径行走方案编号作为离散变量,设计了3条有代表性的路径方案;并联合ABB IRB 1410机械臂中开放可调权限的速度、加速度和加加速度等3个连续变量参与优化实验;选择了若干指标并对其进行归一化,通过实验验证了随机森林贝叶斯优化方法在机械臂运动参数优化场景的有效性;与经典的基于高斯过程的贝叶斯优化进行了对比,证明了随机森林贝叶斯优化方法的经济性。
贝叶斯定理可以通过某一现象的先验分布与似然概率(观测样本的概率分布)来推测后验分布,形式上:设P(A)是没有数据支持下,A事件发生的概率,即先验概率;P(A|B)在事件B发生下,A事件发生的概率,即后验概率;P(B|A)是A事件发生后B事件发生的概率分布,即似然函数;P(B)是边际似然概率,有
利用贝叶斯原理,在解决黑盒全局参数优化问题上,可以基于已有的实验结果,依次进行迭代,不断提高黑盒函数的置信度,得到下一组具有较高事件概率的参数组合。因此式(1)可以改写为[7]
式中:f表示为未知的黑盒函数;D1:t={(x1,y1),(x2,y2),···,(xt,yt)}表 示 已 有 的 实 验 结 果 ,xt=(x1,t,x2,t,···,xn,t)T表 示 第t次 实 验 下 的n维 参 数 组 合 ,yt=f(xt)+εt表示在xt参数组合下,得到的实验结果,εt为观测误差;P(D1:t|f)表示在f下的y的似然分布;P(f)表示f的先验概率分布,即对黑盒函数的假设;P(D1:t)表示f的边际似然分布;P(f|D1:t)表示f的后验概率分布,描述已有实验的数据对先验进行修正后的黑盒函数的置信度。
贝叶斯优化的难点就在于如何选取合适的概率代理函数和采集函数,以提高优化的速度,减少不必要的参数组合选取。
概率代理模型是贝叶斯优化框架的一个核心部分,用于处理优化变量与优化目标之间复杂的映射关系。不同的概率代理模型在贝叶斯优化的应用中往往发挥着不同的效果,这里常用的代理模型有高斯过程、随机森林等。
1.1.1 高斯过程
高斯过程是多元高斯分布向无穷维度的扩展,是一系列服从正态分布的随机变量在一定数据集z内的组合,由均值函数m(z)和协方差函数k(z,z′)确定为
实验数据 D1:t中xt=(x1,t,x2,t,···,xn,t)T可以看作定义在样本空间上的实值函数,为n维随机变量。在贝叶斯优化中,需要假设优化函数的先验为高斯过程为[8]
基于式(4),f1:t(x)的高斯分布参数为
对于下一个点xt+1,满足正态分布N(µt+1,Σt+1),x1:t与xt+1的联合分布参数为
因此ft+1(x)的概率密度函数为
ft+1(x)与f1:t(x)的联合概率密度为
进一步得到条件概率密度函数为
将f1:t,t+1(x)、f1:t(x)代入求解,得到
式(15)中: ︿x1:t=x1:t-µ1:t,方差 Σt+1|1:t可以通过计算xi之 间协 方差 得到 ,均 值 µt+1|1:t中 的 µt+1可以 使用µ1:t近似计算。
因此,之后每个未知点的95%置信区间均可以通过公式(17)求解得到为
1.1.2 随机森林
随机森林本质是多个分类与回归树。树的内部结点的右分支是取值为“是”,左分支是取值为“否”。递归地二分每个特征,将输入空间(特征空间)划分为有限个单元,并在这些单元上确定预测的概率分布,即在输入给定的条件下输出条件概率分布。
回归树主要用于处理数值类或者连续型变量;分类树易于处理一些非数值类或者离散型变量。
针对回归树生成,给定数据集 D1:t,将第i个特征输入空间 (xi,1,xi,2,···,xi,t)T,以及对应的输出空间(y1,y2,···,yt)T,划分为M个单元 (R1,R2,···,RM)。每个单元Rm上的固定输出为Cm;I为0-1函数;回归树模型可表示为[9]
针对固定输出为Cm,使用平方误差(mean square error, MSE)表征回归树对于训练数据的预测误差,基于MSE最小原则求解模型损失函数式(19)的极小值为
MSE直接对Cm求导并让导数为0,Rm空间上样本数为Nm,有
对于特征输入空间的划分,采用最小二乘法,递归地将每个区域划分两个子域与,为
通过求解式(24)得
在数据集 D1:t第i个特征组中,遍历所有的xi,j,求解出全局最小值的解集。之后对R1与R2继续划分,直到满足迭代终止条件。
分类树生成原理与回归树相似,但一般采用基尼指数最小原则进行空间划分。假设有K个类别,样本点属于第k类的概率pk,基尼指数定义为[9]
当样本集合 D1:t根据C特征划分为 D1与 D2两部分时,得
基尼指数越大,则表示样本集合不确定性越大。因此在分类树划分空间时候,需要遵循基尼指数最小原则。
在随机森林算法中,需要根据具体数据类型,确定回归树与分类树数目与最大迭代次数,对每个特征进行训练,去拟合黑盒函数f,该过程可以类比高斯回归过程中n个点所构成的多维正态分布。对于新加入的点,在随机森林的每棵树上有一个预测结果,即条件概率分布。每棵树相互独立,预测结果均值与方差为
95%置信区间为
采集函数是根据后验概率分布构造的,由于原始函数f(x)的计算成本较高,因此构建一个更便宜的函数,称为采集函数 α (x),由代理模型确定下一个待评价点。在参数空间X中选择下一个点xt+1将满足
这能同时能够平衡“探索”以及“精进”两种机制,通过采集函数最值化,来选择下一个最具潜力的评估点,并避免不必要的采样。
常用的采集函数有置信边界策略等。置信边界策略包括置信下界(lower confidence bound, LCB)和置信上界(upper confidence bound, UCB)两种,都为均值函数和协方差函数组成,UCB 多应用于望大优化,LCB 多应用于望小优化。
本文针对机械臂的调优实验指标均为望小指标,采用LCB策略,得
式中: µ (x)与 σ (x)分别为均值函数与方差函数,k为迭代递减步长,用于平衡“探索”与“精进”机制,保证调优前期尽可能进行宽度搜索,调优后期进行深度搜索。
第1步,先根据初始化数据:D1:t={(x1,y1),(x2,y2),···,(xt,yt)},xt为初始化实验中,输入参数组合,yt为在对应参数组合下的归一化指标,获得一个概率代理模型的先验分布;
第2步,选择LCB作为采集函数,计算 α (x);
第3步,得到使采集函数最大的xt+1值;
第4步,以xt+1为下一次实验参数组合,得到该组实验下的归一指标yt+1;
第5步,将{xt+1,yt+1}作为下一次迭代的训练值,根据1.1节所述各公式,计算、更新概率代理模型的先验分布;
重复第2步到第5步,直到达到迭代次数预设值,或者出现重复推荐数据xt+a=xt+b。
流程图如图1所示。
图1 调优流程图
为简化子指标数据的输入,需要对数据归一化处理,本次研究选择层次分析法。计算步骤如下[10]:
(1)建立层次结构模型,将决策目标、决策准则和决策对象按它们之间的相互关系分为最高层、中间层和最低层,绘出层次结构图。
(2)构造判断矩阵,不把所有因素放在一起比较,而是两两相互比较,并按表1对其重要性程度评定等级。其中判断矩阵具有如下性质
表1 因素重要性程度评定等级表
(3)求解判断矩阵最大特征根的特征向量,进行层次单排序,并进行一致性检验。
(4)层次总排序及其一致性检验,将最后因素层加权求解到最高层。
本研究以ABB IRB 1410机械臂为调优对象,机械臂基坐标系遵循右手系原则,Y轴正方向与Z轴正方向如图2所示。选择机械臂安装燃料电池极板为工作场景,该场景要求电池极板在X-Y平面安装的位置精度有着较高要求,极板中心点位置度应处于0.05 mm公差带内,以保证后续极板点胶工艺顺利进行;同时希望保证安装过程平稳,冲击较小。
图2 ABB IRB 1 410
2.2.1 调优参数介绍
机械臂的调优参数从量化角度可以分为连续型变量与离散型变量,连续型调优变量为以下3个:
(1)机械臂工具中心点(tool center point, TCP)线性速度(记作:速度v),为机械臂末端工具坐标系原点相对于机械臂的基坐标系原点的线性运动速率,单位为 m m/s。可以通过机械臂的示教器,在程序数据模块进行设置,更改范围为0~2 100 m m/s。
(2)机械臂加速度百分率(记作:加速度a),对机械臂TCP线性速度进行调节,实现机械臂末端的加减速过程,是机械臂系统内部额定值的百分比,单位为%。可以通过机械臂示教器,在运动指令-AccSet中进行设置,更改范围为20%~100%。
(3)机械臂加速度坡度(记作:加加速度a’),用于调节机械臂加速度曲线,改善机械臂末端加减速过程的平滑性,是机械臂系统内部额定值的百分比,单位为%。可以通过机械臂示教器,在运动指令-AccSet中进行设置,更改范围为20%~100%。
因为燃料电池极板在被吸取时即位于机械臂末端,因而选择调整工具坐标系本身的运动学参数可以更直观地对实验结果进行采集和反馈。
离散型调优变量为机械臂末端空间路径方案p。机械臂末端空间路径是指机械臂末端工具坐标系原点在空间的线性运动轨迹。如表2所示,3条路径选取相同的起点终点,各自选取一个不同空间位置点作为中间插值点用以区分路径。
表2 路径点空间信息表 mm
路径的起点、终点和中间插值点3点构成了2条路径线段,计算路径线段长度如表3所示。
表3 路径线段长度信息 mm
对于插值点—终点这段路径,必然存在着机械臂减速阶段,不可避免地会引起振动。对减速过程的惯性力以机械臂基坐标系按式(33)进行正交分解:
式中:L为插值点—终点路径长度;i为某一正交方向;li为L的在i上正交投影长度; θi为惯性力在i正交方向上分量。并着重关注在X-Y平面的惯性力分量。得到的详细信息如表4所示。
表4 路径的正交分解
路径耗时上,路径2最短、路径1次之、路径3耗时最长;但相应的,从平稳性角度出发,路径3的TCP坐标X-Y平面下的惯性力分量小于路径1、小于路径2,即路径3最平稳、路径1次之、路径2的振动最大。换言之,路径2耗时短但不平稳、路径3平稳但耗时长,路径1的两方面都较折中。这3条路径实际上对应了运动特点的3种不同程度,因而设计此3条路径并将之作为备选路径方案参与优化过程是具代表性及合理性的。
2.2.2 采样空间划分
贝叶斯调优来需要给定调优参数的采样空间以及采样间隔,采样空间由实验设备的参数上下界或者待研究的参数范围确定。采样间隔应该设置合理,避免出现过大过小情况。
对于高斯过程来说,需要将路径方案p进行映射。规定路径1的方案编号为1,路径2的方案编号为2,路径3的方案编号为3。
因此设置如表5调优参数信息,前述数据集D1:t中xt=(v,a,a′,p)T。
表5 优化参数信息
机械臂在燃料电池极板安装过程中,需要考虑装配效率与精度。
机械臂的运动总时间会直接影响生产效率。在实验中,通过机械臂自带的计时器得到单次运行时间T,其分辨率为0.001 s。
机械臂末端残余振动是由于运输过程中的快速定位运动会激发机器人的结构模态,引起末端执行器的惯性振动[11]。在实验中,通过在机械臂末端布置加速度传感器,采集运动结束前后0.4 s的时间段内振动信号,计算振动脉冲指数A,公式为
式中:A表示时间段内信号绝对值最大值与绝对平均值之商,振动冲击越大,A越大。单位为无量纲数,计算后小数位数保持为3位。
机械臂的绝对定位精度是用来衡量机械臂末端工具坐标系原点到达空间目标点准确度的能力,分为位置精度与方向精度,其通常由静态因素和动态因素等共同作用;其中,动态因素包括机械臂在运动过程中受惯性力等因素造成关节微变形与机械抖动等情况。根据前述实验场景,电池极板在X-Y平面安装的位置精度有着较高要求。在实验中,以30 mm×30 mm白色方片纸为定位基准,通过1 200 W像素摄像头采集每次运动结束后白色方片纸中心与参考位置中心的像素偏差,并计算转化为长度单位,计为机械臂的位置精度 µ,计算公式为
单方向上采集分辨率为0.01 mm。传感器设备图3所示。
图3 振动传感器(左)摄像头传感器(右)
选择单次运动的总时间T,振动脉冲指数A以及其末端的位置精度 µ当作判别指标。三者皆为负向指标,即在实验过程中越小越好。层次分析法的判断矩阵赋权原则:相同指标同样重要(量化值为1);时间较于振动脉冲因子同样重要(量化值为1);时间较于绝对定位偏差介于稍微重要和明显重要之间(量化值为2);振动脉冲因子较于绝对定位偏差介于稍微重要和明显重要之间(量化值为2)。按照1.4节所述,计算得到3个子指标的判断矩阵如表6。计算特征向量与一致性,得到时间权重为0.4;振动脉冲因子权重为0.4;绝对定位偏差权重为0.2,如表7所示。最后,总归一指标为各归一化后的子指标加权之和。
表6 判断矩阵
表7 子指标信息
实验结果的归一指标yt的有效位数与子指标保持一致,计算公式如下。
按照前述流程,随机森林与高斯过程的贝叶斯优化各选择3组相同的初始化参数开始调优,其中随机森林的树个数为20,贝叶斯优化的采集函数为LCB(保证95%的置信率);将得到的第一次运动参数x3+1当作第一轮实验数据,计算第一次归一指标y3+1,并将该组数据插入初始化数据集返回计算第二次运动参数;依次迭代计算。得到高斯过程实验结果如表8所示,随机森林实验结果如表9所示。
表8 高斯过程贝叶斯实验结果
表9 随机森林贝叶斯实验结果
如表10所示,将二者调优结果的归一指标与实验次数进行对比分析,规定归一指标下降幅度除以实验次数为单位收益,用于衡量调优算法综合性能。可以得到:二者均收敛到路径1,印证了路径1长度较短,X-Y平面惯性力分量较小的性质,表明了在处理离散参数二者均有着较好性能;随机森林算法的综合调优性能较高斯过程提高了15%。证明了随机森林作为概率代理模型在贝叶斯算法中有着较好的调优性能,之后也可以推广到存在离散与连续参数的工业调优场景中。
表10 随机森林与高斯过程对比实验结果
笔者从理论上分析了贝叶斯优化方法概率代理模型与采集函数,阐述了贝叶斯优化方法在ABB IRB 1410机械臂运动参数优化中的可行性。选择了随机森林与高斯过程两种概率代理模型进行实验对比:在选择的调优参数个数较少情况下,随机森林贝叶斯的优化效果比高斯贝叶斯高15%。
考虑到贝叶斯优化的概率代理模型还有很多种,因此对其他概率模型代理下的贝叶斯优算法在机械臂运动参数中的性能是下一步研究的重点。