2水平3因子设计

本文是实验设计与分析（第6版，Montgomery著傅珏生译)第6章2k析因设计第6.3节的python解决方案。本文尽量避免重复书中的理论，着于提供python解决方案，并与原书的运算结果进行对比。您可以从Detail 下载实验设计与分析（第6版，Montgomery著傅珏生译)电子版。本文假定您已具备python基础，如果您还没有python的基础，可以从Detail 下载相关资料进行学习。

设有3个因子，A，B，C，每一因子有两个水平。此设计叫做23析因设计。现在，如图6.4a所示，在一个立方体上显示出8个处理组合。用“-”和“+”的记号分别代表因子的低水平和高水平，我们在图6.5b中列出23设计的8个实验。有时我们称为设计矩阵。推广6.2节所用的记号，将处理组合依标准顺序写为（1），a，b，ab，c，ac，bc，abc。这些符号亦代表在所指定的处理组

合上的n个观测值的总和。

实际上有3种不同的记号广泛应用于2k设计的实验中。先是“+”与“-”记号，常称之为几何记号（geometric notation）。其次是用小写字母表示处理组合。最后，用记号1和0分别表示因子的高水平和低水平，以代替“+”和“-”，对于23设计，这些不同的记号说明如下：

23设计的8个处理组合间有7个自由度。3个自由度和A，B，C的主效应有关。4个自由度与交互作用有关：AB，AC，BC各有1个，ABC有1个。

考虑估计主效应。首先考虑估计主效应A。当B与C处于低水平时，A的效应是[a-（1）]/n。同理，当B处于高水平、C处于低水平时，A的效应是[ab-c]/n。当C处于高水平，B处于低水平时，A的效应是[ac-bc]/n。最后，当B和C都处于高水平时，A的效应是[abc-bc]/n。这样一来，A的平均效应正是这4个效应的平均值，即

此式也可以用图6.5a的立方体右边一面的4个处理组合（其中A处于高水平）和左边一面的4个处理组合（其中A处于低水平）之间的对照推导出来。也就是说，A效应恰好是A处于高水平时4个试验的平均值（`yA+）减去A处于低水平时4个试验的平均值（`yA-），即

此式可以重新排为

它与（6.11）式相同。

同理，B的效应是立方体前边一面的4个处理组合的平均值和后边一面的4个处理组合的平均值之差，即

C的效应是立方体上边一面的4个处理组合的平均值和下边一面的4个处理组合的平均值之差，即

二因子交互作用效应容易计算出来。AB交互作用的一种度鼠是A在B的两个水平上的平均效应之差。为方便起见，称此差的一半为AB交互作用。用符号表示为，

因为AB交互作用是此差的一半，所以

（6.14）式可写为：

从这种形式容易看出，AB交互作用就是图6.5b的立方体中在两个对角面上的试验的平均值之差。同理并参考图6.5b，得AC和BC交互作用分别是

ABC交互作用定义为AB在C的两个不同水平上的交互作用之差的平均值.于是

和前面一样，可以把ABC交互作用看作为两个平均值之差。如果把两个平均值中的试验分离为两组，那么，它们就是图6.5c的立方体中的两个四面体的顶点。

在（6.11）式至（6.17）式中，方括号中的量是处理组合的对照。一张加减符号表可以从这些对照中导出来，如表6.3所示。主效应的符号是：高水平取加号，低水平取减号。一旦主效应的符号确定，其余各列的符号可以逐行地用前面恰当的列的符号相乘而得。例如，AB列的符号就是A列的符号与B列的符号逐行的乘积。任一效应的对照容易由此表求得

表6.3有几个有趣的性质：（1）除列I之外，每列加号的个数与减号的个数相等。（2）任意两列符号乘积的和为零。（3）列I与任一列相乘，该列的符号不变。也就是说，I是一个单位元素（identity element），（4）任意两列相乘，得出表中的一列。例如，A×B=AB，以及AB×B=AB2=A

我们看到，乘积中的指数可以利用模为2的算术运算法求得。（也就是说，指数仅是0或1：如果大于1，就将它减少2的某个倍数，直至它为0或1为止）所有这些性质都包含在用来估计效应的对照的正交性中。

因为每一效应都有一个相应的单自由度对照，所以效应的平方和是容易计算的。在有n次重复的23设计中，任一效应的平方和是

例6.1用23析因设计在单晶片等离子蚀刻工具中开发氮化物蚀刻工序。设计因子包括电极间的间隙、气流（C2F6作为反应物气体）和应用在阴极的RF功率（见图3.1的晶片蚀刻工具示意图）。每个因子都有两个水平，设计重复两次。响应变量是硅氮化物的蚀刻率（A/m）。表6.4中显示的是蚀刻岸数据。图6.6几何地表示了该设计

用表6.4的处理组合的总和，估计因子效应如下：

最大的效应是功率（C=306.125）、间隙（A=-101.625），功率-间交互作用（AC=-153，625）。

由（6.18）式计算得平方和如下：

总平方和是SST=531420.9375，用减法得SSE=18020.50。表6.5中汇总了效应估计及平方和。标为“百分比贡献率”的列度量了模型的每一项对总的平方和的贡献百分比。贡献百分比通常是粗略的，但它是模型每一项相对重要性的有效指示。注意，主效应C（功率）事实上主导了该工序，因为它占了总体变异性的70％多，面主效应A（间隙）和AC的交互作用各占大约8％和18％

表6.6中的方差分析可用于确定这些效应的大小。它表明间隙和功率的主效应是高度显著的（两个P值都很小）。AC的交互作用也是高度显著的，因此，间家和功率之间有强烈的交互作用.

import numpy as np

import matplotlib.pyplot as plt

import statsmodels.formula.api as smf

import pandas as pd

import statsmodels.api as sm

from sklearn import linear_model

import seaborn as sns

import matplotlib.pyplot as plt

import mistat

from scipy import special

from scipy.stats import ncx2

x1 = [-1,1, -1,1, -1,1, -1,1,  -1,1, -1,1, -1,1, -1,1 ]

x2 = [-1,-1,1,1, -1, -1,1,1, -1,-1,1,1, -1, -1,1,1]

x3 = [-1, -1,-1,-1,1,1,1,1, -1, -1,-1,-1,1,1,1,1 ]

y1 = [550,669,633,642,1037,749,1075,729,604,650,601,635,1052,868,1063,860]

data = {"y1":y1,"x1":x1,"x2":x2,"x3":x3}

data=pd.DataFrame(data)

model1 = smf.ols(formula='data.y1 ~ C(data.x1)*C(data.x2)*C(data.x3)', data=data).fit()

print(model1.summary2())                 

anova_table = sm.stats.anova_lm(model1, typ=2)

anova_table

>>> anova_table

                                       sum_sq      df           F    PR(>F)

C(data.x1)                         41310.5625        1.0   18.339364   0.002679

C(data.x2)                           217.5625        1.0    0.096584  0.763911

C(data.x3)                        374850.0625        1.0    166.410505  0.000001

C(data.x1):C(data.x2)               2475.0625          1.0      1.098776  0.325168

C(data.x1):C(data.x3)              94402.5625          1.0      41.908965  0.000193

C(data.x2):C(data.x3)                 18.0625          1.0      0.008019  0.930849

C(data.x1):C(data.x2):C(data.x3)     126.5625            1.0        0.056186  0.818586

Residual                           18020.5000       8.0         NaN       NaN