CN112286169B - 一种工业机器人故障检测方法 - Google Patents

一种工业机器人故障检测方法 Download PDF

Info

Publication number
CN112286169B
CN112286169B CN202011127074.3A CN202011127074A CN112286169B CN 112286169 B CN112286169 B CN 112286169B CN 202011127074 A CN202011127074 A CN 202011127074A CN 112286169 B CN112286169 B CN 112286169B
Authority
CN
China
Prior art keywords
model
variable
spe
industrial robot
variables
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202011127074.3A
Other languages
English (en)
Other versions
CN112286169A (zh
Inventor
应泽
何雨辰
项剑
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Metrology University
ZHEJIANG QIANJIANG ROBOT Co.,Ltd.
Original Assignee
Zhejiang Qianjiang Robot Co ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang Qianjiang Robot Co ltd filed Critical Zhejiang Qianjiang Robot Co ltd
Priority to CN202011127074.3A priority Critical patent/CN112286169B/zh
Publication of CN112286169A publication Critical patent/CN112286169A/zh
Application granted granted Critical
Publication of CN112286169B publication Critical patent/CN112286169B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G05CONTROLLING; REGULATING
    • G05BCONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B23/00Testing or monitoring of control systems or parts thereof
    • G05B23/02Electric testing or monitoring
    • G05B23/0205Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults
    • G05B23/0218Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterised by the fault detection method dealing with either existing or incipient faults
    • G05B23/0243Electric testing or monitoring by means of a monitoring system capable of detecting and responding to faults characterised by the fault detection method dealing with either existing or incipient faults model based detection method, e.g. first-principles knowledge model

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Automation & Control Theory (AREA)
  • Feedback Control In General (AREA)
  • Manipulator (AREA)
  • Testing And Monitoring For Control Systems (AREA)

Abstract

本发明提供了一种工业机器人故障检测方法,属于工业机器人过程监控与故障检测领域。它解决了工业机器人故障检测准确性低的问题。本发明对预处理后的数据组采用设定不同马尔可夫链构建包含两类隐变量的动态潜在参照模型,能同时提取质量相关的过程变量和质量无关的过程变量信息,即将过程全部信息作为模型故障检测的监督项。统计量可以根据两个隐变量的分布计算得到,用卡方分布计算得到统计量T2和SPE的统计阈值
Figure DDA0002733965140000011
和SPElim,采集工业机器人工作过程中运行的数据并按照训练完成的动态潜在参照模型得到工作过程中的统计量
Figure DDA0002733965140000012
和SPEtest,将统计量
Figure DDA0002733965140000013
和SPEtest和统计阈值
Figure DDA0002733965140000014
和SPElim对比,可以判断工业机器人工作过程是否出现故障。本发明建立动态双隐变量的模型,使得双隐变量和动态性能够相结合,充分提取过程中与质量相关和无关的动态特性,提高了故障检测的准确性。

Description

一种工业机器人故障检测方法
技术领域
本发明属于工业机器人过程监控与故障检测领域,特别涉及一种工业机器人故障检测方法。
背景技术
工业机器人是广泛用于工业领域的多关节机械手或多自由度的机器装置,具有一定的自动性,可依靠自身的动力能源和控制能力实现各种工业加工制造功能。工业机器人被广泛应用于电子、物流、化工等各个工业领域之中。
在工业机器人运行过程中,如何有效的对系统运行过程进行监视,对帮助从业人员及时了解过程状态,消除异常行为,防止灾难性事故发生是非常重要的。其中,基于数据的多元统计过程监测方法为提高过程安全和生产效率提供了有效的建模框架,如质量相关动态模型,质量相关指的是同时包含过程变量和质量变量且过程变量和质量变量的相关性较大。
现有的工业机器人故障检测方法采用的是单一隐变量构建的模型即质量相关动态模型来进行计算和判断,包含单一隐变量质量相关动态模型只提取了过程变量中与质量变量相关的特征信息,而与质量变量无关的过程变量没有及时保留。过程信息的丢失严重限制了单一隐变量的特征表示,从而影响模型故障检测的性能。反映工业运行状态的测量数据来源不同,一些基本的过程变量可以通过比较常规的硬件传感器或者系统内部反馈快速获取,但是质量变量却受到实验成本等各方面的限制较难实现高速采样。由此看来,质量相关过程监测方法在引入质量变量时还造成了数据多采样率问题。
如中国专利申请号为201911377117.0专利公开了基于半监督自回归动态隐变量模型的故障检测方法,该方法首先收集正常工况下包含过程数据和关键质量的数据集,按数据是否被质量标记将数据集划分为过程数据和关键质量数据同时存在的有标签样本以及只有过程数据而缺少该时刻所对应的关键质量数据的无标签样本,并进一步建模。
该篇专利单一特征隐变量只保留过程变量与质量变量相关的信息,并没有保留过程变量中与质量变量无关的信息,无法充分的提取工业机器人运行状态的过程数据特征,并且该专利的采样率只有二种,对于现有工业机器人中采样率多的状况没有更深的进行考虑,从而对于工业机器人故障检测效果差,导致难以有效的对系统运行过程进行监视。
发明内容
本发明的目的是针对现有技术存在的上述问题,提出了一种工业机器人故障检测方法,其所要解决的技术问题是:如何提高工业机器人故障检测准确性。
本发明通过下列技术方案来实现:一种工业机器人的故障检测方法,其特征在于,包括以下步骤:
A、获取工业机器人在正常运行过程中的多种不同采样率的数据,并组成多采样率训练样本集,对多采样率训练样本集进行预处理得到多采样率的参照训练样本集;
B、根据参照训练样本集构建包含双隐变量的动态潜在参照模型,明确所需计算的动态潜在参照模型的模型参数,通过估计上述模型参数和双隐变量分布实现对动态潜在参照模型的训练,根据训练完成的动态潜在参照模型的双隐变量分布和模型参数构建相应的统计量T2和SPE,用卡方分布计算得到统计量T2和SPE的统计阈值
Figure GDA0003356727190000021
和SPElim
C、在统计阈值
Figure GDA0003356727190000022
和SPElim建立后,采集工业机器人运行过程中的数据,使用训练完成的动态潜在参照模型计算工业机器人运行过程中的检测统计量
Figure GDA0003356727190000023
和SPEtest,将得到的检测统计量
Figure GDA0003356727190000024
和SPEtest分别与对应的统计阈值
Figure GDA0003356727190000025
和SPElim相比,判断当前工业机器人是否存在故障。
本发明获取机器人正常运行过程中来自不同传感器与流程设备的多组数据,每组数据的采样率不同,需要提取多组数据中的特征,对于数据组中特征的提取需要先将数据组进行预处理后构建动态潜在参照模型,明确动态潜在参照模型的模型参数,对动态潜在参照模型进行训练,在动态潜在参照模型训练完成后得到模型的参数和隐变量的分布。本发明在预处理数据组后,通过设定不同的马尔可夫链构建包含两类隐变量的动态潜在参照模型,能同时提取质量相关的过程变量和质量无关的过程变量信息,即将过程全部信息作为模型故障检测的监督项。统计量可以根据两个隐变量的分布计算得到,同时用卡方分布计算得到统计量T2和SPE的统计阈值
Figure GDA0003356727190000031
和SPElim,随后采集工业机器人工作过程中运行的数据并按照训练完成的动态潜在参照模型得到工作过程中的统计量
Figure GDA0003356727190000032
和SPEtest,将统计量
Figure GDA0003356727190000033
和SPEtest和统计阈值
Figure GDA0003356727190000034
和SPElim对比,可以判断工业机器人工作过程是否出现故障。
现有大多数技术中对于工业机器人的故障检测采用的是包含单隐变量的模型,对工业机器人数据特征提取不够完整,故障检测准确性不高,虽然双隐变量的模型已经存在,但包含双隐变量的模型只限于静态和确定性假设的模型,目前现代流程工业大规模使用闭环控制策略,控制器会根据之前时刻采集到的过程变量信息来调节系统的运行状态,过程变量之间存在着样本与样本之间的相关关系,也就是过程的动态性。因此双隐变量和动态过程的动态性对于工业过程监控十分重要。本发明为了同时考虑过程数据样本之间的自相关关系和样本变量之间的互相关关系,将双隐变量策略应用于动态过程监测,建立了一种动态双隐变量的模型,该模型能使双隐变量策略和动态性能够相结合,充分提取过程中与质量相关和无关的动态特性,提高了故障检测的准确性。另外本发明是从实际问题出发做到在算法上的创新,不局限于模型或者算法之间简单的组合。
在上述的工业机器人故障检测方法中,在所述的步骤A中,获取工业机器人在正常运行过程中三种不同采样率的数据并进行预处理,得到参照训练样本集记作为G,参照训练样本集满足0均值单位协方差,参照训练样本集表示为:
Figure GDA0003356727190000041
Figure GDA0003356727190000043
Figure GDA0003356727190000044
Figure GDA0003356727190000045
表示采样率较高的两类过程变量,
Figure GDA0003356727190000046
表示采样率最低的质量相关变量;gt表示观测样本集;M1、M2、M3为不同采样率数据所包含的变量个数,T1、T2、T3为样本个数,且满足M=M1+M2+M3;T=T1>T2>T3。本方案中对获取的数据进行预处理,参照训练样本集中的数据满足0均值单位协方差有利于后续双隐变量和模型参数的求解,减少求解过程中不必要的麻烦,X(1),X(2)属于采样率较高的两类过程变量,表示可以由常规的传感器获得,X(3)表示采样率最低的质量相关变量所对应传感器采集数据过程中采样周期长。
在上述的工业机器人故障检测方法中,在上述的步骤B中,根据参照训练样本集构建包含双隐变量的动态潜在参照模型为:
qt=Fqt-1+fq
ht=Wht-1+fh
Figure GDA0003356727190000047
Figure GDA0003356727190000048
Figure GDA0003356727190000049
Figure GDA00033567271900000410
Figure GDA00033567271900000411
为模型的两类隐变量,其中qt用来描述质量变量和过程变量之间的数学关系,ht用来保留与质量变量无关的过程变量剩余信息,l,s为隐变量的维度,l,s的和远小于M,M用来表示观测维度;
Figure GDA00033567271900000412
Figure GDA00033567271900000413
为模型的载荷矩阵,
Figure GDA00033567271900000414
Figure GDA00033567271900000415
Figure GDA00033567271900000416
为不同采样率数据对应的发射矩阵,Mn,Nn分别表示不同采样率n下变量集的维数;
Figure GDA00033567271900000417
为模型误差项;模型误差项以及双隐变量初始值h1,q1满足
Figure GDA0003356727190000051
其中,~N表示服从高斯分布,Σq为特征隐变量q的协方差矩阵,Σh为特征隐变量h的协方差矩阵,
Figure GDA0003356727190000052
为第n个采样率下观测变量集,
Figure GDA0003356727190000053
为特征隐变量q的初始协方差矩阵,
Figure GDA0003356727190000054
为特征隐变量h的初始协方差矩阵,
Figure GDA0003356727190000055
为隐变量的初始化均值;通过上述动态潜在参照模型明确所需计算的模型参数为:
Figure GDA0003356727190000056
本发明的动态潜在参照模型提出包含了数据动态特性,在模型建立后,需要对双隐变量分布和模型参数的进行求解,其中隐变量qt用来描述质量变量和过程变量之间的数学关系,隐变量ht用来保留与质量变量无关的过程变量剩余信息,质量变量一般表示采样率较低的变量,相反过程变量为采样率较高的,本发明建立动态双隐变量的模型能够解决特征提取不完整和多采样率所带来的问题,提高了机器人故障检测准确性,此时动态潜在参照模型还未开始训练。
在上述的工业机器人故障检测方法中,在步骤B中,建立动态潜在参照模型的概率框架,所述的概率框架是动态潜在模型通过概率转换关系得到,表示为:
Figure GDA0003356727190000057
Figure GDA0003356727190000058
本发明中的动态潜在参照模型双隐变量的分布和模型参数的估计是在概率框架中进行的,因为在采集到的传感器数据中往往伴随着往往伴随着噪声和不确定性波动,如果采用确定性假设的方式对模型的双隐变量的分布和模型参数值进行估计不能将噪声和波动考虑进去,那么将会忽略掉数据噪声带来的模型不稳定因素,所以需要引入动态潜在参照模型的概率框架,并且在概率框架下更有利于对双隐变量的分布和模型参数的进行估计。
在上述的工业机器人故障检测方法中,建立所述概率框架的对数似然函数,采用EM算法和卡尔曼滤波数据融合算法得到对数似然函数收敛时的双隐变量的分布和动态潜在参照模型的模型参数值。
通过对数似然函数和EM算法精确估计动态潜在参照模型的模型参数和双隐变量的分布,EM算法中的E步用于估计双隐变量的分布,E步估计双隐变量的分布时需要利用卡尔曼数据融合技术,在传统单卡尔曼滤波器的基础上对滤波参数进行调整,将其扩展为多个具有不同操作时间的相关滤波器,实现模型估计与多采样率数据相融合,M步用于估计模型参数,实际上就是用最近一次E步得到双隐变量的分布求解对数似然函数关于模型参数的导函数,将此时模型得到的模型参数再代入E步中,估计当前模型参数下的双隐变量的分布,E步和M步迭代运算,直到对数似然函数收敛,最后一次迭代得到的参数就认为是最接近真实的模型参数值,两个隐变量的分布也就可以确定了。
在上述的工业机器人故障检测方法中,在步骤B中,卡尔曼滤波数据融合算法中需要构建状态指示器
Figure GDA0003356727190000062
所述的状态指示器
Figure GDA0003356727190000063
为一个非线性分段函数;
Figure GDA0003356727190000064
表示第n种采样率的数据X(n)在t时刻未能采集到变量,
Figure GDA0003356727190000065
表示第n种采样率的数据X(n)在t时刻为有效采集,n=1,2,3。通过状态指示器的用来指示多采样率测量数据中变量的缺失情况,状态指示器对于时变参数的确定有指导作用。
在上述的工业机器人故障检测方法中,在步骤B中所述构建的对数似然函数,表示为
Figure GDA0003356727190000061
其中,g1:T为1时刻到T时刻的观测样本集;q1:T为1时刻到T时刻的隐变量;h1:T为1时刻到T时刻的隐变量;
在EM算法中,E步根据当前的模型参数值,利用卡尔曼滤波数据融合算法估计双隐变量的分布,公式为:
Figure GDA0003356727190000071
Figure GDA0003356727190000072
Figure GDA0003356727190000073
Figure GDA0003356727190000074
Figure GDA0003356727190000075
Figure GDA0003356727190000076
分别表示联合隐变量
Figure GDA0003356727190000077
的均值和方差;Ψt,Rt是由模型载荷矩阵和误差项协方差矩阵组合而成的时变参数,由状态指示器
Figure GDA0003356727190000078
确定Ψt,Rt的取值,双隐变量的后验分布最终表示为:
Figure GDA0003356727190000079
Figure GDA00033567271900000710
Figure GDA00033567271900000711
其中,E<>表示关于隐变量的期望。
本文中,
Figure GDA00033567271900000712
Figure GDA00033567271900000713
分别表示联合隐变量
Figure GDA00033567271900000714
的均值和方差(即隐变量联合分布),联合隐变量的分布已知,那么两个隐变量qt和ht各自的边缘分布也可以求出。另外,对数似然函数可以用于判断双隐变量的分布和模型的参数值是否符合要求,而时变参数用于扩展卡尔曼滤波数据融合算法,将单卡尔曼滤波器扩展为多个具有不同操作时间的时变滤波器,这是解决多采样率问题的关键,其中E步得到双隐变量的联合分布及其两个隐变量各自的边缘分布将转换为隐变量的一阶矩和二阶矩(如上式所示),并用于M步中模型参数的计算,在第一次E步中,模型参数可以为满足模型矩阵大小设定的任意值,即参数的初始化。
在上述的工业机器人故障检测建模中的EM算法中,M步会根据E步得到双隐变量后验期望最大化对数似然函数的期望来计算模型参数,将后验期望代入EM算法的M步中,更新动态潜在参照模型的模型参数,公式中new表示数值的更新,
Figure GDA0003356727190000081
Figure GDA0003356727190000082
Figure GDA0003356727190000083
Figure GDA0003356727190000084
Figure GDA0003356727190000085
Figure GDA0003356727190000086
Figure GDA0003356727190000087
Figure GDA0003356727190000088
Figure GDA0003356727190000089
Figure GDA00033567271900000810
Figure GDA00033567271900000811
Figure GDA00033567271900000812
Figure GDA0003356727190000091
Figure GDA0003356727190000092
在M步算法中,得到的模型参数值后需要将模型的参数值重新代入对数似然函数中,判断对数似然函数是否收敛,如果对数似然函数不收敛需要再根据此时的模型参数代入对数似然函数估计双隐变量的分布,即重复E步得到双隐变量的分布,E步得到双隐变量的分布需要再次转化为隐变量一阶矩和二阶矩来最大化对数似然函数的期望,并用于M步中模型参数的计算,如果对数似然函数收敛,则迭代运算结束,得到对数似然函数收敛时双隐变量的分布和模型的参数,此时动态潜在参照模型训练已经完成,即通过E步和M步的迭代运算,直到对数似然函数收敛,当对数似然函数收敛时,得到训练模型参数,模型训练完成。利用本发明双隐变量的分布和模型的参数值所构建的训练样本集的统计量和统计阈值有利于提高工业机器人故障检测的准确性。
在上述的工业机器人故障检测方法中,根据对数似然函数收敛时的模型参数值得到双隐变量最终的分布均值
Figure GDA0003356727190000093
和协方差矩阵
Figure GDA0003356727190000094
统计量T2和SPE根据双隐变量的分布和模型参数构建为:
Figure GDA0003356727190000095
Figure GDA0003356727190000096
Figure GDA0003356727190000101
用卡方分布计算得到统计量
Figure GDA0003356727190000102
和SPEtest的统计阈值
Figure GDA0003356727190000103
和SPElim,表示为:
Figure GDA0003356727190000104
l,s,d分别为隐变量和观测变量的维度,d和M相等,α,β为显著水平,这里公式中last表示最终在训练阶段确定的模型参数值、双隐变量的分布和时变参数值。通过上述方式将由双隐变量和模型参数值确定的统计量转换为用于和检测统计量
Figure GDA0003356727190000105
和SPEtest比较的统计阈值
Figure GDA0003356727190000106
SPElim,本方案中的模型包含了质量无关和质量相关的特征,使得检测统计量比对的准确性大大提高,通过与统计阈值
Figure GDA0003356727190000107
SPElim的比较就能够准确得到工业机器人的数据是否存在故障。
在上述的步骤C中,采集工业机器人运行过程中与参照训练集样本对应的采样率数据作为数据样本集,并对数据样本集进行预处理得到测试样本集,使用训练完成的动态潜在参照模型和模型参数估计测试样本集对应的双隐变量分布,根据模型参数和测试样本集对应的双隐变量分布构建相应的检测统计量
Figure GDA0003356727190000108
和SPEtest,将动态潜在测试模型中构建的检测统计量
Figure GDA0003356727190000109
和SPEtest与参照训练样本集的统计阈值
Figure GDA00033567271900001010
和SPElim相比,如果
Figure GDA00033567271900001011
且SPEtest>(SPE)lim,判断测试样本集为故障数据。本方案使用工业机器人运行的正常数据构成参照训练样本集,经数据预处理后建立参照训练样本集的模型即动态潜在参照模型,其中测试样本集和参照训练样本集模型结构相同,差别只在于处理数据对象的不同,我们使用训练完成的动态潜在参照模型及其模型参数估计测试样本集对应的双隐变量分布,即成功提取了测试数据的完整特征,对于测试样本集对应双隐变量的分布估计具体步骤和参照训练样本集模型中的E步相同(进行一次E步即可)。所以我们利用已经训练完成的动态潜在参照模型及其模型参数估计测试样本集对应的双隐变量分布能完整的提取测试样本集中包含的特征,并有效判断该组数据是否存在故障,以及故障出现的时间,提高了机器人故障检测的准确性。
与现有技术相比,本发明具有以下优点:
1.本发明将静态模型扩展为动态模型,静态过程监测方法假设样本之间相互独立,未能考虑数据之间的动态关系。因此,对于普遍采用闭环控制策略的现代流程工业而言,静态模型往往无法达到预期的监测效果。动态模型在隐空间中对过程的动态特性加以描述,更适合描述实际的工业过程。在保留描述质量相关隐变量的同时,引入了另一类动态隐变量来保留过程变量中的剩余信息,模型对过程特征信息的提取不再受限于单个隐变量,储存在过程变量中的特征信息也得到了完整保留。
2.本发明利用卡尔曼数据融合技术则根据每一时刻测量变量的维度进行滤波参数的调整,模型的估计会在不同采样率之间交替完成。将多采样数据提前映射到统一维度的隐变量空间中,并提取出用于故障检测的完整特征信息,同时描述单个采样率以及不同采样率之间的相关关系。
附图说明
图1是本发明实施例的流程示意图;
图2是本发明实施例多采样率数据的结构图;
图3是本发明实施例模型联合对数似然函数的收敛趋势曲线;
图4是时本发明实施例的模型对工业机器人故障的检测结果图。
具体实施方式
以下是本发明的具体实施例,并结合附图对本发明的技术方案作进一步的描述,但本发明并不限于这些实施例。
如图1所示,一种工业机器人故障检测方法,包括:
步骤A:获取工业机器人在正常运行过程中的多种不同采样率的数据,并组成构建模型用的多采样率训练样本集,这里不同的采样率采集到的数据包括但不限于,关节水平夹角、位置反馈设备测量的关节实际位置值、速度传感器获取的关节实际速度值、电流有效值、机械臂末端点在水平x轴方向的瞬时加速度、机械臂末端点在水平y轴方向的瞬时加速度等,其中机械臂末端点在水平x轴方向的瞬时加速度、机械臂末端点在水平y轴方向的瞬时加速度为两种采样率较低的数据;对多采样率训练样本集进行预处理得到多采样率的参照训练样本集,流程如下:从工业机器人正常工序中采集到包含三种采样率的4800组样本数据,用于建模和模型参数训练;经过预处理操作后,使4800组训练样本中的变量满足0均值单位协方差;参照训练样本集可记作G,具体表示为:
Figure GDA0003356727190000121
Figure GDA0003356727190000123
步骤B中根据参照训练样本集构建包含双隐变量的动态潜在参照模型,明确所需计算的动态潜在参照模型的模型参数,通过估计上述模型参数和双隐变量分布实现对动态潜在参照模型训练,根据训练完成的动态潜在参照模型的双隐变量分布和模型参数构建相应的的统计量T2和SPE,用卡方分布计算得到统计量T2和SPE的统计阈值
Figure GDA0003356727190000124
和SPElim
本实施例中给出的动态潜在参照模型根据训练样本集的数据及其特征,采用设定不同的马尔可夫链将过程动态性和双隐变量策略相结合造型得到,马尔可夫链算法是建立参照模型的重要环节,通过设定不同的马尔可夫链算法建立的模型使得在隐空间中对过程的动态特性加以描述,更适合描述实际的工业过程。在保留描述质量相关隐变量的同时,引入了另一类动态隐变量来保留过程变量中的剩余信息,模型对过程特征信息的提取不再受限于单个隐变量,储存在过程变量中的特征信息也得到了完整保留,动态潜在参照模型的结构公式表示为:
qt=Fqt-1+fq
ht=Wht-1+fh
Figure GDA0003356727190000125
Figure GDA0003356727190000126
Figure GDA0003356727190000127
其中,
Figure GDA0003356727190000131
Figure GDA0003356727190000132
为模型的两类隐变量,qt用来描述质量变量和过程变量之间的数学关系,ht用来保留与质量变量无关的过程变量剩余信息,l,s为两类隐变量各自的维度;
Figure GDA0003356727190000133
Figure GDA0003356727190000134
为模型的载荷矩阵,
Figure GDA0003356727190000135
Figure GDA0003356727190000136
Figure GDA0003356727190000137
为不同采样率数据对应的发射矩阵,Mn,Nn分别表示不同采样率n下变量集的维数;
Figure GDA0003356727190000138
为模型误差项;隐变量初始值服从高斯分布
Figure GDA0003356727190000139
其中,~N表示服从高斯分布,Σq为特征隐变量q的协方差矩阵,Σh为特征隐变量h的协方差矩阵,
Figure GDA00033567271900001310
为第n个采样率下观测变量集,
Figure GDA00033567271900001311
为特征隐变量q的初始协方差矩阵,
Figure GDA00033567271900001312
为特征隐变量h的初始协方差矩阵,
Figure GDA00033567271900001313
为隐变量的初始化均值;
Figure GDA00033567271900001314
构成所述的模型参数。
模型的概率框架通过模型的概率转换的得到,表示为:
Figure GDA00033567271900001315
Figure GDA00033567271900001316
构建一个状态指示器
Figure GDA00033567271900001317
用来指示多采样率测量数据中变量的缺失情况;
Figure GDA00033567271900001318
表示第n个采样率数据X(n)在t时刻未能采集到数据,
Figure GDA00033567271900001319
表示第n个采样率数据X(n)在t时刻为有效采集,n=1,2,3;通过设定采样指示器的取值便可以描述任意形式的多采样率数据;图2给出了工业机器人运行过程中采集到的含三种采样率数据的数据结构示意图;其中,M1、M2、M3为不同采样率数据所包含的变量个数,本实施例图2中,表格行数不代表变量个数,空白区域表示该采样时刻未能采集到的数据,对应状态指示器为0;画斜线区域表示该采样时刻采集到数据,对应状态指示器为1。
为精确估计出动态潜在参照模型的模型参数,根据动态潜在参照模型模型的概率框架构建对数似然函数,即:
Figure GDA0003356727190000141
其中,g1:T为1时刻到T时刻的观测样本集;q1:T为1时刻到T时刻的隐变量;h1:T为1时刻到T时刻的隐变量;
E步根据当前的模型参数值,估计双隐变量分布,最开始模型参数值只需要符合本方案中模型参数的格式要求,取值可以为任意值,同时双隐变量的分布的估计需要利用卡尔曼数据融合技术和引入时变参数,在传统的单卡尔曼滤波器的基础上对滤波参数进行调整,将其扩展为多个具有不同操作时间的相关滤波器,实现模型估计与多采样率数据相融合,公式表示如下:
Figure GDA0003356727190000142
Figure GDA0003356727190000143
Figure GDA0003356727190000144
Figure GDA0003356727190000145
Figure GDA0003356727190000146
分别为联合隐变量
Figure GDA0003356727190000147
的均值和方差;Ψt,Rt是由模型载荷矩阵和误差项协方差矩阵组合而成的时变参数;Ψt,Rt的具体取值是根据状态指示器的设定而定的;工业机器人运行过程的多采样率数据结构,对Ψt,Rt有了具体设定:
当采样时刻满足
Figure GDA0003356727190000148
时;
Figure GDA0003356727190000149
当采样时刻满足
Figure GDA00033567271900001412
时:
Figure GDA0003356727190000151
当采样时刻满足
Figure GDA0003356727190000154
时:
Figure GDA0003356727190000155
特别的,在初始时刻t=1,状态指示器被设定为
Figure GDA0003356727190000156
对应双隐变量的估计需要单独定义,表示如下:
Figure GDA0003356727190000157
Figure GDA0003356727190000158
Figure GDA0003356727190000159
接着,作为优选,利用卡尔曼后向平滑算法对双隐变量的分布做进一步完善表示如下:
Figure GDA00033567271900001510
Figure GDA00033567271900001511
Figure GDA00033567271900001512
Figure GDA00033567271900001513
Figure GDA00033567271900001514
Kt和Jt表示中间变量,又称卡尔曼增益。采用卡尔曼后向平滑算法完善双隐变量的分布而得到的后验期望有利于M步中模型参数的计算。综合以上推导过程,双隐变量的后验分布最终表示为:
Figure GDA0003356727190000161
Figure GDA0003356727190000162
Figure GDA0003356727190000163
其中,E<>表示关于隐变量的期望。在M步中,根据E步更新的结果,更新模型参数,公式中new表示数值的更新,
Figure GDA0003356727190000164
Figure GDA0003356727190000165
Figure GDA0003356727190000166
Figure GDA0003356727190000167
Figure GDA0003356727190000168
Figure GDA0003356727190000169
Figure GDA00033567271900001610
Figure GDA00033567271900001611
Figure GDA00033567271900001612
Figure GDA00033567271900001613
Figure GDA00033567271900001614
Figure GDA0003356727190000171
Figure GDA0003356727190000172
Figure GDA0003356727190000173
在M步算法中,得到的模型参数值后需要将模型的参数值重新代入对数似然函数中,判断对数似然函数是否收敛,如果对数似然函数不收敛需要再根据此时的模型参数代入对数似然函数估计双隐变量的分布,即重复E步得到双隐变量的分布,E步得到双隐变量的分布需要再次转化为双隐变量的一阶矩阵和二阶矩期望来最大化对数似然函数的期望,并用于M步中模型参数的计算,如果对数似然函数收敛,则迭代运算结束,得到对数似然函数收敛时双隐变量的分布和模型的参数,此时动态潜在参照模型训练已经完成,通过E步和M步的不断迭代运算。如图3所示,E步和M经过25次迭代后对数似然函数基本趋于收敛;最后一次迭代得到的模型参数就认为是最符合要求的模型参数值。双隐变量最终的分布均值
Figure GDA0003356727190000174
和协方差矩阵
Figure GDA0003356727190000175
也可以确定,从而实现了对样本空间的划分;隐变量q对应的是质量相关子空间,变量h对应的是独立子空间,而误差项f对应的是残差空间。统计量T2和SPE的构建如下:
Figure GDA0003356727190000181
Figure GDA0003356727190000182
Figure GDA0003356727190000183
统计阈值
Figure GDA0003356727190000189
SPElim反映的是系统正常运行时的统计分布极限,可以由统计量T2和SPE和χ2分布确定,即
Figure GDA0003356727190000184
其中,l,s,d分别为隐变量和观测变量的维度,d和M相等,α,β为显著水平,这里公式中last表示最终在训练阶段确定的模型参数值、双隐变量的分布和时变参数值。
在步骤C中在统计阈值
Figure GDA00033567271900001810
和SPElim建立后,采集工业机器人运行过程中的数据,使用训练完成的动态潜在参照模型计算工业机器人运行过程中的检测统计量
Figure GDA00033567271900001811
和SPEtest,将得到的检测统计量
Figure GDA00033567271900001813
和SPEtest分别与对应的统计阈值
Figure GDA00033567271900001812
和SPElim相比,判断当前工业机器人是否存在故障。
在线收集工业机器人运行过程中与参照训练样本集对应的不同采样率数据。作为优选,采取的数据为机器人中途超载情况下采集的一批故障数据,一共收集了4400组样本并在t=2376时引入故障;接着,对采取的数据进行预处理得到测试样本集,记为Gtest,表示如下:
Figure GDA0003356727190000185
Figure GDA0003356727190000187
利用动态潜在参照模型得到对数似然函数收敛时(即动态潜在参照模型已经在正常数据中训练完成)的模型参数,得到测试样本集对应的双隐变量的分布表示为
Figure GDA0003356727190000188
Figure GDA0003356727190000191
Figure GDA0003356727190000192
Kk=Pk-1k)TkPk-1k)T+Rk]-1
Figure GDA0003356727190000193
公式中last表示最终在训练阶段确定的模型参数值和双隐变量的分布,根据模型参数和测试样本集对应的双隐变量的分布构建相应的检测统计量
Figure GDA0003356727190000194
和SPEtest,测试样本集对应的检测统计量构建与参照训练样本集保持一致。将测试样本集的检测统计量
Figure GDA0003356727190000195
和SPEtest与统计阈值
Figure GDA0003356727190000196
和SPElim相比,一旦满足
Figure GDA0003356727190000197
且SPEtest>(SPE)lim就判定测试集样本集为故障数据;如图4所示,故障数据的检测统计量
Figure GDA0003356727190000198
和SPEtest曲线在t=2376时有明显的断层发生,并超过了统计阈值线,成功检测出故障;T2统计量的故障诊断准确率达到93.23%,SPE统计量的故障诊断准确率达到98.67%,图中对于检测统计量
Figure GDA0003356727190000199
和SPEtest以log10(SPEtest)和
Figure GDA00033567271900001910
出现,主要是为了在图中更好的表示。
本发明在工业机器人故障检测中,并通过设定不同马尔可夫链将静态模型扩展为动态模型,在保留描述质量相关隐变量的同时,引入了另一类质量无关的隐变量来保留过程变量中的剩余信息,能完整的保留数据的特征。现有技术中虽然融入了质量变量(质量相关),但也遗弃了大量和质量变量无关的过程变量信息,导致质量相关下特征不完整提取。信息的丢失就严重限制了单个隐变量的特征表示,本发明建立的模型对过程特征信息的提取不再受限于单个隐变量,引入双隐变量的动态模型使得储存在过程变量中的特征信息得到完整的保留。本发明还利用卡尔曼数据融合技术根据每一时刻测量变量的维度进行滤波参数的调整,模型的估计会在不同采样率之间交替完成。将多采样数据提前映射到统一维度的双隐变量空间中,同时描述单个采样率以及不同采样率之间的相关关系。本发明能够构建包含双隐变量的动态模型,同时解决多采样率的问题,完整提取出用于故障检测的完整特征信息,提高了机器人故障检测准确性。
本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代,但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (9)

1.一种工业机器人故障检测方法,其特征在于,包括以下步骤:
A、获取工业机器人在正常运行过程中的多种不同采样率的数据,并组成多采样率训练样本集,对多采样率训练样本集进行预处理得到多采样率的参照训练样本集;
B、根据参照训练样本集构建包含双隐变量的动态潜在参照模型,所述动态潜在参照模型为:
qt=Fqt-1+fq
ht=Wht-1+fh
Figure FDA0003356727180000011
Figure FDA0003356727180000012
Figure FDA0003356727180000013
Figure FDA0003356727180000014
为模型的两类隐变量,其中qt用来描述质量变量和过程变量之间的数学关系,ht用来保留与质量变量无关的过程变量剩余信息,l,s为隐变量的维度,l,s的和远小于M,M用来表示观测维度;
Figure FDA0003356727180000015
Figure FDA0003356727180000016
为模型的载荷矩阵,
Figure FDA0003356727180000017
为不同采样率数据对应的发射矩阵,Mn,Nn分别表示不同采样率n下变量集的维数;fq,fh,fx (n)为模型误差项;模型误差项以及双隐变量初始值h1,q1满足fq~N(0,∑q),fh~N(0,∑h),
Figure FDA0003356727180000018
其中,~N表示服从高斯分布,Σq为特征隐变量q的协方差矩阵,Σh为特征隐变量h的协方差矩阵,
Figure FDA0003356727180000019
为第n个采样率下观测变量集,
Figure FDA00033567271800000110
为特征隐变量q的初始协方差矩阵,
Figure FDA00033567271800000111
为特征隐变量h的初始协方差矩阵,
Figure FDA00033567271800000112
为隐变量的初始化均值;通过上述动态潜在参照模型明确所需计算的模型参数为:
Figure FDA00033567271800000113
通过估计上述模型参数和双隐变量分布实现对动态潜在参照模型的训练,根据训练完成的动态潜在参照模型的双隐变量分布和模型参数构建相应的统计量T2和SPE,用卡方分布计算得到统计量T2和SPE的统计阈值
Figure FDA00033567271800000114
和SPElim
C、在统计阈值
Figure FDA0003356727180000021
和SPElim建立后,采集工业机器人运行过程中的数据,使用训练完成的动态潜在参照模型计算工业机器人运行过程中的检测统计量
Figure FDA0003356727180000022
和SPEtest,将得到的检测统计量
Figure FDA0003356727180000023
和SPEtest分别与对应的统计阈值
Figure FDA0003356727180000024
和SPElim相比,判断当前工业机器人是否存在故障。
2.根据权利要求1所述的工业机器人故障检测方法,其特征在于,在所述的步骤A中,获取工业机器人在正常运行过程中三种不同采样率的数据并进行预处理,得到参照训练样本集记作为G,参照训练样本集满足0均值单位协方差,参照训练样本集表示为:
Figure FDA0003356727180000025
Figure FDA0003356727180000026
Figure FDA0003356727180000027
表示采样率较高的两类过程变量,
Figure FDA0003356727180000028
表示采样率最低的质量相关变量;gt表示观测样本集;M1、M2、M3为不同采样率数据所包含的变量个数,T1、T2、T3为样本个数,且满足M=M1+M2+M3;T=T1>T2>T3
3.根据权利要求2所述的工业机器人故障检测方法,其特征在于,在步骤B中,建立动态潜在参照模型的概率框架,所述的概率框架是动态潜在参照模型通过概率转换关系得到,表示为:
Figure FDA0003356727180000029
4.根据权利要求3所述的工业机器人故障检测方法,其特征在于,对动态潜在参照模型开始进行训练,建立所述概率框架的对数似然函数,采用EM算法和卡尔曼滤波数据融合算法得到对数似然函数收敛时的双隐变量的分布和动态潜在参照模型的模型参数值。
5.根据权利要求4所述的工业机器人故障检测方法,其特征在于,卡尔曼滤波数据融合算法中需要构建状态指示器
Figure FDA0003356727180000031
所述的状态指示器
Figure FDA0003356727180000032
为一个非线性分段函数;
Figure FDA0003356727180000033
表示第n种采样率的数据X(n)在t时刻未能采集到变量,
Figure FDA0003356727180000034
表示第n种采样率的数据X(n)在t时刻为有效采集,n=1,2,3。
6.根据权利要求5所述的工业机器人故障检测方法,其特征在于,所述构建的对数似然函数表示为:
Figure FDA0003356727180000035
其中,g1:T为1时刻到T时刻的观测样本集;q1:T为1时刻到T时刻的隐变量;h1:T为1时刻到T时刻的隐变量;在EM算法中,E步根据当前的模型参数值,利用卡尔曼滤波数据融合算法估计双隐变量的分布,公式为:
Figure FDA0003356727180000036
Figure FDA0003356727180000037
Kt=Pt-1Ψt TtPt-1Ψt T+Rt]-1
Figure FDA0003356727180000038
Figure FDA0003356727180000039
分别表示联合隐变量
Figure FDA00033567271800000310
的均值和协方差;Ψt,Rt是由模型载荷矩阵和误差项协方差矩阵组合而成的时变参数,由状态指示器
Figure FDA00033567271800000311
确定Ψt,Rt的取值,双隐变量的后验分布最终表示为:
Figure FDA0003356727180000041
其中,E<>表示关于隐变量的期望。
7.根据权利要求6所述的工业机器人故障检测方法,其特征在于,在EM算法中,M步根据E步得到双隐变量后验期望最大化对数似然函数的期望来计算模型参数,将后验期望代入EM算法的M步中,求解动态潜在参照模型的模型参数,公式中new表示数值的更新,
Figure FDA0003356727180000042
Figure FDA0003356727180000043
Figure FDA0003356727180000044
Figure FDA0003356727180000045
Figure FDA0003356727180000046
Figure FDA0003356727180000047
Figure FDA0003356727180000048
Figure FDA0003356727180000049
Figure FDA00033567271800000410
Figure FDA00033567271800000411
Figure FDA0003356727180000051
Figure FDA0003356727180000052
Figure FDA0003356727180000053
Figure FDA0003356727180000054
通过E步和M步的迭代运算,直到对数似然函数收敛,当对数似然函数收敛时,得到构建统计量T2和SPE所需的模型参数值。
8.根据权利要求7所述的工业机器人故障检测方法,其特征在于,根据对数似然函数收敛时的模型参数值得到双隐变量最终的分布均值
Figure FDA0003356727180000055
和协方差矩阵
Figure FDA0003356727180000056
统计量T2和SPE根据双隐变量的分布和模型参数构建为:
Figure FDA0003356727180000057
Figure FDA0003356727180000058
Figure FDA0003356727180000059
用卡方分布计算得到统计量T2和SPE的统计阈值
Figure FDA0003356727180000061
和SPElim,表示为:
Figure FDA0003356727180000062
l,s,d分别为隐变量和观测变量的维度,d和M相等,α,β为显著水平,公式中last表示最终在训练阶段确定的模型参数值、双隐变量的分布和时变参数值。
9.根据权利要求1所述的工业机器人故障检测方法,其特征在于,在上述的步骤C中,采集工业机器人运行过程中与参照训练集样本对应的采样率数据作为数据样本集,并对数据样本集进行预处理得到测试样本集,使用训练完成的动态潜在参照模型和模型参数估计测试样本集对应的双隐变量分布,根据模型参数和测试样本集对应的双隐变量分布构建相应的检测统计量
Figure FDA0003356727180000063
和SPEtest,将动态潜在测试模型中构建的检测统计量
Figure FDA0003356727180000064
和SPEtest与参照训练样本集的统计阈值
Figure FDA0003356727180000065
和SPElim相比,如果
Figure FDA0003356727180000066
且SPEtest>(SPE)lim,判断测试样本集为故障数据。
CN202011127074.3A 2020-10-20 2020-10-20 一种工业机器人故障检测方法 Active CN112286169B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011127074.3A CN112286169B (zh) 2020-10-20 2020-10-20 一种工业机器人故障检测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011127074.3A CN112286169B (zh) 2020-10-20 2020-10-20 一种工业机器人故障检测方法

Publications (2)

Publication Number Publication Date
CN112286169A CN112286169A (zh) 2021-01-29
CN112286169B true CN112286169B (zh) 2022-02-01

Family

ID=74423923

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011127074.3A Active CN112286169B (zh) 2020-10-20 2020-10-20 一种工业机器人故障检测方法

Country Status (1)

Country Link
CN (1) CN112286169B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115248558A (zh) * 2022-08-15 2022-10-28 浙江大学 一种基于变分贝叶斯概率潜变量模型的动态过程故障检测方法
CN116339275A (zh) * 2023-02-06 2023-06-27 浙江科技学院 基于全结构动态自回归隐变量模型的多尺度过程故障检测方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105334823A (zh) * 2015-11-05 2016-02-17 浙江大学 基于有监督的线性动态系统模型的工业过程故障检测方法
CN105759787A (zh) * 2016-03-19 2016-07-13 浙江大学 一种基于切换有监督线性动态系统模型的故障诊断方法
WO2017094267A1 (ja) * 2015-12-01 2017-06-08 株式会社Preferred Networks 異常検出システム、異常検出方法、異常検出プログラム及び学習済モデル生成方法
CN109085805A (zh) * 2018-07-24 2018-12-25 浙江科技学院 一种基于多采样率因子分析模型的工业过程故障检测方法
CN109917777A (zh) * 2019-04-16 2019-06-21 浙江科技学院 基于混合多采样率概率主成分分析模型的故障检测方法
CN111142501A (zh) * 2019-12-27 2020-05-12 浙江科技学院 基于半监督自回归动态隐变量模型的故障检测方法
CN111507003A (zh) * 2020-04-20 2020-08-07 中国计量大学 基于质量相关动态特性提取的脱丁烷塔关键变量预测方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
TWI399660B (zh) * 2008-07-09 2013-06-21 Inotera Memories Inc 偵測半導體製程變異之方法
CN108267172B (zh) * 2018-01-25 2020-06-30 神华宁夏煤业集团有限责任公司 矿用智能机器人巡检系统
CN110209145B (zh) * 2019-05-16 2020-09-11 浙江大学 一种基于核矩阵近似的二氧化碳吸收塔故障诊断方法

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105334823A (zh) * 2015-11-05 2016-02-17 浙江大学 基于有监督的线性动态系统模型的工业过程故障检测方法
WO2017094267A1 (ja) * 2015-12-01 2017-06-08 株式会社Preferred Networks 異常検出システム、異常検出方法、異常検出プログラム及び学習済モデル生成方法
CN105759787A (zh) * 2016-03-19 2016-07-13 浙江大学 一种基于切换有监督线性动态系统模型的故障诊断方法
CN109085805A (zh) * 2018-07-24 2018-12-25 浙江科技学院 一种基于多采样率因子分析模型的工业过程故障检测方法
CN109917777A (zh) * 2019-04-16 2019-06-21 浙江科技学院 基于混合多采样率概率主成分分析模型的故障检测方法
CN111142501A (zh) * 2019-12-27 2020-05-12 浙江科技学院 基于半监督自回归动态隐变量模型的故障检测方法
CN111507003A (zh) * 2020-04-20 2020-08-07 中国计量大学 基于质量相关动态特性提取的脱丁烷塔关键变量预测方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
Double Layer Distributed Process Monitoring;YUCHEN HE;《IEEE Access》;20190214;全文 *

Also Published As

Publication number Publication date
CN112286169A (zh) 2021-01-29

Similar Documents

Publication Publication Date Title
CN107480440B (zh) 一种基于两阶段随机退化建模的剩余寿命预测方法
CN111275288B (zh) 基于XGBoost的多维数据异常检测方法与装置
CN109657945B (zh) 一种基于数据驱动的工业生产过程故障诊断方法
CN105700518B (zh) 一种工业过程故障诊断方法
CN109146246B (zh) 一种基于自动编码器和贝叶斯网络的故障检测方法
CN112286169B (zh) 一种工业机器人故障检测方法
CN111913803B (zh) 一种基于akx混合模型的服务负载细粒度预测方法
CN109062189B (zh) 一种用于复杂故障的工业过程故障诊断方法
CN109917777B (zh) 基于混合多采样率概率主成分分析模型的故障检测方法
CN107861492A (zh) 一种基于裕度统计量的广义非负矩阵分解故障监测方法
CN111142501A (zh) 基于半监督自回归动态隐变量模型的故障检测方法
CN108960329B (zh) 一种包含缺失数据的化工过程故障检测方法
CN105607631B (zh) 间歇过程弱故障模型控制限建立方法及弱故障监测方法
CN111241744A (zh) 一种基于双向lstm的低压铸造机时间序列数据异常检测方法
CN109298633A (zh) 基于自适应分块非负矩阵分解的化工生产过程故障监测方法
CN112666918A (zh) 一种基于在线压缩keca自适应工业过程故障检测方法
CN114962390A (zh) 液压系统故障诊断方法、系统及作业机械
CN105629959A (zh) 一种工业过程故障检测方法
CN116467653A (zh) 一种基于概率分布和XGBoost决策算法的织机异常数据处理方法
CN110084301B (zh) 一种基于隐马尔可夫模型的多工况过程工况辨识方法
CN114418420A (zh) 基于因果推断的竞争风险生存分析方法
CN118172910A (zh) 基于卷积神经网络的煤矿架空乘人装置安全预警系统
CN111045415A (zh) 一种基于局部概率密度双子空间的多模态过程故障检测方法
CN111414943A (zh) 一种基于混合隐朴素贝叶斯模型的异常检测方法
CN111062848A (zh) 一种消防工程异常状态监测的智能监控方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
TR01 Transfer of patent right

Effective date of registration: 20220302

Address after: 317500 south side of Wuyang Road, Wenqiao Town Industrial Park, Wenling City, Taizhou City, Zhejiang Province

Patentee after: ZHEJIANG QIANJIANG ROBOT Co.,Ltd.

Patentee after: China Metrology University

Address before: 317500 south side of Wuyang Road, Wenqiao Town Industrial Park, Wenling City, Taizhou City, Zhejiang Province

Patentee before: ZHEJIANG QIANJIANG ROBOT Co.,Ltd.

TR01 Transfer of patent right