题目
- 来自《试验设计与建模》方开泰 刘民千 周永道 主编 高等教育出版社 2 课后习题3.4
试验结果的直观分析
正交试验表如下：
| No. | A | B | C | 4 | y |
|---|---|---|---|---|---|
| 1 | 1 | 1 | 1 | 1 | 32 |
| 2 | 1 | 2 | 2 | 2 | 55 |
| 3 | 1 | 3 | 3 | 3 | 39 |
| 4 | 2 | 1 | 2 | 3 | 53 |
| 5 | 2 | 2 | 3 | 1 | 49 |
| 6 | 2 | 3 | 1 | 2 | 42 |
| 7 | 3 | 1 | 3 | 2 | 56 |
| 8 | 3 | 2 | 1 | 3 | 61 |
| 9 | 3 | 3 | 2 | 1 | 63 |
| \(T_1\) | 126 | 141 | 135 | ||
| \(T_2\) | 144 | 165 | 171 | ||
| \(T_3\) | 180 | 144 | 144 | ||
| \(m_1\) | 42 | 47 | 45 | ||
| \(m_2\) | 48 | 55 | 57 | ||
| \(m_3\) | 60 | 48 | 48 | ||
| \(R\) | 18 | 8 | 12 |
- 通过直接比较转化率，可以看到最佳的水平组合来自第九组试验:\((A_3,B_3,C_2)=(90^oC,150min,6\%)\).
- 计算各因子在每个水平下的转化率之和\(T_i\)和平均转化率\(m_i\)，最佳的组合均为\((A_3,B_2,C_2)\).
- 通过计算三个因子三个水平下的平均转化率的极差\(R\),得到因子对响应的影响程度如下：\(A>C>B\).
试验结果的方差分析
给定统计模型如下：
\[\left\{\begin{array}{l}
y_{i j k}=\mu+\alpha_{i}+\beta_{j}+\gamma_{k}+\varepsilon_{i j k}, i, j, k=1,2,3, \\
\alpha_{1}+\alpha_{2}+\alpha_{3}=\beta_{1}+\beta_{2}+\beta_{3}=\gamma_{1}+\gamma_{2}+\gamma_{3}=0, \\
\varepsilon_{i j k} \sim N\left(0, \sigma^{2}\right),iid .
\end{array}\right.\]
通过R得到的方差分析表如下：
- 从方差分析得到的结果可以看出，因素A是最显著的，其次是因素C，因素B是最不显著的。
得到对各个因素各个水平主效应的估计如下：
试验结果的回归分析
在正式建模之前，画出三个因素与响应的散点图如下，先看看它们之间的关系是否为线性：
- 通过最后一行可以看出，因素\(A\)与响应\(y\)大致是线性的，因素\(B,C\)可能存在非线性关系。
然后建立\(A,B,C\)三个因素关于响应\(y\)的线性模型，得到的结果如下：
- 可以看到，只有因素\(A\)是比较显著的，模型拟合的并不优秀
考虑加入二次项\(A^2,B^2,C^2\)和交叉项\(AB,AC,BC\),继续拟合。并使用逐步回归进行变量筛选得到的最优回归结果如下，不管是向前回归还是向后回归还是逐步回归最后筛选出的结果都相同，均为\(y\sim B^2+C^2+AC+BC\)
- 此时，因素\(AC\)是显著的，模型的决定系数\(R^2\)也提高到0.7445，这与之前两种分析得到的结论相同。
结论
通过以上三种分析，可以知道转化率最高的水平组合为\(A_3B_2C_2\)或\(A_3B_3C_2\),因素\(A\)反应温度对转化率的影响是最显著的，其次是因素\(C\)用碱量，因素\(B\)反应时间对转化率的影响很小。
代码
y=c(32,55,39,53,49,42,56,61,63)
tag=matrix(c(1,1,1,2,2,2,3,3,3,
1,2,3,1,2,3,1,2,3,
1,2,3,2,3,1,3,1,2,
1,2,3,3,1,2,2,3,1,
y),ncol = 5)
colnames(tag)=c("A","B","C","","y")
dat=as.data.frame(tag)
sum(dat$y[dat$A==1])
sum(dat$y[dat$A==2])
sum(dat$y[dat$A==3])
sum(dat$y[dat$B==1])
sum(dat$y[dat$B==2])
sum(dat$y[dat$B==3])
sum(dat$y[dat$C==1])
sum(dat$y[dat$C==2])
sum(dat$y[dat$C==3])
sum(dat$y[dat$A==1])/3
sum(dat$y[dat$A==2])/3
sum(dat$y[dat$A==3])/3
sum(dat$y[dat$B==1])/3
sum(dat$y[dat$B==2])/3
sum(dat$y[dat$B==3])/3
sum(dat$y[dat$C==1])/3
sum(dat$y[dat$C==2])/3
sum(dat$y[dat$C==3])/3
m=c(42,48,60,47,55,48,45,57,48)
x=c(80,85,90,90,120,150,5,6,7)
plot(m)
plot(x[1:3],m[1:3],ylim = c(40,65),type = "b",main = "A",xlab = "A",ylab = "Y")
plot(x[4:6],m[4:6],ylim = c(40,65),type = "b",main = "B",xlab="B",ylab = "Y")
plot(x[7:9],m[7:9],ylim = c(40,65),type = "b",main = "C",xlab="C",ylab = "Y")
#方差分析
aov1=aov(lm(y~factor(A)+factor(B)+factor(C),data = dat))
summary(aov1)
model.tables(aov1)
#回归分析
M=matrix(c(80,80,80,85,85,85,90,90,90,
90,120,150,90,120,150,90,120,150,
5,6,7,6,7,5,7,5,6,
1,2,3,3,1,2,2,3,1,
y),ncol = 5)
colnames(M)=c("A","B","C","","y")
dat0=data.frame(M)
pairs(dat0)
lm0=lm(y~A+B+C,data = dat0)
summary(lm0)
dat0$A2=dat0$A^2
dat0$B2=dat0$B^2
dat0$C2=dat0$C^2
dat0$BC=dat0$B*dat0$C
dat0$AB=dat0$A*dat0$B
dat0$AC=dat0$A*dat0$C
lm1=lm(y~A+AC++C2,data = dat0)
summary(lm1)
lm2=lm(y~A+B2+C2+BC,data = dat0)
summary(lm2)
lm=lm(y~A+B2+C2+AB+AC+BC,data = dat0)
summary(lm)
step(lm,method="forward")
step(lm,method="both")
step(lm,method="backward")
lm3=lm(y~B2+C2+AC+BC,data = dat0)
summary(lm3)