Minitab实验之试验设计
实验目的：
Minitab_专业论坛
www.mimitab.cn
本实验主要引导学生利用Minitab统计软件进行试验设计分析，包括全因子设计、部分因子设计、响应曲面设计、混料设计、田口设计以及响应优化，并能够对结果做出解释。
实验仪器：Minitab软件、计算机
实验原理：
“全因子试验设计”（fullfactorialdesign）的定义是：所有因子的所有水平的所有组合都至少要进行一次试验的设计。由于包含了所有的组合，全因子试验所需试验的总次数会比较多，但它的优点是可以估计出所有的主效应和所有的各阶交互效应。所以在因子个数不太多，而且确实需要考察较多的交互作用时，常常选用全因子设计。一般情况下，当因子水平超过2时，由于试验次数随着因子个数的增长而呈现指数速度增长，因而通常只作2水平的全因子试验。
进行2水平全因子设计时，全因子试验的总试验次数将随着因子个数的增加而急剧增加，例如，6个因子就需要64次试验。但是仔细分析所获得的结果可以看出，建立的6因子回归方程包括下列一些项：常数项、主效应项有6项、二阶交互作用项15项、三阶交互项20项，…，6阶交互项1项，除了常数项、主效应项和二阶交互项以外，共有42项是3阶以及3阶以上的交互作用项，而这些项实际上已无具体的意义了。部分因子试验就是在这种思想下诞生的，它可以使用在因子个数较多，但只需要分析各因子和2阶交互效应是否显著，并不需要考虑高阶的交互效应，这使得试验次数大大减少。
在实际工作中，常常要研究响应变量Y是如何依赖于自变量，进而能找到自变量的设置使得响应变量得到最佳值（望大、望小或望目）。如果自变量的个数较少（通常不超过3个），则响应曲面方法（responsesurfacemethodology，RSM）是最好的方法之一，本方法特别适合于响应变量望大或望小的情形。通常的做法是：先用2水平因子试验的数据，拟合一个线性回归方程（可以包含交叉乘积项），如果发现有弯曲的趋势，则希望拟合一个含二次项的回归方程。其一般模型是（以两个自变量为例）：
这些项比因子设计的模型增加了各自的变量的平方项。由于要估计这些项的回归系数，原来因子设计所安排的一些设计点就不够用了，需要再增补一些试验点。这种先后分两阶段完成全部试验的策略就是“序贯试验”的策略。适用于这种策略的方法有很多种，其中最常用的就是中心复合设计（centralcompositedesign，CCD）。
稳健参数设计（robustparameterdesign）（也称健壮设计、鲁棒设计，简称参数设计）是工程实际问题中很有价值的统计方法。它通过选择可控因子的水平组合来减少一个系统对噪声变化的敏感性，从而达到减小此系统性能波动的目的。过程的输入变量有两类：可控因子和参数因子。可控因子是指一旦选定就保持不变的变量，它包括产品或生产过程设计中的设计参数，而噪声因子是在正常条件下难以控制的变量。在做参数设计时，就是把可控因子的设计当做研究的主要对象，与此同时让噪声因子按照设定的计划从而系统改变其水平的方法来表示正常条件下的变化，最终按照我们预定的望大、望小或望目地目标选出最佳设置。田口玄一博士在参数设计方法方面贡献非常突出，他在设计中引进信噪比的概念，并以此作为评价参数组合优劣的一种测度，因此很多文献和软件都把稳健参数设计方法称为田口方法（Taguchidesign）。
在实际工作中，常常需要研究一些配方配比试验问题。这种问题常出现在橡胶、化工、制药、冶金等课题中。例如不锈钢由铁、镍、铜和铬4种元素组成；闪光剂由镁、硝酸钠、硝酸锶及固定剂组成；复合燃料、复合塑料、混纺纤维、混泥土、粘结剂、药品、饲料等都是由多种成分按相应比例而不是其绝对数值；而且显然所有分量之和总是为1的。对于这种分量之和总是为1的试验设计，称为混料设计（mixturedesign）。
实验内容和步骤：
实验之一：全因子试验设计
:例：改进热处理工艺提高钢板断裂强度问题。合金钢板经热处理后将提高其断裂其抗断裂性能，但工艺参数的选择是个复杂的问题。我们希望考虑可能影响断裂强度的4个因子，确认哪些因子影响确实是显著的，进而确定出最佳工艺条件。这几个因子及其试验水平如下：
A：加热温度，低水平：820，高水平：860（摄氏度）
B：加热时间，低水平：2，高水平：3（分钟）
C：转换时间，低水平：1.4，高水平：1.6（分钟）
D：保温时间，低水平：50，高水平：60（分钟）
由于要细致考虑各因子及其交互作用，决定采用全因子试验，并在中心点处进行3次试验，一共19次试验。
步骤1：全因子设计的计划（创建）
选择[统计]=>[DOE]=>[因子]=>[创建因子设计]，单击打开创建因子设计对话框。
,
选择两水平因子（默认生成元），在因子数中选择4，单击“设计”选项，弹出“设计”选项对话框。选择“全因子”试验次数为16的那行，并在“每个区组的中心点数”中选择3，其他项保持默认（本例中没有分区组，各试验点皆不需要完全复制）。单击确定。
单击“因子”选项打开，分别填写四个因子的名称及相应的低水平和高水平的设置。单击确定。
“选项”选项可以使用折叠设计（这是一种减少混杂的方法）、指定部分（用于设计生成）、使设计随机化以及在工作表中存储设计等；“结果”选项用于控制会话窗口中显示的输出。本例中这两项保持默认。单击确定，计算机会自动对于试验顺序进行随机化，然后形成下列表格。在表的最后一列，写上响应变量名（强度），这就完成了全部试验的计划阶段的工作。
步骤2：拟合选定模型
按照上图的试验计划进行试验，将结果填入上表的最后一列，则可以得到试验的结果数据（数据文件：DOE_热处理（全因）），如下：
拟合选定模型的主要任务是根据整个试验的目的，选定一个数学模型。通常首先可以选定“全模型”，就是在模型中包含全部因子的主效应及全部因子的二阶交互效应。在经过细致的分析之后，如果发现某些主效应和二阶交互效应不显著，则在下次选定模型的时候，应该将不显著的主效应和二阶交互效应删除。
选择[统计]=>[DOE]=>[因子]=>[分析因子设计]，打开分析因子设计对话框。
点击“项”选项后，在“模型中包含项的阶数”中选择2（表示模型中只包含2阶交互作用和主效应项，三阶以上交互作用不考虑），对默认的“在模型中包括中心点”保持不选。单击确定。
在“图形”选项中，“效应图”中选择“正态”和“Pareto”，“图中的标准差”中选择“正规”，“残差图”中选择“四合一”，在“残差与变量”图中将“加热温度”、“加热时间”、“转换时间”和“保温时间”选入，单击确定。
在“存储”选项中，在“拟合值与残差”中选定“拟合值”和“残差”，在“模型信息”中选定“设计矩阵”。单击确定。
结果如下：
拟合因子:强度与加热温度,加热时间,转换时间,保温时间
强度的估计效应和系数（已编码单位）
系数标
项效应系数准误TP
常量541.6321.377393.390.000
加热温度20.03810.0191.5006.680.000
加热时间16.8878.4441.5005.630.000
转换时间3.8131.9061.5001.270.240
保温时间11.1135.5561.5003.700.006
加热温度加热时间0.7370.3691.5000.250.812
加热温度转换时间-0.487-0.2441.500-0.160.875
加热温度保温时间3.0621.5311.5001.020.337
加热时间转换时间1.2630.6311.5000.420.685
加热时间保温时间7.1133.5561.5002.370.045
转换时间保温时间0.8370.4191.5000.280.787
S=6.00146PRESS=1778.45
R-Sq=92.49%R-Sq（预测）=53.68%R-Sq（调整）=83.11%
强度的方差分析（已编码单位）
来源自由度SeqSSAdjSSAdjMSFP
主效应43298.853298.85824.7122.900.000
2因子交互作用6252.17252.1742.031.170.408
残差误差8288.14288.1436.02
弯曲19.929.929.920.250.633
失拟5169.72169.7233.940.630.709
纯误差2108.50108.5054.25
合计183839.16
强度的估计系数（使用未编码单位的数据）
项系数
常量932.26
加热温度-0.25063
加热时间-111.262
转换时间43.812
保温时间-16.5637
加热温度加热时间0.036875
加热温度转换时间-0.121875
加热温度保温时间0.0153125
加热时间转换时间12.6250
加热时间保温时间1.42250
转换时间保温时间0.83750
结果分析：
分析要点一：分析评估回归的显著性。包含三点：
（1）看方差分析表中的总效果。方差分析表中，主效应对应的概率P值为0.000小于显著性水平0.05，拒绝原假设，认为回归总效果是显著的。
（2）看方差分析表中的失拟现象。方差分析表中，失拟项的P值为0.709，无法拒绝原假设，认为回归方程并没有因为漏掉高阶交互作用项而产生失拟现象。
（3）看方差分析表中的弯曲项。方差分析表中，弯曲项对应的概率P值0.633，表明无法拒绝原假设，说明本模型中没有弯曲现象。
分析要点二：分析评估回归的总效果
（1）两个确定系数R-Sq与R-Sq（调整），计算结果显示，这两个值分别为92.49%和83.11%，二者的差距比较大，说明模型还有待改进的余地。
（2）对于预测结果的整体估计。计算结果显示R-Sq和R-Sq（预测）分别为92.49%和53.68%，二者差距比较大；残差误差的SSE为288.14，PRESS为1778.45，两者差距也比较大；说明在本例中，如果使用现在的模型，则有较多的点与模型差距较大，模型应该进一步改进。
分析要点三：分析评估各项效应的显著性。计算结果显示，4个主效应中，加热温度、加热时间和保温时间是显著的，只有转换时间不显著；6个2因子水平交互效应中，只有加热时间保温时间是显著的。说明本例中还有不显著的自变量和2因子交互作用，改进模型时应该将这些主效应和交互作用删除。
对于各项效应的显著性，计算机还输出了一些辅助图形来帮助我们判断和理解有关结论。
Pareto图是将各效应的t检验的t值的绝对值作为纵坐标，按照绝对值的大小排列起来，根据选定的显著性水平，给出t值的临界值，绝对值超过临界值的效应将被选中，说明这些效应是显著的。从图中可以看到，加热温度、加热时间、保温时间以及加热时间保温时间是显著的。
正态效应图，凡是因子效应离直线不远者，就表明这些效应是不显著的；反之，则是显著的。从图中可以看到，加热温度、加热时间、保温时间以及加热时间保温时间是显著的。
步骤3：残差诊断
残差诊断的主要目的是基于残差的状况来诊断模型是否与数据拟合得比较好。如果数据和模型拟合得比较好，则残差应该是正常的。残差分析包括四个步骤：
（1）在“四合一”图的右下角图中，观察残差对于以观测值顺序为横轴的散点图，重点考察此散点图中，各点是否随机地在水平轴上下无规则的波动着。
（2）在“四合一”图的右上角图中，观察残差对于以响应变量拟合预测值为横轴的散点图，重点考察此散点图中，残差是否保持等方差性，即是否有“漏斗型”或“喇叭型”。
（3）在“四合一”图的左上角正态概率图（或右下角的直方图）中，观察残差的正态检验图，看残差是否服从正态分布。
（4）观察残差对于以各自变量为横轴的散点图，重点观察此散点图中是否有弯曲趋势。
从上面这些图可以看到，这些图形都显示残差是正常的。
步骤4：判断模型是否需要改进
这一步需要综合前面的分析：包括残差诊断和显著性分析。从上面的分析我们得知，在模型中包含不显著项，应该予以删除，所以需要建立新的模型。
选择[统计]=>[DOE]=>[因子]=>[分析因子设计]，打开分析因子设计对话框。主要是修改“项”选项中的设置，在选取的项中将加热温度、加热时间和保温时间保留，其他项皆删去，操作中的其余各项都保持不变。单节确定。
结果如下：
拟合因子:强度与加热温度,加热时间,保温时间
强度的估计效应和系数（已编码单位）
系数标
项效应系数准误TP
常量541.3191.363397.270.000
加热温度20.03810.0191.3637.350.000
加热时间16.8878.4441.3636.200.000
保温时间11.1125.5561.3634.080.001
加热时间保温时间7.1133.5561.3632.610.022
CtPt1.9813.4290.580.573
S=5.45038PRESS=724.350
R-Sq=89.94%R-Sq（预测）=81.13%R-Sq（调整）=86.07%
强度的方差分析（已编码单位）
来源自由度SeqSSAdjSSAdjMSFP
主效应33240.713240.711080.2436.360.000
2因子交互作用1202.35202.35202.356.810.022
弯曲19.929.929.920.330.573
残差误差13386.19386.1929.71
失拟3151.52151.5250.512.150.157
纯误差10234.67234.6723.47
合计183839.16
强度的估计系数（使用未编码单位的数据）
项系数
常量212.788
加热温度0.500938
加热时间-61.3500
保温时间-2.44500
加热时间保温时间1.42250
CtPt1.98125
结果分析：
从方差分析表中可以看到，主效应和2阶交互作用对应的概率都小于显著性水平0.05，应该拒绝原假设，认为本，本模型总的来说是有效的；失拟值和弯曲对应的概率分别为0.157和0.573，都大于显著性水平，不应拒绝原假设，说明本模型删除了很多项之后，并没有造成失拟的现象。
再看删减后的模型是否比原来的有所改进。从上述表中，可以看到，由于模型的项数减少了6项，R-Sq通常都会有微小的降低（本例由0.9249降到0.8968），但关键还是要看调整的R-Sq（调整）是否有所提高，本例中，该值从0.8311提高到0.8673，可见删除不显著的效应之后，回归的效果明显好了；而s的值有6.00146降为5.31913，PRESS由1778.45降到704.408，再次证明删除不明显的主效应和交互效应后，回归的结果更好了。
步骤5：对选定的模型进行分析解释
经过前三步的多次反复以后，我们可以获得一个满意的回归方程：
对选定的模型进行分析，主要是在拟合选定模型后输出更多的图形和信息，并做出有意义的解释。主要包括下面四个方面：
（1）再次进行残差诊断。
具体做法是：选择[统计]=>[DOE]=>[因子]=>[分析因子设计]，打开分析因子设计对话框。点击“图形”窗口后，在“图中的残差”中选择“标准化”，在“残差图”中，在单独视图下选择“直方图”，单击确定。点击“存储”窗口后，在拟合值与残差中，选择“标准化残差”和“删后”。单击确定。
结果如下：
从得出的直方图可知，残差及所有残差数据都是正常的。
（2）确认主效应及交互作用的显著性，并考虑最优设置
通过输出各因子的主效应图和交互效应图来判定。具体做法是：选择[统计]=>[DOE]=>[因子]=>[因子图]，打开因子图对话框。选定“主效应图”和“交互作用图”，在图中使用的均值类型中选择“数据均值”。在主效应图的设置中，将“强度”选入到响应中，将可用中的所有项选入所选中；在交互作用图的设置中，重复前面主效应图设置的步骤。单击确定。
结果如下：
从主效应图中可以看到，加热温度、加热时间和保温时间三者的回归线比较陡，顾主效应影响确实显著，而转换时间的回归线较平，故主效应影响不显著；为了使断裂强度达到最大，三因子都是取值越大越好，即加热温度应取上限860摄氏度，加热时间应取上限3分钟，保温时间应取上限60分钟。从交互作用图可以看出，只有加热时间和保温时间二者效应线明显不平行，说明二者交互作用显著。
（3）输出等值线图、响应曲面图等以确认最佳设置
本例中，只有加热时间和保温时间的交互作用显著，因此绘制这组等值线图和响应曲面图，而设定另一个影响显著的变量（加热温度）为最佳设置。具体操作为：选择[统计]=>[DOE]=>[因子]=>[等值线/曲面图]，打开等值线/曲面图对话框。选定“等值线图”和“曲面图”。在等值线图设置中，在因子中，X轴选为加热时间，Y轴选为保温时间，在设置中，选择保留附加因子在高设置，并在加热时间中设置860，单击确定；在曲面图设置中，X轴中选择加热时间，Y轴中选择保温时间，单击确定。
结果如下：
从等值线图和曲面图可以看出，断裂强度的最大值确实在加热时间为3分钟，保温时间为60分钟，加热温度固定在860摄氏度时达到最大。
（4）实现最优化
Minitab软件中有专门的响应变量优化器窗口。具体做法：[统计]=>[DOE]=>[因子]=>[响应优化器]，打开响应优化器对话框。将“可用项”中的强度选入到“所选项”中；点击“设置”窗口，根据本例的要求，在“目标”中选择“望大”，在“下限”中填入560（这个值是在做过的试验中已经实现了的），在“望目”中填入600（这个值是在做过的试验中未能达到的，是较高理想），上限留为空白。
结果如下：
这个图中共有3列，分别为选中的自变量。最上端列出各变量的名称、取值范围以及最优设置，上半图是合意值d的取值情况，下半图是最优化结果：最大值在加热温度取860摄氏度、加热时间取3分钟、保温时间取60分钟达到，断裂强度最终可以达到569.2066。合意度d为0.23016。
步骤6：进行验证试验
通常的做法是在先算出在最佳点的观测值的预测值及其变动范围，然后再最佳点做若干次验证试验，如果验证试验结果的平均值落在事先计算好的范围内，则说明一切正常，模型是正确的，预测结果可信；否则就要进一步分析发生错误的原因，改进模型，再重新验证，以求得符合实际数据的统计模型。具体做法是：选择[统计]=>[DOE]=>[因子]=>[分析因子设计]，打开分析因子设计对话框。在前面建立的模型的基础上，即在“项”中已经将最终选定的模型中包括了加热温度、加热时间、保温时间以及加热时间和保温时间的交互作用项。再打开“预测”窗口，在“因子”中按顺序设定各个主效应的最优值，分别为860360。单击确定。
结果如下：
根据该模型在新设计点处对强度的预测响应
拟合值
点拟合值标准误95%置信区间95%预测区间
1569.2072.926(562.931,575.483)(556.186,582.227)
结果解释：最左侧给出的拟合预测值是569.207，就是将自变量值代入回归方程所得的结果，这与最优值的预测是一致的。拟合值标准误为2.926，是拟合值的标准差，此值在作进一步计算时还有用。预测值平均值置信区间的结果是（562.931,575.438），具体的理解可以是：当加热温度取860摄氏度，加热时间取3分钟，保温时间取60分钟时，我们有95%的把握断言，断裂强度平均值将落入（562.931,575.438）之内。95%的预测区间是将来一次验证试验时将要落入的范围，可供做验证试验时使用，具体的理解是：当加热温度取860摄氏度，加热时间取3分钟，保温时间取60分钟时，我们有95%的把握断言，任何一块钢板的断裂强度将落入（556.186,582.227）之内。
试验之二：部分因子试验设计
部分因子试验设计与全因子试验设计的不同之处在于大大减少了试验的次数，具体表现在试验设计创建阶段的不一致，下面主要就部分因子试验设计的创建进行讲述。
步骤1：部分因子试验的计划（创建）——默认生成元的计划
例：用自动刨床刨制工作台平面的工艺条件试验。在用刨床刨制工作台平面试验中，考察影响其工作台平面光洁度的因子，并求出使光洁度达到最高的工艺条件。
共考察6个因子：
A因子：进刀速度，低水平1.2，高水平1.4（单位：mm/刀）
B因子：切屑角度，低水平10，高水平12（单位：度）
C因子：吃刀深度，低水平0.6，高水平0.8（单位：mm）
D因子：刀后背角，低水平70，高水平76（单位：度）
E因子：刀前槽深度，低水平1.4，高水平1.6（单位：mm）
F因子：润滑油进给量，低水平6，高水平8（单位：毫升/分钟）
要求：连中心点在内，不超过20次试验，考察各因子主效应和2阶交互效应AB、AC、CF、DE是否显著。由于试验次数的限制，我们在因子点上只能做试验16次，另4次取中心点，这就是的试验，通过查部分因子试验分辨度表可知，可达分辨度为Ⅳ的设计。具体操作为：选择[统计]=>[DOE]=>[因子]=>[创建因子设计]，单击打开创建因子设计对话框。在“设计类型”中选择默认2水平因子（默认生成元），在“因子数”中选定6。
单击“显示可用设计”就可以看到下图的界面，可以确认：用16次试验能够达到分辨度为Ⅳ的设计。
单击“设计”选项，选定1/4部分实施，在每个区组的中心点数中设定为4，其他的不进行设定，单击确定。
单击“因子”选项，设定各个因子的名称，并设定高、低水平值。点击确定。
再点击确定后，就可以得到试验计划表，如下：
与全因子设计不同的是，我们不能肯定这个试验计划表一定能满足要求，因为部分因子试验中一定会出现混杂，这些混杂如果破坏了试验要求，则必须重新进行设计，从运行窗中可以看到下列结果：
设计生成元:E=ABC,F=BCD
别名结构
I+ABCE+ADEF+BCDF
A+BCE+DEF+ABCDF
B+ACE+CDF+ABDEF
C+ABE+BDF+ACDEF
D+AEF+BCF+ABCDE
E+ABC+ADF+BCDEF
F+ADE+BCD+ABCEF
AB+CE+ACDF+BDEF
AC+BE+ABDF+CDEF
AD+EF+ABCF+BCDE
AE+BC+DF+ABCDEF
AF+DE+ABCD+BCEF
BD+CF+ABEF+ACDE
BF+CD+ABDE+ACEF
ABD+ACF+BEF+CDE
ABF+ACD+BDE+CEF
从此表得知，计算机自己选择的生成元是：E=ABC,F=BCD。后面的别名结构中列出了交互作用项的混杂情况，即每列中互为别名的因子有哪些；从上表可以看出，主效应与三阶及四阶交互作用混杂，二阶交互作用与四阶交互作用混杂，三阶交互作用与四阶交互作用混杂；关键是要检查一下题目所要求的2阶交互作用情况，将3阶以上的交互作用忽略不计，混杂的情况有：AB=CE,AC=BE,AD=EF,AF=DE,AE=BC=DF,BD=CF,BF=CD。本例中所要求的4个2阶交互作用是AB,AC,CF,DE，显然可以看到，这四个2阶交互作用均没有混杂。因此可以看到此试验计划是可行的。
步骤2：指定生成元的部分因子试验计划
例：和前面的例子是一样的，考察的是各因子主效应和2阶交互效应AB,AC,CE和DE是否显著。
从上例的别名结构表中可以看出，AB与CE是相互混杂，因此用默认的生成元构造的试验计划是不能满足要求的。
指定生成元的步骤：由要求条件可知，AB,AC,CE和DE不能混杂，这相当于AB≠CE,AB≠DE,AC≠DE,运用移项法则，变形后可知，即E≠ABC,E≠ABD,E≠ACD.对于分辨度为Ⅳ的设计生成元中，只能含3个字母。而试验次数为16的的各列中，字母个数为3的项只有4个：ABC,ABD,ACD以及BCD。既然给定条件中有3个选择不可接受，因此，生成元只能选择E=BCD，试验计划对于F没有要求，因此F可以任选，取F=ABC。
具体操作为：选择[统计]=>[DOE]=>[因子]=>[创建因子设计]，单击打开创建因子设计对话框。在“设计类型”中选择2水平因子（指定生成元），在“因子数”中选定4（这是基本设计的因子数，其他两个因子是通过指定生成元加入的）。
打开“因子”对话框，选定全因子，并在“每个区组的中心点数”中选择4。打开“生成元”选项，在“通过列出生成元将因子添加到基本设计中”中填写生成元：E=BCDF=ABC,单击确定。
单击确定后，得到的结果如下：
设计生成元:E=BCD,F=ABC
别名结构（直到4阶项）
I+ABCF+ADEF+BCDE
A+BCF+DEF+ABCDE
B+ACF+CDE+ABDEF
C+ABF+BDE+ACDEF
D+AEF+BCE+ABCDF
E+ADF+BCD+ABCEF
F+ABC+ADE+BCDEF
AB+CF+ACDE+BDEF
AC+BF+ABDE+CDEF
AD+EF+ABCE+BCDF
AE+DF+ABCD+BCEF
AF+BC+DE+ABCDEF
BD+CE+ABEF+ACDF
BE+CD+ABDF+ACEF
ABD+ACE+BEF+CDF
ABE+ACD+BDF+CEF
从上面的结果可以看出，AB,AC,CE和DE均没有相互混杂，此设计满足原定的要求。
部分因子试验的分析步骤总体来说与全因子试验设计是一致的。但是有一个要注意的地方：在第一步选定模型中显著的主效应和2阶交互作用时，当某些2阶交互作用效用显著时，不能仅从表面上的结果来定取舍，要仔细分析混杂结构，查看在结构表中，此显著项是与哪个（或哪些）2阶交互作用效应相混杂的，再根据背景材料予以判断，最终决定入选。比如：数据显示B,C,D以及AD是显著的，但是背景材料又说明A和D没有交互作用，而AD与BC是相混杂的，这个时候，应该是B,C,D以及BC是显著的。
实验之三：响应曲面设计
响应曲面设计包括两种方法：中心复合设计和Box-Behnken设计。在中心复合设计中，整个试验由下面三部分试验点组成：
（1）立方体点[或称角点]，各坐标皆为1或-1，这是因子试验的组成部分；
（2）中心点，各点坐标皆为0；
（3）星号点[或轴点]，除了一个自变量的坐标为±ɑ外，其余自变量皆为0，在k个因子的情况下，共有2k个星号点。
中心复合设计，包括三种设计：
中心复合序贯设计，当“水平定义”栏中选定“立方点”时，表示这时设定的水平作为立方点，星号点将超出立方体。
中心复合有界设计，当“水平定义”栏中选定“轴点”时，表示这时设定的水平作为轴上的星号点，立方点将向内收缩。
中心复合表面设计，意味着将星号点的位置向中心收缩而设定在立方体的表面上。
Box-Behnken设计，这种设计是将各试验点取在立方体的棱的中点上，除非极端重视试验次数，否则通常不采用这种设计。
步骤1：响应曲面的计划
例：提高密封胶条黏合力试验。影响黏合力的3个因子是：A:烘烤温度（220-240摄氏度）、B:烘烤时间（6-10秒）、C:黏合压力（100-140帕）。在因子设计中，分别取下列条件，安排了全因子试验：
A:烘烤温度,低水平220，高水平240（摄氏度）
B:烘烤时间，低水平7，高水平9（秒）
C:黏合压力，低水平110，高水平130（帕）
试验后发现，数据明显呈现弯曲状况，希望进一步安排些实验以拟合响应曲面方程。由于要进行序贯试验，最好选中心复合设计。
具体做法是：选择[统计]>[DOE]>[响应曲面]>[创建响应曲面设计]，打开创建响应曲面设计对话框。在“设计类型”中选择“中心复合”，在“因子数”中设定为3。
打开“显示可用设计”对话框，可以看到未划分区组时试验次数为20。打开“设计”后，本例中需要的试验次数为20次，这是可行的，因此不必修改，中心点数也不用另设；但是选取哪种中心复合设计，需要考虑更多条件，由于在烘烤温度上，原来的试验温度条件已经取在边界上了，不允许再超界因而不能使用中心复合序贯设计，但是又考虑到要保持序贯性，只能放弃中心复合边界设计（没有序贯性），因而选用中心复合表面设计，即在Alpha值中选择表面中心；在“因子”中，选择“立方点”，并填写各因子的名称及水平；在“选项”中，为了看清楚结构，暂时先删除随机化。单节确定。
结果如下：
在表中的20次试验中，第1至第8号因子点以及第15至17号中心点，已经在因子设计阶段获得了数据，只要将这些结果填在后面第8列上，然后再补充其他9个试验，及可以完成全部响应曲面的试验任务。
步骤2：响应曲面设计的分析
例：提高烧碱纯度问题。在烧碱生产过程中，经过因子的筛选，最后得知反应炉内压力及温度是两个关键因子。在改进阶段进行全因子试验，因子A压力的低水平和高水平分别取为50帕和60帕，因子B反应温度的低水平和高水平分别取为260及320摄氏度，在中心点处也作了3次试验，试验结果在数据文件：DOE_烧碱纯度（响应1）。
对于这批数据按全因子试验进行分析，具体操作为：选择[统计]=>[DOE]=>[因子]=>[分析因子设计]，打开分析因子设计对话框。首先将全部备选项列入模型，删除在模型中包括中心点，在“图形”中的残差与变量下将压力和温度选入进去。得到的结果如下：
纯度的效应和系数的估计（已编码单位）
项效应系数系数标准误TP
常量96.9610.4150233.630.000
压力-2.665-1.3320.5490-2.430.094
温度-0.765-0.3820.5490-0.700.536
压力温度0.0350.0180.54900.030.977
S=1.09803PRESS=134.203
R-Sq=68.01%R-Sq（预测）=0.00%R-Sq（调整）=36.01%
对于纯度方差分析（已编码单位）
来源自由度SeqSSAdjSSAdjMSFP
主效应27.68747.687453.843723.190.181
2因子交互作用10.00120.001230.001230.000.977
残差误差33.61703.617011.20567
弯曲13.51783.517813.5178170.920.014
纯误差20.09920.099200.04960
合计611.3057
从上述表中可以看到，主效应和2因子交互作用对应的概率P值均大于0.1，说明模型的总效应不显著，而且弯曲对应的概率P值为0.014，拒绝原假设，认为存在明显的弯曲趋势；R-Sq和R-Sq（预测）的值都比较小，说明了模型的总效果不显著。
从残差与各变量的图也验证了存在严重的弯曲现象。这些都表明，对响应变量单纯地拟合一阶线性方程已经不够了，需要再补充些“星号点”，构成一个完整的响应曲面设计，拟合一个含二阶项的方程就可能问题了。补充的4个星号点的实验结果见数据表：DOE_烧碱纯度（响应2）。
下面对全部11个点构成的中心复合序贯设计进行分析，拟合一个完整的响应曲面模型。分析如下：
第一步：拟合选定模型。
选择[统计]>[DOE]>[响应曲面]>[分析响应曲面设计]，打开分析响应曲面设计对话框。点击窗口“项”以后，可以看到模型中将全部备选项都列入了模型，包括A(压力)、B（温度）以及它们的平方项AA、BB和交互作用项AB；打开“图形”窗口，选定“正规”、“四合一”以及残差与变量，并将压力和温度都选入残差与变量中；打开“储存”窗口，选定“拟合值”、“残差”以及“设计矩阵”。单击确定。
得到的结果如下：
纯度的估计回归系数
项系数系数标准误TP
常量97.78040.10502931.0660.000
压力-1.89110.09114-20.7500.000
温度-0.60530.09092-6.6570.001
压力压力-2.58220.15339-16.8350.000
温度温度-0.46150.15314-3.0140.030
压力温度0.03510.182530.1920.855
S=0.181900PRESS=0.693667
R-Sq=99.35%R-Sq（预测）=97.27%R-Sq（调整）=98.70%
对于纯度的方差分析
来源自由度SeqSSAdjSSAdjMSFP
回归525.231025.23105.04620152.510.000
线性215.712715.71277.85635237.440.000
平方29.51719.51714.75853143.820.000
交互作用10.00120.00120.001230.040.855
残差误差50.16540.16540.03309
失拟30.06620.06620.022080.450.747
纯误差20.09920.09920.04960
合计1025.3964
结果解释：
（1）看方差分析表中的总效果。在本例中，回归项的P值为0.000，表明应该拒绝原假设，认为本模型总的来说是有效的。
看方差分析表中的失拟现象，本例中，失拟项对应的P值为0.747，明显大于显著性水平0.05，接受原假设，认为本模型中不存在失拟现象。
（2）看拟合的总效果。本例中，R-Sq与R-Sq（调整）比较接近，认为模型的拟合效果比较好；R-Sq（预测）比较接近于R-Sq值且这个值比较大，说明将来用这个模型进行预测的效果比较可信。
（3）各效应的显著性。从表中可以看到，压力、温度以及它们的平方项对应的概率值都小于显著性水平，说明这些效应都是显著的；而压力和温度的交互效应项对应的概率值为0.855，显然大于显著性水平，认为该效应项是不显著的。
第二步：进行残差诊断
利用自动输出的残差图来进行残差诊断。
从上述残差图中可以看出，残差的状况是正常的。
第三步：判断模型是否需要改进。
根据第一步的分析，我们得知压力和温度的交互作用项是不显著的，应该予以剔除，因此需要重新拟合新的模型，使得新的模型中不包含交互作用项。
具体实现步骤是：在项中将交互作用项剔除，在结果中输出标准化残差和删后残差。
得到的结果为：
纯度的估计回归系数
项系数系数标准误TP
常量97.78040.096221016.1770.000
压力-1.89110.08350-22.6470.000
温度-0.60530.08331-7.2650.000
压力压力-2.58220.14054-18.3730.000
温度温度-0.46150.14031-3.2890.017
S=0.166665PRESS=0.546550
R-Sq=99.34%R-Sq（预测）=97.85%R-Sq（调整）=98.91%
对于纯度的方差分析
来源自由度SeqSSAdjSSAdjMSFP
回归425.229825.22986.30744227.070.000
线性215.712715.71277.85635282.830.000
平方29.51719.51714.75853171.310.000
残差误差60.16670.16670.02778
失拟40.06750.06750.016870.340.836
纯误差20.09920.09920.04960
合计1025.3964
纯度的估计回归系数，使用未编码单位的数据
项系数
常量-59.9731
压力5.36834
温度0.134611
压力压力-0.0512244
温度温度-2.56700E-04
结果解释：
（1）先看方差分析表中的总效果。回归项对应的P值为0.000，拒绝原假设，说明回归模型总的来说是有效的；看方差分析表中的失拟现象，可以看到失拟对应的P值为0.836，大于0.05，接受原假设，即可以判定，本模型删去了一项，但没有造成失拟现象。
（2）看删减后的模型是否比原来的有所改进。
全模型 变化 删减模型 R-Sq 99.35% 减小 99.34% R-Sq（调整） 98.70% 增大 98.91% S 0.181900 减小 0.166665 R-Sq（预测） 97.27% 增大 97.85% PRESS 0.693677 减小 0.546550 由于模型项缺少了一项，R-Sq通常会有所降低，但关键要看调整的R-Sq（调整）是否有所提高，s值是否有所降低，预测残差平方和PRESS是否有所降低，R-Sq（预测）是否有所提高。从表中来看，均符合上述要求，表明删除了不显著的交互作用后，回归的效果更好了。
此外，我们还可以得到最后确定的回归方程：
从标准化残差以及删后残差的结果分析表中，可以看到这些值都小于2，因此认为新的模型的残差没有发现任何不正常的情况。
第四步：对选定的模型进行分析解释。
通过前面得到的回归方程，运用数学方法我们可以得到使得纯度最大的A和B分别取什么值，但是不能保证该最大值就一定落在试验范围之内。在求解前，先看一下等值线图和曲面图，具体实现：[统计]>[DOE]>[响应曲面]>[等值线图/曲面图]。从图中可以看到，在原试验范围内确实有个最大值。
运用人工解方程的方法，可以得到当压力=52.4、温度=262.2时所获得的纯度最高。
与此同时，计算机也提供了响应优化器，可以直接获得最优解。具体操作是：[统计]>[DOE]>[响应曲面]>[响应优化器]，打开响应优化器对话框。将“纯度“选入到所选项中，打开“设置”对话框，在“目标”中选择望大，在下界和望目中分别输入80和100。单击确定。
得到的结果为：
从上图中也可以看到，在压力=52.3465、温度=262.1616时，纯度达到最大值为98.3249，与我们手算的结果是一样的。
为了获得置信区间，从“统计]>[DOE]>[响应曲面]>[分析响应曲面设计]”入口，选定“响应”为纯度，在“预测”中，在自变量设置处，填写“52.4,262.2”则可以得到如下结果：
使用纯度模型的新设计点数的预测响应
点拟合值拟合值标准误95%置信区间95%预测区间
198.32500.0859139(98.1148,98.5353)(97.8662,98.7839)
从结果中可以看到，预测结果的值与我们最优化的值是一样的，说明预测结果是可信的。前一个置信区间表明的是回归方程上的点的置信区间，此值可以作为改进的结果的预报写在总结报告中；后一个置信区间表明的是以上述回归方程上的预测值的置信区间为基础，加上观测值固有的波动所给出的置信区间，这就是将来做一次验证试验时将要落入的范围，可供做验证试验时使用。
Minitab_专业论坛www.mimitab.cn
|
|
