一种基于D-最优内表设计的田口试验设计方法

摘要

本发明提供一种基于D-最优内表设计的田口试验设计方法，该方法的具体步骤是：1.根据试验目的、条件及工程经验，给出试验响应、因子、因子间约束等基本试验信息；2.根据噪声因子的个数以及拟采用的水平个数，选择合适的均匀设计表；3.确定响应与因子之间的关系模型以及内表试验次数；4.在设计区域内，采用D-最优设计方法给出给定试验次数的内表设计方案，完成内表设计；5.采用基于经验分布函数等分位点方法，完成外表设计；6.根据试验方案进行试验，根据试验得到性能指标计算质量特性信噪比，计算得到使信噪比最大化的可控因子最佳水平；7.基于灵敏度分析法确定稳定因子最佳水平；8.变动调整因子，确定最佳因子水平组合。

说明书

技术领域

本发明涉及一种基于D-最优内表设计(即最以大化试验点行列式 (Determinant)为准则的设计方法)的田口试验设计方法，不仅可以有效解决传 统田口设计难以处理的设计区域不规则的情况，而且可有效缓解田口设计采用直 积表形式导致试验设计次数剧增的问题，适用于产品设计、生产制造、质量控制 以及工艺优化等相关技术领域。

背景技术

一个合理的试验设计方案可以在有限的时间、成本等资源条件下，充分揭示 各个影响因子对目标响应的影响，在有限的试验次数下获得最优的目标响应结果， 为提高产品质量寻求最优的工艺参数组合等。

为了提高产品性能的稳定性，田口玄一博士提出了田口设计方法，采用内表、 外表直积的形式，在内表中安排可控因子的各个水平，在外表中安排噪声因子进 行试验，而所谓直积法，是指对内表中可控因子的每个搭配，使用外表的误差因 子模拟各种干扰，计算该搭配的抗干扰能力，即信噪比(SN比)，从而通过对内表 设计各种方案的比较，利用SN比寻找最佳的可控因子搭配。

广义上说，所有的试验设计本身都是资源受限的设计，试验设计的目的就是 在尽量节约资源的情况下，获取尽量多的系统或者过程的信息，优化设计和工艺。 在传统的田口设计中，内表、外表均使用正交设计进行，通过方差分析进行数据 分析，其设计区域要求是规则的超立方体，一般不考虑因子之间存在约束而导致 设计区域不规则的情况。在实际工程中，由于影响因子之间经常是不独立的、存 在一定的相互作用，从而产生设计区域不规则的情况，这时，使用正交设计进行 内表设计，需切割或填补试验区域进行反复的折衷、对比，不仅设计过程复杂， 而且由于降维过多，导致设计效率较低。同时，田口设计的直积表形式，使得在 噪声因子较多的情况下，试验次数剧增，经常超出试验资源允许的范围，导致企 业难以承受；另外，由于传统的田口内外表设计使用正交设计进行，而正交设计 是基于设计表的设计，导致田口设计试验次数不能根据实际需求灵活调整，容易 造成资源不能充分利用或资源不足的情况。因此，研究在影响因子存在约束的情 况下，合理地选择试验点，从而科学地收集和分析数据，获得更优的参数组合， 优化设计和工艺，具有重要的理论意义和迫切的现实需求。

为此，本发明给出了一种基于D-最优内表设计的田口试验设计方法。

发明内容

(1)本发明的目的：本发明针对传统的田口设计在有限的试验次数内，难以 灵活高效解决试验区域不规则时的试验设计问题，提供一种基于D-最优内表设计 的田口试验设计方法，以筛选出更具有代表性的试验点，从而进行精确、高效、 灵活的试验方案设计。

D-最优设计是一种基于数学模型的设计方法，可以根据模型参数个数，在给定 的设计区域内，给出满足一定条件范围内任意试验次数的最优设计；本发明在田 口设计的框架下，利用D-最优设计进行内表设计，可以在因子间存在约束的情况 下，给出更加灵活、高效的试验设计方案，建立更加贴近产品加工或生产过程实 际的数学模型，进行更准确的数据分析和预测，从而确定最佳的因子水平组合。

(2)技术方案：

本发明给出了一种基于D-最优内表设计的田口试验设计方法。

田口设计运用内表和外表相结合的方法进行稳健设计，并使用信噪比(SN) 作为质量评价指标，以寻求最稳定的参数组合。

传统的田口设计，内表采用正交设计安排可控因子的试验，外表采用正交设 计安排噪声因子的试验，使用信噪比(SN)度量指标波动的大小。由于正交设计 的局限性，当设计区域不规则时，传统的田口设计方法无法给出高效的设计方案。 另外，正交试验次数是由设计方法确定的，而非试验精度，在田口设计直积表的 作用下，正交设计“整齐可比”的性质使得试验次数剧增，并且调整非常不灵活。 然而，在信噪比计算过程中，并不要求试验方案中各误差因素“整齐可比”，因 此，可以放松“整齐可比”的要求，探索更为有效的试验设计方法。

本发明采用D-最优设计进行内表设计，采用基于经验分布函数等分位点的均 匀设计进行外表设计。由一般性等价定理，对于有p个参数的线性模型，存在由 n^*(p≤n^*≤p(p+1)/2)个试验点组成的D-最优设计ξ^*。即D-最优设计可用在一定范 围内，给出任意设计区域、任意试验次数的内表试验设计方案，且随着试验因子 个数的增加，D-最优设计所需要的最小试验次数将远远小于其他设计。在满足试 验精度的前提下，使用D-最优设计进行内表设计，并结合均匀设计进行外表设计， 不仅可以有效解决传统田口设计难以处理的因子间存在约束的情况，而且可以在 较少的试验次数下，给出高效的试验设计方案。

基于上述理论和思路，本发明一种基于D-最优内表设计的田口试验设计方法， 具体实施步骤如下：

步骤一：根据试验目的、条件以及工程经验，给出试验基本试验信息，包括试 验响应(即所关注的质量特性)、影响因子(包括可控因子和噪声因子)及其取值 范围，各因子之间的约束关系、以及试验资源允许的试验总次数n^total。

步骤二：根据噪声因子的个数以及拟采用的水平个数m，选择合适的均匀设计 表，并以此确定外表试验次数n^outer。选择方法简介如下：

假设有噪声因子l个，拟采用水平数m，则从均匀设计表库中选择均匀表 同时外表试验次数确定为n^outer＝ω。

步骤三：确定响应与因子之间的关系模型以及内表试验次数n^inner。

该步骤中响应与因子之间的关系模型具体确定方法如下：

基于步骤一中的试验基本信息，采用数学模型的形式表达出响应与因子之间的 关系模型，记为y＝Xβ+ε。其中,为已知的n^inner×p维参数模 型矩阵；f(x_i)为关于x_i的已知函数，反应全部可控因子及其之间的约束关系； β＝(β₀,β₁,β₂,...,β_p)^T为p个待估参数，反应各个因子与响应之间的影响关系。

内表试验次数n^inner的确定方法如下：

由一般性等价定理，对于有p个参数的线性模型，存在由n^*(p≤n^*≤p(p+1)/2)个 试验点组成的D-最优设计ξ^*。所以，根据总试验次数n^total与外表试验次数n^outer，不 难给出D-最优内表试验次数范围[p,min(p(p+1)/2,n^total/n^outer)]。

因此，根据实际试验条件，选择任意n^inner∈[p,min(p(p+1)/2,n^total/n^outer)]即可。

步骤四：在设计区域内，采用D-最优设计方法给出给定试验次数n^inner的内表设 计方案，完成内表设计，其具体设计方法如下：

首先，将各因子取值范围规范至[-1，1]区间范围内，根据实际因子间的约束 关系，确定设计区域；然后，在该设计区域内，选择满足给定次数的设计ξ，使得 信息矩阵行列式最大，至此完成D-最优内表设计。进一步， 使用G-最优效率衡量该设计的优良性。G-最优效率的定义为：

$G_{eff}=\frac{p}{\underset{x\in X}{max}d(x,\xi )}$

这里，G-最优效率值G_eff越接近于1，表示设计ξ越好。

步骤五：采用基于经验分布函数等分位点方法，针对内表中每个外表试验点， 确定噪声因子水平，完成外表设计。

该步骤中噪声因子水平的确定方法如下：

记x_[i]为变量x的基于经验分布函数第i等分位点，其定义如下：

其中，F(x)为随机变量x的分布函数，m为噪声因子水平数；argF(·)表示函数F(·) 的反函数。换言之，上式表示，当时，所对应的随机变量x的值记为x_[i]，x_[i]即 为因子x的第i个水平值。

步骤六：根据试验方案进行具体试验，根据试验得到的产品性能指标计算质量 特性信噪比，拟合得到信噪比与因子之间的回归函数η(x_i)，计算得到使信噪比最 大化的可控因子最佳水平。

步骤七:基于灵敏度分析方法确定稳定因子及其最佳水平。

信噪比对影响因子x_i,i＝1,2,...,n的灵敏度定义如下：

S_ij的数值表示对变量x_i的敏感程度，S_ij的符号表示对变量x_i的单调 性。数值较大的S_ij对应的因子为稳定因子，取其最优水平组合以达到系统的稳定 性，使系统输出的波动减少到最低程度。其余因子为调整因子，调整系统的输出 值使之达到或接近于目标值，从而在尽量减少波动的情况下，得到使得质量特性 最接近其目标值的最佳参数组合。

步骤八:回归拟合求得关系模型y中未知参数，并依此对调整因子进行调整， 确定最佳因子水平组合，使得在该试验方案下质量特性估计值最接近目标值。

(3)优点和功效：

本发明提供一种基于D-最优内表设计的田口试验设计方法，其优点为：

①本发利用D-最优设计进行内表设计，可以在试验因子之间存在约束条件 下，给出高效合理的田口试验设计方案。

②本发明有效缓解了田口设计直积表形式导致的试验次数剧增问题，并且 内表可以给出满足一定条件下任意次数的试验，使得试验方案次数可以灵活调整。

附图说明

图1是本发明方法流程图。

具体实施方式

下面以小球氮化过程为例，结合附图，对本发明做进一步详细说明。

反馈弹簧是电液伺服阀的关键部件，其反馈杆上的小球与阀芯上的反馈槽要 求处于零间隙配合状态。在反馈弹簧的使用过程中，若小球硬度表面不够，容易 引发磨损现象，造成空载流量曲线突变，从而导致错误伺服信号，使电液伺服阀 产生致命性故障。为了提高小球表面硬度，需要对小球进行离子软氮化处理。氮 化层深度越深，则小球的耐磨性越好，但同时氮化层深度受到小球直径变化的制 约。为了保证小球的直径在允许范围内的同时，提高小球表面硬度，设定氮化层 深度目标值为0.05mm。见图1，基于本发明具体实施步骤如下：

步骤一：确定试验响应(即所关注的质量特性)、影响因子(包括可控因子和 噪声因子)及其取值范围，各因子之间的约束关系、以及试验资源允许的试验总 次数n^total。这里，考虑试验资源和进度要求，总试验次数要求不超过150次。

将氮化层的深度作为响应，在离子软氮化的过程中，每次都要经历一个溅射阶 段和两个氮化阶段，表1汇总了各个影响因子及其取值范围：

表1影响因子及其取值范围

根据离子软氮化理论，为保证氮化效果，设定约束条件：x₅-x₂≥20℃， x₉-x₅≥10℃。另外，Ar,H₂,N₂通气量(L/min)难以精确控制，存在波动，为噪声因 子，其容差分别为${\Delta }_{x_3}=0.08L/min,{\Delta }_{x_6}=0.1L/\min ,{\Delta }_{x_7}=0.02L/min.$

步骤二：根据噪声因子的个数以及拟采用的水平个数m，选择合适的均匀设计 表，并以此确定外表试验次数n^outer。选择方法简介如下：

假设有噪声因子l个，拟采用水平数m，则从均匀设计表库中，依据偏差选择均 匀表或且外表试验次数确定为n^outer＝ω。

假设有噪声因子l个，拟采用水平数m，则从均匀设计表库中选择均匀表 U_ω(m^l)，且外表试验次数确定为n^outer＝ω。

在本例中，共有3个噪声因子，为在保证精度的前提下减少试验次数，对每个 因子选取6个水平进行试验，选用均匀设计表进行外表设计，外表试验次数 n^outer＝6。

步骤三：确定响应与因子之间的关系模型以及内表试验次数n^inner。

该步骤中响应与因子之间的关系模型具体确定方法如下：

基于步骤一中的试验基本信息，采用数学模型的形式表达出响应与因子之间的 关系模型，记为y＝Xβ+ε。其中,为已知的n^inner×p维参数模 型矩阵；f(x_i)为关于x_i的已知函数，反应全部可控因子及其之间的约束关系； β＝(β₀,β₁,β₂,...,β_p)^T为p个待估参数，反应各个因子与响应之间的影响关系。

内表试验次数n^inner的确定方法如下：

由一般性等价定理，对于有p个参数的线性模型，存在由n^*(p≤n^*≤p(p+1)/2)个 试验点组成的D-最优设计ξ^*。所以，根据总试验次数n^total与外表试验次数n^outer，不 难给出D-最优内表试验次数范围[p,min(p(p+1)/2,n^total/n^outer)]。因此，根据实际试验 条件，选择任意n^inner∈[p,min(p(p+1)/2,n^total/n^outer)]即可。

在本实施例中，根据工程经验，选取如下模型：

y＝β₀+x₁β₁+x₂β₂+xβ₃+x₄β₄+x₅β₅+x₆β₆+x₇β₇+x₈β₈+x₉β₉+x₂x₃β₂₃+

x₂x₅β₂₅+x₂x₈β₂₈+x₃x₅β₃₅+x₃x₆β₃₆+x₃x₇β₃₇+x₃x₈β₃₈+x₃x₉β₃₉+x₆x₉β₆₉+e

其中，误差e～N(0,σ²)。

该模型中共有19个未知参数，即p＝19。结合内表所需次数，可确定外表试验 次数的取值范围19≤n^inner≤min(19(19+1)/2,n^total/6)。考虑总试验资源，选用内表试验 次数n^inner＝19。

步骤四：在设计区域内，采用D-最优设计方法给出给定试验次数n^inner的内 表设计方案，完成内表设计。具体设计方法如下：

首先，将各因子取值范围规范至[-1，1]区间范围内，根据实际因子间的约 束关系，确定设计区域；然后，在该设计区域内，选择满足给定次数的设计ξ，使 得信息矩阵行列式最大，其中，f(x_i)为关于x_i的已知函数， 至此完成D-最优内表设计。进一步，使用G-最优效率衡量该设计的优良性。G-最 优效率的定义为：

$G_{eff}=\frac{p}{\underset{x\in X}{max}d(x,\xi )}$

其中，d(x,ξ)＝f^T(x)M^-1(ξ)f(x)。G-最优效率值G_eff越接近于1，则表示设计ξ越好。

在本实施例中，D-最优内表设计表如下所示：

表2基于D-最优设计的内表试验方案(编码)

序号 x₁x₂x₃x₄x₅x₆x₇x₈x₉1 -1 -1 -1 -1 -1 -1 1 -0.38418 1 2 -1 -0.12675 1 -1 0.66667 1 1 -1 1 3 0.02211 -0.19135 0.001343 0.021786 -0.1083 -0.01131 0.050151 -0.03109 1 4 1 1 0.022368 1 0.66667 -0.00126 1 -0.01406 1 5 -1 1 1 -1 -0.33333 -1 -1 1 1 6 1 0.31204 0.036341 -1 -1 1 1 -1 -0.66667 7 -1 1 0.073784 1 0.66667 1 -1 1 1 8 1 -1 -1 -0.15154 -1 1 1 0.025892 1 9 0.072806 -1 -1 1 -1 -0.03687 -1 -1 -0.66667 10 1 0.18165 1 1 -1 -1 1 -1 -0.66667 11 -1 1 -1 -1 0.66667 -1 -0.10588 -0.045 1 12 -1 -1 0.020447 1 -0.09342 -1 -1 -1 0.23992 13 -1 -0.03265 1 1 -1 1 1 1 1 14 1 -1 1 -1 -1 -1 1 -0.00697 0.27705 15 1 -1 -0.26793 -1 -1 -1 -1 -1 1 16 -1 -1 1 1 0.66667 -1 -1 1 1

17 1 0.020626 0.074386 1 -1 -1 -1 0.055004 1 18 -1 0.33333 1 -0.04891 -1 1 1 1 -0.66667 19 -0.08383 1 1 0.048094 0.66667 -1 -1 -1 1

该设计的G效率可见其设计效率较高。

步骤五：采用基于经验分布函数等分位点方法，针对内表中每个外表试验点， 确定噪声因子水平，完成外表设计。

该步骤中噪声因子水平的确定方法如下：

记x_[i]为变量x的基于经验分布函数第i等分位点，其定义如下：

其中，F(x)为随机变量x的分布函数，m为噪声因子水平数；argF(·)表示函数F(·) 的反函数。换言之，上式表示，当时，所对应的随机变量x的值记为x_[i]，x_[i]即 为因子x的第i个水平值。

本例中，假定噪声因子均服从正态分布N(0,σ²)，容差Δ与标准差σ之间存在关 系Δ＝3σ，可得三个噪声因子的方差分别为： ${\sigma }_{x_3}^2={(0.08/3)}^2,{\sigma }_{x_6}^2={(0.1/3)}^2,{\sigma }_{x_7}^2={(0.02/3)}^2,$均值μ则根据相应的内表确定。以表2 中2号内表为例，该条件下确定的三个噪声因子水平数见表3。

表32号内表对应噪声因子水平数设定

则对应的外表设计为：

表42号内表对应的外表设计

步骤六：根据试验方案进行具体试验，根据试验得到的产品性能指标计算 质量特性信噪比，拟合得到信噪比与因子之间的回归函数η(x_i)，计算得到使信噪 比最大化的可控因子最佳水平。

该例中，所关注的响应氮化层深度呈现望目特性，根据所得试验数据，计 算各个条件下的信噪比见表5

表5试验设计方案(解码)及信噪比

利用上述数据，对信噪比η与各因子进行模型参数估计，采用逐步回归法可 得到以下模型

η＝3465.311+0.033862x₁-1.34615x₂+16.60486x₃-0.1968x₄-6.99205x₅-384.59x₆

+84.20613x₇+0.3755x₈-5.12495x₉+0.148254x₁x₃+0.002458x₂x₅+0.328684x₄x₆

+0.255114x₄x₇+2.575394x₆x₈+0.88726x₇x₈+0.010374x₅x₉

利用非线性求解法，求出在满足试验约束的情况下使得信噪比最大的因子 水平组合为：

x^*＝[30,500,0.6,120,520,0.07,0.475,150,530]

步骤七:基于灵敏度分析方法确定稳定因子及其最佳水平。

信噪比对影响因子x_i,i＝1,2,...,n的灵敏度定义如下：

S_ij的数值表示对变量x_i的敏感程度，S_ij的符号表示对变量x_i的单 调性。数值较大的S_ij对应的因子为稳定因子，取其最优水平组合以达到系统的稳 定性，使系统输出的波动减少到最低程度。其余因子为调整因子，调整系统的输 出值使之达到或接近于目标值，从而在尽量减少波动的情况下，得到使得质量特 性最接近其目标值的最佳参数组合。

本实施例中，目标函数η对各个影响因子的灵敏度计算如下：

目标函数η对所有影响因子的灵敏度自大至小汇总为：

表6η对各影响因子灵敏度汇总表

显然，目标函数η对N₂通气量最为敏感，其次是H₂通气量和Ar通气量。因 此，确定N₂、H₂Ar通气量三种因子为稳定因子，选取其对信噪比的最佳水平； 确定三个阶段通气时间和通气温度为调整因子，利用调整因子将设计方案的质量 特性调整到目标值。

步骤八:回归拟合求得关系模型y中未知参数，并依此对调整因子进行调 整，确定最佳因子水平组合，使得在该试验方案下质量特性估计值最接近目标值。

本实施例中，利用最小二乘法，求得关系模型y中未知参数值为：

β＝(β₀,β₁,...,β₆₉)

＝(1.9666,0.0003,-0.086,-0.7368,-0.0002,-0.001,-0.0019,0.1384,0.78,

-0.0083,0,0572,-0.0006,0.0028,-0.0006,0.0011,-0.1771,0.001,0.003,0.231)

根据关系模型，调动调整因子，确定最终的最优方案为：

x^*＝[30,440,0.6,120,350,0.07,0.475,120,500]

该参数组合下，小球氮化层深度预测值为5.01mm，达到目标值要求，圆满完成 了小球表面氮化试验设计的目标。