生产机台组件的多变量预测保养方法与其计算机程序产品与流程

1.本发明是有关于一种生产机台组件的多变量预测保养方法与其计算机程序产品，且特别是有关于一种基于剩余使用寿命(remaining useful life；rul)预测的生产机台组件的多变量预测保养方法与其计算机程序产品。


背景技术：

2.生产机台是任何制造厂不可缺少的部分。生产机台中的组件、模块或装置(例如：加热器、压力模块和节流阀(throttle valve)等)的失效会引起生产异常，导致不良的产品品质和/或降低产能，因而造成重大损失。
3.一般，解决上述问题最常用的方法是定期的预防保养(preventive maintenance；pm)。即，在预设时间间隔下执行保养相关作业。此预设时间间隔基本上是根据标的装置(td)的平均故障时间间隔(mean time between failure)来决定。因此，如何安排适当的pm计划通常是工厂的关键议题。一个不当的定期pm计划会增加维修成本或降低产能。
4.预测性维护目的在找出标的设备(即生产机台的组件)何时生病并在标的设备死亡发生之前以进行即时性维护，以避免意外的标的设备停机时间。通过这种方式，不仅提高了生产机台的稼动率和制造品质，而且还可以降低预防性保养中的过度维护的额外成本。
5.为改善机台保养计划以增加晶圆厂的绩效，国际半导体技术制造协会(international sematech manufacturing initiative；ismi)提出一种预测性和预防性保养(predictive and preventive maintenance；ppm)的指标。如ismi所定义，ppm包含预防保养(pm)、基于条件的保养(condition-based maintenance；cbm)、预测保养(predictive maintenance；pdm)和故障后维修(breakdown maintenance；bdm)。其中，ismi主张cbm和pdm的技术应被发展，并以单一模块或多个模块的型式被使用，使得终端使用者能有效率地使用这些技术。cbm的定义为：“在指出机台将要失效或机台的性能正在恶化的一或多个指标出现后进行保养”。错误侦测及分类(fault detection and classification；fdc)是一种与cbm相关的方法，其定义为：“监控机台与工厂数据以评估机台的健康，并在侦测到错误时发出警报和/或关闭机台”。另一方面，pdm是一种应用预测模型的技术，找出设备状态信息与保养信息间的关联，来预测机台或标的装置(td)的剩余寿命(remaining useful life；rul)，以达到减少非计划性停机的保养事件的目标。
6.在一些已知技术中，只用单一个特征来预测设备的老化特征，进而判断设备的剩余使用寿命，如何增加预测的准确度，为此领域技术人员所关心的议题。


技术实现要素：

7.由于已知演算法的限制，当标的设备即将死亡时，若标的设备的老化特征突然上升或变得平滑，指数模型可能无法跟上即时预测甚至错误地预测标的设备的rul。
8.本发明的一目的是在提供一种生产机台组件的多变量预测保养方法与其计算机程序产品，借以即时并准确地预测生产机台组件的rul，而可及时地进行生产机台组件的维
修。
9.本发明的另一目的是在提供一种生产机台组件的多变量预测保养方法与其计算机程序产品，通过提出了预警机制和死亡相关指数(death correlation index；dci)，以便在生产机台组件可能很快死亡状态时立即进行维护，并以数据化的方式来呈现生产机台组件进入死亡状态的可能性。
10.根据本发明的一态样，提供一种生产机台组件的多变量预测保养方法。首先，获得生产机台的组件依序处理多个工件时所使用的多组制程数据，其中每一个组制程数据包含多个参数的数值。然后，获得分别对应至此些组制程数据的多个事件指示值，其中此些事件指示值指出当此组件在处理每一个工件时该组件是否发生异常事件。接着，分别使用多个演算法将每一组制程数据的此些参数的数值转换成多个参数指标的数值。然后，对此些组制程数据中的每一个参数指标与此些事件指示值进行一相关性分析，而获得分别对应至此些参数指标的多个相关系数。接着，选取对应至此些相关系数中一最大者的参数指标为一老化特征，并设定一辅助老化特征。接着，进行第一判断步骤，以根据每一个工件对应的此老化特征的数值，来判断此组件在处理此些工件时是否处于一生病状态，其中一旦此组件在处理该些工件的一者时是处于生病状态时，则将此工件设定为一样本选取点。然后，进行一多变量建模步骤，此多变量建模步骤包含：使用此样本选取点前n个工件所对应的n组制程数据中对应至此老化特征的n个数值以及对应辅助老化特征的n个数值为一组建模样本数据，其中n为正整数。接着，对老化特征与辅助老化特征执行一格兰杰因果关系检验(granger causality test)以判断老化特征与辅助老化特征之间的相关性，若不相关则从组建模样本数据中删除辅助老化特征。接着，以使用此组建模样本数据并根据一多变量时间序列预测演算法来建立一老化特征预测模型，而获得依工件生产次序排列的此老化特征的多个预测数值。接着，使用每一个工件的一工件处理时间和此组件无法使用时的此老化特征的一死亡规格值，来将此些预测数值转换成的多个剩余使用寿命(remaining useful life；rul)预测值(rul
t
)，其中t代表第t个工件，t为整数。然后，进行第二判断步骤，以根据此些剩余使用寿命预测值(rul
t
)来判断此组件是否需要更换或维修。
11.在一些实施例中，在前述的第一判断步骤中，首先以一组转换公式分别将每一组制程数据的上述老化特征的数值转换成分别对应至工件的多个装置健康指数(device health index；dhi)，此组转换公式为：
12.当
13.当
14.当
15.当
16.当
17.当
18.其中为此些组制程数据的该老化特征的数值的平均值，其中为所对应的转换值；
19.max y
t
为此些组制程数据的该老化特征的最大数值，max y
t_mapping
为max y
t
所对应的转换值；
20.min y
t
为此些组制程数据的此老化特征的最小数值，min y
t_mapping
为min y
t
所对应的转换值；
21.lsl为规格下限；lcl为管制下限；ucl为管制上限；usl为规格上限；lsl
_mapping
为lsl所对应的转换值；lcl
_mapping
为lcl所对应的转换值；ucl
_mapping
为ucl所对应的转换值；usl
_mapping
为usl所对应的转换值。
22.然后，依序判断此些装置健康指数是否大于或等于一门槛值，并将此些装置健康指数中最先大于或等于一门槛值的一者所应的工件设定为前述的样本选取点。
23.在一些实施例中，在前述的多变量建模步骤中，首先使用一向量自回归模型(vector autoregression model；var)为前述的多变量时间序列预测演算法，来建立前述的老化特征预测模型。使用偏自相关函数(partial autocorrelation function；pacf)选出此向量自回归模型的最大落后期数。然后，对建模样本数据中的数值进行一白噪音检定，其中当此些数值是白噪音时，则加入前述的样本选取点前n+1个工件所对应的又一组制程数据中对应至此老化特征的数值至建模样本数据。接着，使用此向量自回归模型的最大落后期数，来建立多个向量自回归模型组合。然后，使用一讯息准则演算法，来计算出每一个向量自回归模型组合的讯息量。接着，选出此些向量自回归模型组合中具有最大讯息量的一者为一最佳模型。
24.在一些实施例中，前述的讯息准则演算法为贝氏讯息准则(bayesian information criteria；bic)。
25.在一些实施例中，前述的多变量建模步骤还包含：判断建模样本中的数值的变异数是否会随着时间而越来越大，其中当此些数值的变异数随着时间而越来越大时，对建模样本数据的每一个数值进行对数转换；对此些数值进行一单根检定(unit root test)，以确认依序排列的此些数值是否为稳态状态，其中当此些数值不是稳态状态时，对建模样本数据的每一个数值进行差分转换。
26.在一些实施例中，前述的单根检定为扩充迪基-福勒(augmented dickey-fuller；adf)检验或kwiatkowski-phillips-schmidt-shin(kpss)检定。
27.在一些实施例中，前述的第二判断步骤包含：
28.判断(rul
t
─
rul
t-1
)/rul
t-1
是否大于或等于一门槛值，而获得第一结果，其中t-1代表第t-1个工件；判断rul
t
是否小于一维修缓冲时间，而获得第二结果，其中前述的组件必须在此维修缓冲时间进行维修；当第一结果和第二结果均为否时，此组件处于生病状态但未急速恶化，不需进行维修；当第一结果为否而第二结果为是时，此组件未急速恶化但其剩余使用寿命不足，需进行维修；当第一结果为是而第二结果为否时，此组件急速恶化，若处理连续第t+i个工件的每一者的第一结果为是而第二结果为否，则需检查或维修此组件，其中i为正整数；以及当第一结果和第二结果均为是时，此组件需进行维修。
29.在一些实施例中，前述的第二判断步骤包含：
30.以一组转换公式分别将每一组制程数据的前述老化特征的数值转换成分别对应至工件的前述组件的多个死亡相关指数(death correlation index；dci)，此组转换公式为：
[0031][0032]
其中y
death
为此组件在死亡状态时所对应的此老化特征的数值，y
t-1
为此组件在处理第t-1个工件时所对应的此老化特征的数值，conv为共变异数计算，var为变异数计算；
[0033]
当dci
t
大于一门槛值时，代表此组件在处理第t个工件时接近死亡状态，其中此门槛值的计算是根据dci
t
的标准差。
[0034]
在一些实施例中，前述的组件为一加热器、一压力模块、一节流阀、一无油衬套或一轴承，该些参数包含：一轴偏度、一阀开度、一振动振幅、一驱动电压、一驱动电流、一温度和一压力。
[0035]
在一些实施例中，前述的参数指标包含：一转换至频域后的k倍频(其中k大于0)、一整体相似度指标(global similarity index；gsi)、一统计数据分布的峰度(kurtosis)、一统计数据分布的偏度(skewness)、一标准差、一均方根(root mean square)、一平均值、一最大值和一最小值。
[0036]
根据本发明的又一态样，提供一种计算机程序产品，当计算机载入此计算机程序产品并执行后，能够完成前述的生产机台组件的多变量预测保养方法。
[0037]
因此，应用本发明实施例，可即时并准确地预测生产机台组件的rul，而可及时地进行生产机台组件的维修；并可在生产机台组件可能很快死亡状态时立即进行维护，且以数据化的方式来呈现生产机台组件进入死亡状态的可能性。
附图说明
[0038]
为了更完整了解实施例及其优点，现参照结合所附附图所做的下列描述，其中
[0039]
图1a为绘示根据本发明一些实施例的用以说明生产机台组件的多变量预测保养方法的方块示意图；
[0040]
图1b为绘示根据本发明一些实施例的用以说明事件指示值与参数指标(平均开度)的关系；
[0041]
图1c为绘示根据本发明一些实施例的用以说明事件指示值与参数指标(1/4倍频)的关系；
[0042]
图2a和图2b为绘示根据本发明一些实施例的多变量建模步骤的流程示意图；
[0043]
图3a为绘示本发明一实施例的老化特征的偏自相关函数的示意图；
[0044]
图3b为绘示本发明一实施例的辅助老化特征的偏自相关函数的示意图；
[0045]
图4为绘示根据本发明一些实施例的用以说明预警模式的方块示意图；
[0046]
图5为绘示组件的轴偏度的预测结果；
[0047]
图6为绘示组件的rul的预测结果；
[0048]
图7为绘示组件的一dci(tsp)的预测结果；以及
[0049]
图8为绘示组件的另一dci(tsp
mva
)的预测结果。
[0050]
【符号说明】
[0051]
200:制程数据
[0052]
210:数据品质检查
[0053]
220:转换成参数指标
[0054]
222,224,226,228:曲线
[0055]
230:演算法库
[0056]
240:dhi模块
[0057]
250:第一判断步骤(dhi《0.7)
[0058]
260:应用多变量的rul预测建模步骤
[0059]
270:预警模式
[0060]
280:dci模式
[0061]
300:选取建模样本数据
[0062]
302:建模样本中的数值的变异数随时间而越来越大？
[0063]
304:对建模样本数据的数值进行对数转换
[0064]
306:建模样本数据的数值是否为稳态状态
[0065]
308:对建模样本数据的数值进行差分转换
[0066]
309-1:执行格兰杰因果关系检验
[0067]
309-2:从组建模样本数据中删除辅助老化特征(z
t
)
[0068]
310:使用pacf选出var(p)模型的最大落后期数
[0069]
312:对建模样本数据的数值进行白噪音检定
[0070]
314:加入又一老化特征的数值至建模样本数据
[0071]
316:建立var模型组合
[0072]
318:计算出var模型组合的讯息量
[0073]
320:选出最佳var模型
[0074]
324:去除最佳var模型中不显著的预测成分
[0075]
326:对最佳var模型的残差进行检定
[0076]
328:确认最佳var模型
[0077]
400:rul
t
是否下降大于或等于门槛值
[0078]
410:rul
t
是否小于维修缓冲时间
[0079]
420:rul
t
是否小于维修缓冲时间
[0080]
510,520:预测值
[0081]
530:实际值
[0082]
540:死亡规格值
[0083]
550:生病规格值
[0084]
560,570:rul
[0085]
580:dci门槛值
[0086]
dci:死亡相关指数
[0087]
rul:剩余使用寿命
[0088]
tsp:时间序列预测演算法
similarity index；gsi)演算法可将此组阀开度的时序数据转换成一gsi值；演算法库230中的时频转换演算法可将此组振动振幅的时序数据转换成一1/4倍频、一1/2倍频、一2倍频、一4倍频等参数指标的数值。gsi演算法可参照美国专利前案第8095484b2号。以上所述的演算法库230中的演算法仅是举例说明，故本发明实施例并不在此限。
[0096]
另一方面，本方法获得分别对应至此些组制程数据的多个事件指示值，其中这些事件指示值指出当生产机台组件在处理每一个工件时生产机台组件是否发生异常事件。请参照图1b，图1b为绘示根据本发明一些实施例的用以说明事件指示值与参数指标(平均开度)的关系，其中曲线222指出生产机台组件在处理各工件时的平均开度；曲线224指出生产机台组件在处理各工件时是否有异常事件发生。如图1b所示，生产机台组件在处理第289个工件之前，均无异常事件发生，其事件指示值可为例如“0”；生产机台组件在处理第289个工件之后，有异常事件发生，其事件指示值可为例如“1”。请参照图1c，图1c为绘示根据本发明一些实施例的用以说明事件指示值与参数指标(1/4倍频)的关系，其中曲线226指出生产机台组件在处理各工件时在1/4倍频的振动幅度值；曲线228指出生产机台组件在处理各工件时是否有异常事件发生。如图1c所示，生产机台组件在处理第287个工件之前，均无异常事件发生，其事件指示值可为例如“0”；生产机台组件在处理第287个工件之后，有异常事件发生，其事件指示值可为例如“1”。
[0097]
然后，对此些组制程数据中的每一个参数指标与事件指示值进行一相关性分析，而获得分别对应至此些参数指标的多个相关系数，如表一所示。接着，选取对应至此些相关系数中一最大者的参数指标为一老化特征(y
t
)，如表一所示的开度。此外，也可以挑选一辅助老化特征，此辅助老化特征可以与上述的老化特征(y
t
)结合，一起用来预测下个时间点的老化特征。由于老化特征(y
t
)与辅助老化特征合并成为一个向量，在此提出的方法也可称为向量自回归模型(vector autoregression model；var)。辅助老化特征可以是表一除了开度以外的任意其他特征(如轴偏度)，在一些实施例中可以人为的挑选辅助老化特征，例如挑选温度做为辅助老化特征。在一些实施例中也可以用任意的特征选择方法，例如自适应增强(adaptive boosting)，但本揭露并不在此限。
[0098]
演算法相关系数开度0.75171/4倍频0.68451/2倍频-0.41481倍频-0.39242倍频0.35343倍频-0.2734std-0.2725skewness-0.2723kurtosis-0.2618rms0.1767average0.1454max0.0012min0.0012
gsi0.0012
[0099]
表一
[0100]
接着，进行第一判断步骤250，以根据每一个工件对应的此老化特征的数值，来判断生产机台组件在处理此些工件时是否处于一生病状态，其中一旦生产机台组件在处一个工件时是处于生病状态时，则将此工件设定为一样本选取点(t)。以下举例说明第一判断步骤的一种实施方式，然本发明实施例并不在此限。如图1a所示，第一判断步骤中，首先，将此些组制程数据的老化特征(y
t
)的数值输入至装置健康指数(device health index；dhi)模块240，dhi模块240以一组转换公式分别将每一组制程数据的此老化特征(y
t
)的数值转换成分别对应至工件的多个装置健康指数(dhi)。然后，依序判断此些装置健康指数是否大于或等于一门槛值(例如：0.7)，并将此些装置健康指数中最先大于或等于一门槛值的一者所应的工件设定为样本选取点。此组转换公式为：
[0101]
当
[0102]
当
[0103]
当
[0104]
当
[0105]
当
[0106]
当
[0107]
其中为此些组制程数据的此老化特征的数值的平均值，其中为所对应的转换值；
[0108]
max y
t
为此些组制程数据的此老化特征的最大数值，max y
t_mapping
为max y
t
所对应的转换值；
[0109]
min y
t
为此些组制程数据的此老化特征的最小数值，min y
t_mapping
为min y
t
所对应的转换值；
[0110]
lsl为规格下限；lcl为管制下限；ucl为管制上限；usl为规格上限；lsl
_mapping
为lsl所对应的转换值；lcl
_mapping
为lcl所对应的转换值；ucl
_mapping
为ucl所对应的转换值；usl
_mapping
为usl所对应的转换值。dhi演算法参照类似于美国专利前案第10,242,319b2号。本发明的实施例引用此美国专利前案第10,242,319b2号的相关规定(incorporated by reference)。
[0111]
接着，进行一多变量建模步骤260。首先，使用样本选取点(t)前n个工件所对应的n组制程数据中对应至此老化特征(y
t
)的n个数值以及对应辅助老化特征的n个数值做为一组建模样本数据，其中n为正整数。然后，对老化特征(y
t
)与辅助老化特征执行一格兰杰因果关系检验(granger causality test)以判断老化特征(y
t
)与辅助老化特征之间的相关性，若不相关则从此组建模样本数据中删除辅助老化特征。例如：假设老化特征有一个(如
sx1)，而后选的辅助老化特征有两个(如sy1、sz1)，其中辅助老化特征sy1与老化特征(y
t
)有相关，则保留；辅助老化特征sz1与老化特征(y
t
)不相关，则删除，进而使此组建模样本数据对应老化特征sx1与辅助老化特征sy1(如图3a、3b所示)。在一实施例中，老化特征sx1与辅助老化特征sy1、sz1可为轴偏度，轴偏度为工件(如轴承)的轴所偏离的角度。老化特征sx1为工件的轴于x平面所偏离的角度，辅助老化特征sy1为工件的轴于y平面所偏离的角度，辅助老化特征sz1为工件的轴于z平面所偏离的角度，但本发明不以此为限；接着，以使用此组建模样本数据并根据一多变量时间序列预测演算法来建立一老化特征预测模型，而获得依工件生产次序排列的此老化特征(y
t
)的多个预测数值。接着，使用每一个工件的一工件处理时间和生产机台组件无法使用时的此老化特征(y
t
)的一死亡规格值，来将此些预测数值转换成的多个剩余使用寿命(remaining useful life；rul)预测值(rul
t
)，其中t代表第t个工件，t为整数。值得一提的是，本发明使用多变量时间序列预测演算法及多变量建模步骤260的缘由是因已知技术只采用一个变数来做预测，其预测准确度有限。而本发明的多变量的时间序列预测扩充了只能使用一个变数的限制，所预测的老化特征不仅取决于过去的老化特征，也依赖其他特征(即辅助老化特征)，因此采用多变数以提升预测准确度。至于多变量时间序列预测演算法和详细的多变量建模步骤260将于后说明。
[0112]
然后进行一第二判断步骤，以根据此些剩余使用寿命预测值(rul
t
)来判断生产机台组件是否需要更换或维修。在一些实施例中，第二判断步骤包含有预警模式270和dci模式280。至于预警模式270和dci模式280将于后说明。
[0113]
请参照图2a和图2b，图2a和图2b为绘示根据本发明一些实施例的多变量建模步骤260的流程示意图。首先，进行步骤300，以选取使用前述的样本选取点(t)(即dhi》0.7的工件)前n个工件所对应的n组制程数据中对应至此老化特征(y
t
)的n个数值(老化特征实际值)以及对应至辅助老化特征的n个数值做为一组建模样本数据(ym)，其中n为正整数，例如30，即dhi》0.7的工件前30个工件，建模样本数据中的老化特征表示为ym＝{y
t-30
,y
t-29
,...,y
t-2
,y
t-1
}，而辅助老化特征表示为zm＝{z
t-30
，z
t-29
，...，z
t-2
，z
t-1
}。
[0114]
接着，进行步骤302，以判断此组建模样本中老化特征的数值的变异数是否会随着时间而越来越大。换言之，若y
t
＝(1+α)
×yt-1
，其中α大于0，代表y
t
随时间增加，var(y
t
)随α成长，则进行步骤304，否则进行步骤306。在步骤304中，对此组建模样本数据的每一个数值进行对数转换，以强迫数据的增加率分布具有某种程度的规则性，然后进行步骤306。
[0115]
在步骤306中，对此组建模样本中的数值进行一单根检定(unit root test)，以确认此组建模样本中依序排列的数值是否为稳态状态，其中当此些数值不是稳态状态时，对此组建模样本数据的每一个数值进行差分转换(步骤308)。确认时序是否为稳态状态的演算法有扩充迪基-福勒检验(augmented dickey-fuller test；adf)检验或kwiatkowski-phillips-schmidt-shin(kpss)检定等。adf检验和kpss检定的公式与应用方法是习于此技艺的人士所知，故不在此赘述。
[0116]
当依序排列的数值不为稳态状态(步骤306的结果为否)时，进行步骤308，以对此组建模样本数据的每一个数值进行差分转换，以使在依序排列的数值(时序)达到稳态。差分转换的公式为：分转换的公式为：当依序排列的数值为稳态状态(步骤306的结果为是)时，进行步骤309-1，对老化特征(y
t
)与辅助老化特征(z
t
)
执行格兰杰因果关系检验(granger causality test)以判断老化特征(y
t
)与辅助老化特征(z
t
)之间的相关性，若不相关则在步骤309-2从组建模样本数据中删除辅助老化特征(z
t
)，若相关则保留辅助老化特征(z
t
)。本领域具有通常知识者当可理解格兰杰因果关系检验，在此并不详细赘述。
[0117]
接下来进行步骤310，以使用偏自相关函数(partial autocorrelation function；pacf)选出向量自回归(var(p))模型的最大落后期数p。其公式为：
[0118]
b＝arg max(ρρk).
ꢀꢀ
(1)；
[0119][0120]
γk＝cov(y
t
,y
t-k
)＝e(y
t-μ)(y
t-k-μ)
ꢀꢀ
(3-1)；
[0121]
ρρk＝corr(y
t
，y
t-k
|y
t-1
，y
t-2
，...，y
t-k+1
)
ꢀꢀ
(3-2)；
[0122]
其中b为与y
t-1
相关的最大pacf的时间(工件号码)；var(y
t
)＝var(y
t+k
)＝e(y
t-μ)2＝γ0.；e[
·
]为期望值函数；μ为y
t
的平均值(mean)。举例来说，参阅图3a与图3b，其中图3a为绘示本发明一实施例的老化特征sx1的偏自相关函数的示意图；图3b为绘示本发明一实施例的辅助老化特征sy1的偏自相关函数的示意图。老化特征sx1与辅助老化特征sy1为轴偏度。由图3a与图3b可知，老化特征sx1经pacf选出的向量自回归(var(p))模型的最大落后期数p为3，而辅助老化特征sy1经pacf选出的向量自回归(var(p))模型的最大落后期数p为3。借此，本发明透过pacf可决定var(p)模型的最大落后期数p。
[0123]
接着，进行步骤312，以对此组建模样本数据中的各数值(老化特征的数值及辅助老化特征的数值)进行白噪音检定，白噪音检定的演算法有ljung-box检定等。步骤312的主要目的在确认时序是否为白噪音(即状态与状态之间彼此不相关)。当此些数值是白噪音时，此组建模样本数据需要多一点数据，故加入样本选取点前n+1个工件所对应的又一组制程数据中对应至老化特征的数值(老化特征实际值，即y
t-(n+1)
或y
t-31
)与辅助老化特征的数值(辅助老化特征实际值，即z
t-(n+1)
或z
t-31
)至此组建模样本数据(步骤314)。
[0124]
然后，进行步骤316，建立多个var(p)模型组合，并透过var(p)模型的最大落后期数p决定var(p)模型组合的数量。例如：当前述var(p)模型的最大落后期数p为3时，建立3种模型组合：var(1)、var(2)、var(3)。var(1)、var(2)、var(3)分别为var(p)的p等于1、2、3。以下说明var(p)模型。
[0125]
var(p)模型的定义为：
[0126][0127]
其中为在时间点t(第t个工件)的老化特征预测值；b为常数向量；为在时间点(第t个工件)的var(p)模型的最小平方预估系数，i＝1,2,
…
,p；y
t-i
为在时间点t-i(第t-i个工件)的老化特征实际值(对应图3a的老化特征sx1)；z
t-i
为在时间点t-i(第t-i个工件)的辅助老化特征实际值(对应图3b的辅助老化特征sy1)；βi为所要训练(寻找)的系数；ε
t
为在时间点t(第t个工件)的白噪音项，亦即误差向量。如果在步骤309-2中删除了辅助老化特征，则在式子(4)中的系数βi都设为0。
[0128]
然后，进行步骤318，以使用一讯息准则演算法，来计算出每一个var(p)模型组合的讯息量，其中讯息准则演算法为贝氏讯息准则(bayesian information criteria；bic)。
bic演算法的公式如下:
[0129][0130][0131]
其中sse为误差平方的和；m为建模样本数据的大小。
[0132]
表二为以bic演算法来计算出每一个var(p)模型组合var(1)、var(2)、var(3)的讯息量的举例说明。
[0133]
p123bic-106.7993-97.3780-90.0793
[0134]
表二
[0135]
接着，进行步骤320，以选出此些var(p)模型组合var(1)、var(2)、var(3)中具有最大讯息量(即最小bic)的一者为一最佳模型，例如：var(1)。
[0136]
然后，进行步骤322，以去除最佳模型中不显著预测子(参数)。当预测子的预估系数大于95％的信心区间时，此预测子不显著，应予以去除。在正常分布的假设下，95％的信心区间等于1.96，最佳模型中不显著预测子(参数)的判断公式如下:
[0137][0138]
|βi|＞1.96
×
s.e.(βi)
ꢀꢀ
(8)；
[0139]
其中i＝1,2,
…
,p；s.e.(
·
)为系数的标准差。
[0140]
接着，进行步骤324，以在除不显著预测子(参数)后针对var(p)模型重行估计权重。
[0141]
然后，进行步骤326，以对重新估计的模型的残差进行检定，此检定的演算法有ljung-box检定等。重新估计的模型的残差均已被解释时，确认重新估计的模型为最佳模型(步骤328)，并以此获得依工件生产次序排列的老化特征(y
t
)的多个预测数值。接着，进行步骤330，以使用生产机台组件处理每一个工件的工件处理时间(dt)和生产机台组件无法使用时的老化特征(y
t
)的死亡规格值，来将该些预测数值转换成的多个剩余使用寿命(remaining useful life；rul)预测值(rul
t
)，其公式为rul
t
＝k
d-k
t
.，其中t代表第t个工件；k
t
代表第t个工件所对应的时间点(即第7
×
dt)，t为整数kd代表老化特征(y
t
)的死亡规格值所对应的时间点。
[0142]
在获得rul预测值(rul
t
)后，进行第二判断步骤，以根据rul
t
来判断生产机台组件是否需要更换或维修。如图1a所示，在一些实施例中，第二判断步骤包含有预警模式270和dci模式280。
[0143]
当rul预测值(rul
t
)有大幅度的下降或在靠近死亡状态震荡时，使用者难以判断生产机台组件是否需要更换或维修。因此，本发明实施例提出预警模式来解决此问题。请参照图4，图4为绘示根据本发明一些实施例的用以说明预警模式的方块示意图。首先，在第一阶中进行步骤400，以判断目前的rul
t
相较于前一个rul
t-1
是否下降大于或等于一门槛值(例如30％)，即(rul
t-1-rul
t
)/rul
t-1
≥0.3是否成立？当步骤400的结果为是时，在第一阶中进行步骤410或420，以判断rul
t
是否小于一维修缓冲时间(buffer time；bt)，而获得一第
二结果，其中bt是由生产机台组件的原厂提供，当生产机台组件异常时，必须在此维修缓冲时间(bt)对生产机台组件进行维修或更换。当第一结果和第二结果均为否时，生产机台组件处于生病状态但未急速恶化，不需进行维修，而显示例如绿灯。当第一结果为否而第二结果为是时，生产机台组件未急速恶化但其剩余使用寿命不足，需进行维修，而显示例如蓝灯。当第一结果为是而第二结果为否时，生产机台组件急速恶化，而显示例如褐灯。若生产机台组件处理连续第t+i个工件的每一者均显示例如褐灯时，则需检查或维修生产机台组件，其中i为正整数。当第一结果和第二结果均为是时，生产机台组件需进行维修，而显示例如红灯。
[0144]
以下以组件的轴偏度(老化特征sx1与辅助老化特征sy1)、rul预测结果以及dci预测结果来说明本发明实施例。轴偏度代表组件的轴所偏离的角度，老化特征sx1代表组件的轴于x平面所偏离的角度，辅助老化特征sy1代表组件的轴于y平面所偏离的角度，其中辅助老化特征sy1依据前述多变量建模步骤260的格兰杰因果关系检验所检定出来的结果。请参照图5、6、7、8，其中图5为绘示组件的轴偏度的预测结果；图6为绘示组件的rul的预测结果；图7为绘示组件的一dci(tsp)的预测结果；图8为绘示组件的另一dci(tsp
mva
)的预测结果。其中tsp代表已知的单变量时间序列预测演算法(仅考量单个老化特征)；tsp
mva
代表本发明的多变量时间序列预测演算法(考量多个老化特征)。在图5中，横轴为时间，纵轴为轴偏度，曲线510为使用本发明tsp
mva
所获得的组件的轴偏度的预测值；曲线520为使用已知tsp所获得的组件的轴偏度的预测值；点群530(由“*”所组成)为实际值；直线540为组件的轴偏度的死亡规格值；直线550为组件的轴偏度的生病规格值。在图6中，横轴为工件号码，纵轴为rul(天数)，曲线560与曲线570分别为使用已知tsp与本发明tsp
mva
算出的rul。曲线560于工件号码中第106个工件出现预警，曲线570于工件号码中第119个工件出现预警。在图7与图8中，直线580为dci门槛值(dci
t
)。当dci小于dci门槛值时，不需进行维修；当dci大于dci门槛值时，代表组件在处理工件时接近死亡状态，需要维修。图7是用已知tsp算出的dci，图8是用本发明tsp
mva
算出的dci。已知tsp在距离死亡时间为24天前发出预警，而本发明tsp
mva
在距离死亡时间为11天前才发出预警。两者相较之下可知，已知tsp过早发出预警，因此本发明tsp
mva
评估的rul比较准确。
[0145]
可理解的是，本发明的生产机台组件的多变量预测保养方法为以上所述的实施步骤，本发明的内储用于量测抽样的计算机程序产品，是用以完成如上述的量测抽样的方法。上述实施例所说明的各实施步骤的次序可依实际需要而调动、结合或省略。上述实施例可利用计算机程序产品来实现，其可包含储存有多个指令的机器可读取媒体，这些指令可程序化(programming)计算机来进行上述实施例中的步骤。机器可读取媒体可为但不限定于软盘、光盘、只读光盘、磁光盘、只读记忆体、随机存取记忆体、可抹除可程序只读记忆体(eprom)、电子可抹除可程序只读记忆体(eeprom)、光卡(optical card)或磁卡、快闪记忆体、或任何适于储存电子指令的机器可读取媒体。再者，本发明的实施例也可做为计算机程序产品来下载，其可通过使用通讯连接(例如网络连线之类的连接)的数据信号来从远端计算机转移本发明的计算机程序产品至请求计算机。
[0146]
亦可注意的是，本发明亦可描述于一制造系统的语境中。虽然本发明可建置在半导体制作中，但本发明并不限于半导体制作，亦可应用至其他制造工业。制造系统系配置以制造工件或产品，工件或产品包含但不受限于微处理器、记忆体装置、数字信号处理器、专
门应用的电路(asics)或其他类似装置。本发明亦可应用至除半导体装置外的其他工件或产品，如车辆轮框、螺丝。制造系统包含一或多个处理工具，其可用以形成一或多个产品或产品的一部分，在工件(例如：晶圆、玻璃基板)上或中。发明本领域具有通常技艺者应可知，处理工具可为任何数目和任何型式，包含有微影机台、沉积机台、蚀刻机台、研磨机台、退火机台、工具机和类似工具。在实施例中，制造系统亦包含散射仪、椭圆偏光仪、扫描式电子显微镜和类似仪器。
[0147]
综上所述，本发明实施例可利用老化特征与辅助老化特征来预测老化特征，这意味着用多个变数来预测一个变数。已知技术只采用一个变数来做预测，则预测准确度有限，而本发明的多变量的时间序列预测扩充了只能使用一个变数的限制，所预测的老化特征不仅取决于过去的老化特征，也依赖其他特征(辅助老化特征)，因此采用多变数可以提升预测准确度。本发明实施例也可以即时并准确地预测生产机台组件的rul，而可及时地进行生产机台组件的维修；并可在生产机台组件可能很快死亡状态时立即进行维护，且以数据化的方式来呈现生产机台组件进入死亡状态的可能性。
[0148]
虽然本发明已以实施方式揭露如上，然其并非用以限定本发明，任何熟悉此技艺者，在不脱离本发明的精神和范围内，当可作各种的更动与润饰，因此本发明的保护范围当视所附的权利要求书所界定的范围为准。

技术特征：
1.一种生产机台组件的多变量预测保养方法，其特征在于包含：获得一生产机台的一组件依序处理多个工件时所使用的多组制程数据，其中每一该些组制程数据包含多个参数的数值；获得分别对应至该些组制程数据的多个事件指示值，其中该些事件指示值指出当该组件在处理每一该些工件时该组件是否发生异常事件；分别使用多个演算法将每一该些组制程数据的该些参数的数值转换成多个参数指标的数值；对该些组制程数据中的每一该些参数指标与该些事件指示值进行一相关性分析，而获得分别对应至该些参数指标的多个相关系数；选取对应至该些相关系数中一最大者的参数指标为一老化特征，并设定一辅助老化特征；进行一第一判断步骤，以根据每一该些工件对应的该老化特征的数值，来判断该组件在处理该些工件时是否处于一生病状态，其中一旦该组件在处理该些工件的一者时是处于该生病状态时，则将该些工件的该者设定为一样本选取点；进行一多变量建模步骤，该多变量建模步骤包含：使用该样本选取点前n个工件所对应的n组制程数据中对应至该老化特征的n个数值以及对应该辅助老化特征的n个数值为一组建模样本数据，其中n为正整数；对该老化特征与该辅助老化特征执行一格兰杰因果关系检验(granger causality test)以判断该老化特征与该辅助老化特征之间的相关性，若不相关则从该组建模样本数据中删除该辅助老化特征；以使用该组建模样本数据并根据一多变量时间序列预测演算法来建立一老化特征预测模型，而获得依工件生产次序排列的该老化特征的多个预测数值；及使用每一该些工件的一工件处理时间和该组件无法使用时的该老化特征的一死亡规格值，来将该些预测数值转换成的多个剩余使用寿命(remaining useful life；rul)预测值(rul
t
)，其中t代表第t个工件，t为整数；以及进行一第二判断步骤，以根据该些剩余使用寿命预测值(rul
t
)来判断该组件是否需要更换或维修。2.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该第一判断步骤包含：以一组转换公式分别将每一该些组制程数据的该老化特征的数值转换成分别对应至该些工件的多个装置健康指数(device health index；dhi)，该组转换公式为：当当当当
当当其中为该些组制程数据的该老化特征的数值的平均值，其中为所对应的转换值；max y
t
为该些组制程数据的该老化特征的最大数值，max y
t_mapping
为max y
t
所对应的转换值；min y
t
为该些组制程数据的该老化特征的最小数值，min y
t_mapping
为min y
t
所对应的转换值；lsl为规格下限；lcl为管制下限；ucl为管制上限；usl为规格上限；lsl
_mapping
为lsl所对应的转换值；lcl
_mapping
为lcl所对应的转换值；ucl
_mapping
为ucl所对应的转换值；usl
_mapping
为usl所对应的转换值；以及依序判断该些装置健康指数是否大于或等于一门槛值，并将该些装置健康指数中最先大于或等于一门槛值的一者所应的工件设定为该样本选取点。3.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该多变量建模步骤还包含：使用一向量自回归模型(vector autoregression model；var)为该多变量时间序列预测演算法，来建立该老化特征预测模型；使用一偏自相关函数(partial autocorrelation function；pacf)选出该向量自回归模型的最大落后期数；对该组建模样本数据中的该些数值进行一白噪音检定，其中当该些数值是白噪音时，则加入该样本选取点前n+1个工件所对应的又一组制程数据中对应至该老化特征的数值至该组建模样本数据；使用该向量自回归模型的最大落后期数，来建立多个向量自回归模型组合；使用一讯息准则演算法，来计算出每一该些向量自回归模型组合的讯息量；以及选出该些向量自回归模型组合中具有最大讯息量的一者为一最佳模型。4.根据权利要求3所述的生产机台组件的多变量预测保养方法，其特征在于该讯息准则演算法为贝氏讯息准则(bayesian information criteria；bic)。5.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该多变量建模步骤还包含：判断该组建模样本中的该些数值的变异数是否会随着时间而越来越大，其中当该些数值的变异数随着时间而越来越大时，对该组建模样本数据的每一该些数值进行对数转换；及对该些数值进行一单根检定(unit root test)，以确认依序排列的该些数值是否为稳态状态，其中当该些数值不是稳态状态时，对该组建模样本数据的每一该些数值进行差分转换。6.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该第二判断步骤包含：
判断(rul
t-rul
t-1
)/rul
t-1
是否大于或等于一门槛值，而获得一第一结果，其中t-1代表第t-1个工件；判断rul
t
是否小于一维修缓冲时间，而获得一第二结果，其中该组件必须在该维修缓冲时间进行维修；当该第一结果和该第二结果均为否时，该组件处于生病状态但未急速恶化，不需进行维修；当该第一结果为否而该第二结果为是时，该组件未急速恶化但其剩余使用寿命不足，需进行维修；当该第一结果为是而该第二结果为否时，该组件急速恶化，若处理连续第t+i个工件的每一者的该第一结果为是而该第二结果为否，则需检查或维修该组件，其中i为正整数；以及当该第一结果和该第二结果均为是时，该组件需进行维修。7.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该第二判断步骤包含：以一组转换公式分别将每一该些组制程数据的该老化特征的数值转换成分别对应至该些工件的该组件的多个死亡相关指数(death correlation index；dci)，该组转换公式为：其中y
death
为该组件在死亡状态时所对应的该老化特征的数值，y
t-1
为该组件在处理第t-1个工件时所对应的该老化特征的数值，conv为共变异数计算，var为变异数计算；当dci
t
大于一门槛值时，代表该组件在处理第t个工件时接近死亡状态，其中该门槛值的计算是根据dci
t
的标准差。8.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该组件为一加热器、一压力模块、一节流阀、一无油衬套或一轴承，该些参数包含：一轴偏度、一阀开度、一振动振幅、一驱动电压、一驱动电流、一温度和一压力。9.根据权利要求1所述的生产机台组件的多变量预测保养方法，其特征在于该些参数指标包含：一转换至频域后的k倍频(其中k大于0)、一整体相似度指标(global similarity index；gsi)、一统计数据分布的峰度(kurtosis)、一统计数据分布的偏度(skewness)、一标准差、一均方根(root mean square)、一平均值、一最大值和一最小值。10.一种计算机程序产品，其特征在于，当计算机载入此计算机程序产品并执行后，能够完成如权利要求1至9中任一项所述的生产机台组件的多变量预测保养方法。

技术总结
本发明实施例提供一种生产机台组件的多变量预测保养方法与其计算机程序产品，其使用多变量时间序列预测演算法(Multiple-VAriable Time SeriesPrediction；TSP


技术研发人员：林晋逸 谢昱铭 郑芳田 黄宪成
受保护的技术使用者：郑芳田
技术研发日：2022.05.07
技术公布日：2023/10/8