1、试验设计和建模周永道四川大学数学学院试验设计和分析教材教材: : 方开泰、刘民千、周永道方开泰、刘民千、周永道(2011(2011) ),试验设计,试验设计 和建模。和建模。 期末最终成绩构成期末最终成绩构成: 期末考试期末考试: 70% 作业作业: 10% 随堂测试随堂测试: 10% 小课题小课题: 10%软件: MATLAB, SPSS 参考书 方开泰、马长兴,正交与均勻试验设计,科学出版社,2001. Douglas C. Montgomery. Design and Analysis of Experiments, 6th Edition,中国邮电出版社, 2007. Hamada,
2、M. and Wu, Jeff C.F., Experiments: Planning, Analysis, and Parameter Design Optimization , Wiley, 2000. Cornell, J.A. Experiments with Mixtures, 3nd Ed., Wiley, 2002. Fang, K.T., Li, R. and Sudjianto, A. Design and Modeling for Computer Experiments, Chapman & Hall/CRC Press, London, 2005.第一章1.1 1.1
3、科学试验科学试验1.1.1 试验的重要性试验的重要性 科学试验是人们认识自然、了解自然的重要手段。 许多重要的科学规律都通过科学试验发现和证实。 随着科学和技术的发展,试验涉及的因素越来越多,它们之间的关系更加复杂,光凭经验已不能达到预期要求,于是产生了试验设计这门学科。 设计一个试验涉及到试验目的、试验方案、技术保证、分析数据以及有关组织管理等。这些环节有的是属于管理科学,有的是需要数学和统计学的方法来设计试验方案,后者称为统计试验设计统计试验设计, 它是统计学的一个重要分支。 是统计学的重要分支,它能大量节省试验是统计学的重要分支,它能大量节省试验的次数。能将试验数据从随机误差的烟幕的次数
4、。能将试验数据从随机误差的烟幕中去伪存真,抓住事物的规律。中去伪存真,抓住事物的规律。 所以所以一个精心设计的试验是认识世界的一个精心设计的试验是认识世界的有效方法有效方法 (Atkinson and Donev (1992)。 奠定了现代遗传理论的基础例例1.1.1 1孟德尔豌豆实验孟德尔豌豆实验例例1.1.2 2化工试验化工试验在某化工产品的合成工艺中,考虑反应温度(A)、压力(B) 和催化剂用量(C),并选择了试验范围分别为:温度(A): 80oC120oC;压力(B): 46 大气压;催化剂用量(C): 0.5%1.5%; 我们需要选择这三个因素的最佳组合,以达到高产的目的。许多产品都
5、是混合多种成分在一起形成的。许多产品都是混合多种成分在一起形成的。面粉水糖蔬菜汁 椰子汁盐发酵粉乳酸钙 咖啡粉香料色素面包面包怎样确定各种成分的比例呢?怎样确定各种成分的比例呢?经验经验试验试验混料试验混料试验例例1.1.3 3 加工面包试验加工面包试验例例1.41.4环保试验环保试验 在水及食物中的某些化学元素,吃多了对人体是在水及食物中的某些化学元素,吃多了对人体是有害的,为了研究这些元素对人体健康的影响。有害的,为了研究这些元素对人体健康的影响。Cadmium (Cd)镉 Copper (Cu) 铜Zinc (ZN) 锌Nickel (Ni)镍Chromium(Cr) 铬Lead (Pb
6、)铅因素因素0 200 200 200 200 20 0 20 0.01, 0.05, 0.1, 0.2, 0.4, 0.8, 1, 2, 4, 5, 8, 10, 12, 14, 16, 18, 200.01, 0.05, 0.1, 0.2, 0.4, 0.8, 1, 2, 4, 5, 8, 10, 12, 14, 16, 18, 200.01, 0.05, 0.1, 0.2, 0.4, 0.8, 1, 2, 4, 5, 8, 10, 12, 14, 16, 18, 200.01, 0.05, 0.1, 0.2, 0.4, 0.8, 1, 2, 4, 5, 8, 10, 12, 14, 16
7、, 18, 200.01, 0.05, 0.1, 0.2, 0.4, 0.8, 1, 2, 4, 5, 8, 10, 12, 14, 16, 18, 200.01, 0.05, 0.1, 0.2, 0.4, 0.8, 1, 2, 4, 5, 8, 10, 12, 14, 16, 18, 20范围 和 水平 增加产量增加产量 提高质量提高质量 降低成本降低成本 缩短研究时间缩短研究时间科学试验是人类了解自然的手段,通过试验来了解科学试验是人类了解自然的手段,通过试验来了解因素和指标(响应)之间的关系,希望因素和指标(响应)之间的关系,希望一个好的试验设计是用最少的试验次数获得最多的一个好的试验设
8、计是用最少的试验次数获得最多的有用信息。有用信息。 试验设计试验设计的目的的目的v水平组合的比较v建模v参数估计v证实猜想v优化v筛选v发现规律v等等Experiments are performed by investigators in virtually all fields of inquiry, usually to discover something about a particular process or system.Scientific experiments are of essential importance in peoplessurviving and expl
9、oring of nature.A well designed experiment is an efficient method of learning about the worldThe Usefulness of Experimental Design多因素多因素统计模型未知统计模型未知响应曲面多峰响应曲面多峰当代科学当代科学试验的复杂性试验的复杂性非线性非线性响应曲面无解析表达响应曲面无解析表达多峰非线性模型瞎子摸象1.1.2 试验的重要元素试验的重要元素让我们首先通过一个例子来介绍一些重要元素 例例 在一个化工试验中, 试验者希望通过如下 的可控变量来增加产量: x1:原料品种原料
10、品种 m1, m2, m3x2:加酸量加酸量 (ml) 10,28x3:反应时间反应时间 (时) 0.5, 3.5因素(因子) 在试验中可控的并用于考察对试验结果(y)的变量称为 因素因素 或 因子因子 (factor)。如反应温度、压力、催化剂品种、施化肥量、水稻品种等。如反应温度、压力、催化剂品种、施化肥量、水稻品种等。 因素可以是 定量定量 的, 也可以是 定性定性 的。 水平 (level) 因素变化的范围称为试验区域试验区域, 在本例中,试验区域为: m1, m2, m3 x 10,28 x 0.5,3.5.原料品种: m1, m2, m3加酸量:10,19,28反应时间: 0.5,
11、 1.5, 2.5, 3.5水平组合 因素诸水平的组合称为 水平组合水平组合 (level-combination), 如 m2, 10, 2.5, m1, 28, 0.5 。水平组合在文献中又称为 处理组合处理组合 。 一一个个因因子子设计设计 (Factorial design) 是一组水平是一组水平组合组合。 处理, 响应 在试验环境下对确定的水平组合所作的试验称为一个 处理处理 (trial 或或 run) 。 试验的结果称为 响响应应 (response), 响应可以是定性的, 也可是定量的。 不可控的诸微小因素之总和, 称为 随机误差随机误差。 同样条件下的两次试验结果可能不同。随
12、机误差存在于一切试验之中。随机误差存在于一切试验之中。 随机误差 (random error)随机误差随机误差随机误差可假定遵从 正态分布 。 方差 给出随机误差大之度量 。),N(202令 为重复试验之响应值 , y, yn11 iiy , i, , n,me=+=这里, 为真值, m独立同分布,遵从 。 1i, i, , ne=),N(201221/, 1() .1niiniiyynyynms=- 和 的无偏估计为 m2AA1A2yLevely29.532.0Response31.033.030.534.5Total190.5均值 = 190.5/6 = 31.75自由度 : 562211
13、()16.3753.2755iymeandfs=-=y - mean-2.250.25-0.751.25-1.252.750(y - mean)25.06250.06250.56251.56251.56257.562516.375随机误差随机误差: : 部份因子设计部份因子设计设有 s 个因素, 它们分别取个水平。则全部水平组合有 s, q, q1一个一个 水平组合水平组合 可视为可视为 s 维空间的一个点维空间的一个点, 称为称为试验点试验点 。1N siq=个。 例如, 一个六因素, 五水平的全面试验至少需要 次试验。6515625 N = 全面试验若所有的水平组合都作相同重复数的试验,
14、称为全面试验全面试验 。 在农业、生物等试验中,很难做到试验条件完全一样。 区组的概念成为古典试验设计中非常有用的工具,同一区组的试验有十分近似的试验环境。 区组设计可以避免或减少系统误差的干扰,从而大大提高试验结论的可靠性。 在体育比赛中,区组及有关设计已在普遍使用。 区组 试验的环境随着时间的推移,可能有趋势型的变化,如室温渐高、湿度渐小、电压波动加剧等。 为了使试验的结论更加可靠,随机化随机化是用来减少试验误差的重要手段。常用的是对试验次序随机化,哪个试验先做,哪个试验后做,随机决定。 若试验有区组,要根据试验的具体情形采取所有试验的完全随机化,或仅区组内的试验随机化。 随机化 同一个试
15、验重复两次或多次是减少试验误差干扰的一种方法,在传统的计算方法中经常使用。 若 y1, ym 是同一个试验条件下的响应,且 yi 独立同分布,方差为 ,则均值 均值的方差 重复11miiyym=2s( )2/Var yms=传统试验的三个基本原则三个基本原则: 重复性、随机化、分区组重复性、随机化、分区组 针对不同的试验,试验者要选择合适的试验方法,建立相应的统计模型 统计模型试验的组织和管理 一支专业队伍 明确的试验目标 科学的试验方案 试验中,处理可控与不可控因素A. 试验实现方式:试验实现方式:1.1.3 试验的类型试验的类型 传统的试验传统的试验 实验室试验实验室试验 工业试验工业试验
16、 计算机试验计算机试验 计算机模拟计算计算机模拟计算 寻找近似模型寻找近似模型B. 因素约束条件因素约束条件 无约束试验无约束试验 诸因素可以自由的选择试验的值,不受其它因素约束,试验区域是一个超矩形 混料试验混料试验 因素之间的取值会相互影响,例如或 单因素试验单因素试验 水平数可以适当多取,而且可以考虑做重复试验 多因素试验多因素试验 各因素的水平数一般不能取得很大 二水平试验 多水平试验C. 因素个数因素个数D. 响应个数响应个数 单响应试验单响应试验 每次试验只观察一个响应值。如产量 多响应试验多响应试验 每次试验观察多个响应值。如鞋子橡胶底的试验响应:强度、弹性和最大弯曲次数等等 多
17、媒体试验多媒体试验 试验有无穷多个的响应。例如,响应是人的指纹、化学或生物中指纹曲线、声音的曲线、图像的颜色及深浅,等等E. 试验轮次试验轮次 单一试验单一试验 一次试验达到要求一次试验达到要求 序贯试验序贯试验 优选法优选法 响应曲面分析响应曲面分析 均匀序贯试验均匀序贯试验 单区组试验单区组试验 每次试验在相同或十分近似的条件下进行每次试验在相同或十分近似的条件下进行 区组试验区组试验 目的是使得组内的差异比组间差异小目的是使得组内的差异比组间差异小 常见的区组有以日、月、年、批次、双胞胎,常见的区组有以日、月、年、批次、双胞胎,等等等等F. 试验分组试验分组例例1.5. (自由落体运动自
18、由落体运动) 若不计空气阻力,自由落体运动的初始速度为零,记下落时间为x (秒)(s),下落距离y (米)(m),人们发现它们之间有如下规律 g 为重力加速度。设想试验者对关系(1.3) 一无所知,希望通过试验来揭示y 和x 之间的关系 (1.3)试验结果可用二次回归模型拟合 方差分析模型方差分析模型 因子设计,正交设计因子设计,正交设计 参数回归模型参数回归模型 最优设计最优设计 非参数回归模型非参数回归模型 均匀设计均匀设计 稳健回归模型稳健回归模型例例1.6:威布尔生长曲线:威布尔生长曲线方差分析模型方差分析模型1 +, 1 , 1 0ijjijjijjqyj= , , q, i= ,
19、, n ,memaeaa=+=在在0,10中取若干个点作试验,设中取若干个点作试验,设 x1, , xq 为为试验点,试验点,n1, , nq 为其重复数,其统计模型为为其重复数,其统计模型为 二水平试验在西方被广泛推荐二水平试验在西方被广泛推荐 二水平不足以揭示非线性关系二水平不足以揭示非线性关系 多水平试验值得推荐多水平试验值得推荐二水平试验的不足二水平试验的不足试验范围对,试验范围对,但水平不合适但水平不合适试验范围及水平试验范围及水平都对,但不能揭都对,但不能揭示示A和和Y之间更复之间更复杂的关系杂的关系试验范试验范围错围错只能预报四个水平处的响只能预报四个水平处的响应值,进一步采用回
20、归模应值,进一步采用回归模型是有益的。型是有益的。 方差分析模型方差分析模型 因子设计,正交设计因子设计,正交设计 参数回归模型参数回归模型 最优设计最优设计 非参数回归模型非参数回归模型 均匀设计均匀设计 稳健回归模型稳健回归模型根据专业知识根据专业知识, 可选用适当的回归模型可选用适当的回归模型, 比如比如用二次模型用二次模型 ,2210 xxy(x) ,332210 xxxy(x) ,)()(11xfxfy(x)mm回归模型回归模型其中函数其中函数 f1, fm已知已知, 但参数但参数 b1, bm未知。未知。或三次模型或三次模型更一般地更一般地, 给定试验次数给定试验次数 n,希望能获
21、得最精确的,希望能获得最精确的回归系数回归系数 b b0, b, b1, 的估计。的估计。 缺点:对模型的变化缺乏稳健性。缺点:对模型的变化缺乏稳健性。 Kiefer, J.C. (1958), Ann Math.Stat. Kiefer, J.C. (1959), JRSS, B, with discussion Atkinson, A.C. and Donev, A.N. (1992), Optimal Experimental Designs, Clavendon Press, Oxford三次回归模型的三次回归模型的D-最优设计及其拟合最优设计及其拟合如果采用如果采用4 4次多项式模型
22、次多项式模型, ,效果会显著地改进。效果会显著地改进。 方差分析模型方差分析模型 因子设计,正交设计因子设计,正交设计 参数回归模型参数回归模型 最优设计最优设计 非参数回归模型非参数回归模型 均匀设计均匀设计 稳健回归模型稳健回归模型, g(x)y式中函数形式式中函数形式 g(x)未知。希望通过试验求未知。希望通过试验求得得g(x)一个近似模型。这时,一个自然的一个近似模型。这时,一个自然的想法是将试验点在想法是将试验点在0,2上均匀散布,即均上均匀散布,即均匀设计。匀设计。 若试验者对模型未知,这时将面对非参若试验者对模型未知,这时将面对非参数回归模型数回归模型非参数回归模型非参数回归模型
23、均匀设计及其拟合多项式回归均匀设计及其拟合多项式回归均匀设计均匀设计是一种是一种试验设计方试验设计方法法。它可以用较少的试验次数,。它可以用较少的试验次数,安排多因素、多水平的析因试验安排多因素、多水平的析因试验,当试验者对析因试验的统计模,当试验者对析因试验的统计模型未知时,型未知时,均匀设计均匀设计是最好的设是最好的设计方法。计方法。均匀设计均匀设计也是也是仿真试验仿真试验设计设计和和稳健设计稳健设计的重要方法的重要方法。 方差分析模型方差分析模型 因子设计,正交设计因子设计,正交设计 参数回归模型参数回归模型 最优设计最优设计 非参数回归模型非参数回归模型 均匀设计均匀设计 稳健回归模型
24、稳健回归模型 稳健回归模型常用于部分模型已知的情形2012 y( )xxh xbbb例如: )()(xhxfy此时,可用一些 稳健设计稳健设计 或 均匀设计均匀设计.其中 f(x) 为已知函数,h(x) 为 偏离真实函数 的部分。)(xfy( )yg x稳健回归模型稳健回归模型即 f(x) 为参数 的线性函数。b回归模型回归模型:yi, xi1, xi2, , xi,p-1, i = 1, , ny = b0 + b1x1 + b2x2 + + b p-1 x p-1 + E() = 0,Var() = s2 未知或 yi = b0 + b1xi1 + b2xi2 + + b p-1 xi,
25、p-1 + i i,n i.i.d. E(i) = 0, Var(i) = s2.或更一般的, yi= b1 g1(xi) + b2 g2(xi) + + b p gp(xi) + i i,n i.i.d. E(i) = 0, Var(i) = s2.xi=( xi1, xi2, , xi,p-1), i=1, , n 一般回归模型的矩阵表示一般回归模型的矩阵表示:y = Gb b + E( ) = 0,Cov( ) = s2In其中 y : n1, G : np, b b : p1, : n1 其元素 i.i.d.(1.14).,)()()()()()()()()(,121212222111
26、21121npnpnnppnxgxgxgxgxgxgxgxgxgyyybbb b bG Gy线性模型 (1.14) 包括很多有用的模型: 线性模型 通过原点的线性模型 二次模型 中心化二次模型y = b0 + b1x1 + b2x2 + + b p-1 x p-1 + y = b1x1 + b2x2 + + b p-1 x p-1 + bbbkjijiijkiiixxxy110bbbkjijjiiijkiiiixxxxxxy110)()(1.16)对于模型 (1.14) (a) 估计估计 模型 (1.14) 的最小二乘估计为 性质: yGGG)(1b12)()(Cov ,)(GGEsbbb)/
27、()()(22pnGyGysbbs其中 M=GG 为 信息矩阵信息矩阵, 或有时称 M=GG/n 为信息矩阵. s s 2 2 的估计的矩阵表达形式的估计的矩阵表达形式() () (),1Q()yX yXy(IH)yHX(XX) X其中帽子矩阵yHIyQ)( E)(Eb bE(y) = Xb bCov(y) = s2In)(tr)()(2nIHIXHIXsbb0)()tr()tr(tr21222pnnnssssXX(XXXX)X(XHI1yHIyQ)( )(112npnpnb bs为无偏估计. 对于线性模型 (1.16), 在实际中常检验下面的假设: (k=p-1)A. 检验模型是否有意义.
28、H0: b1 = = bk = 0 VS H1: 某 bj 0, 1 j kB. 检验某个变量 xj 是否对模型有显著贡献. H0: bj = 0 VS H1: bj 0(b) 假设检验假设检验 更一般的假设 H0: Cb b = f VS H1: Cb b f情形 A 和 B 都是情形 C 的特殊情形y = Xb b + , Nn (0,s2 In)Ab b = c, A: q p (q p), c: q 1. Rank(A) = q, 在某些限制 Ab b = c 下的估计;没有限制下的估计.22,ssb bb bHH定理定理 1 在有限制的模型下, 我们可得(i)(ii) 其中111()
29、 () ()() () ()( ) ().HHHHQQQQQQX XA A X X AcAAcA X X AAc且bbbbbbbbbbbbbb假设检验理论:定理定理 2 记 如定理 1 所示. 则对于 , 我们有以下结论:(i)(ii)(iii) 当 H0 为真, 检验统计量为( ), (), HHHQQQQ其中和bbbbbbbb)( )()()()( )()(11211cAAAcAcAAAcAb bb bb bb bGGqQQEGGQQHHspnqHHFqsQQpnQqQQF,2)/(/ )(证明可参 Seber (1977).假设检验理论(续):情形情形 A H0: b1 = = bk =
30、 0 VS H1: 某 bj 0, 1 j k记 A = (0, Ip), c = 0, 则:*1, (,)pb,bAc其中bbbbbb从定理 2 可得, 检验统计量为:pnpGGGGFspksF, 12*2*) 1(b bb bb bb bSS1111111111 () ().GGGGnGGGGGGGGxxxG GxG GSSSSAASS其中且情形情形 A (续续) H0: b1 = = bk = 0 VS H1: 某 bj 0, 1 j k 其检验常分解为 ANOVA 表:niiniiiniiyyyyyy121212)()()(SSTotalSSResSSRegRegReg22Total ,(1)FRpsSSSSSS则为确定性系数情形情形 A (续续) ANOVA 表: ( k = p - 1)22221 11 .x-xnpRFRpRFR可简单证明即 和等价,因为为单调函数情形情形 B H0: bj = 0 VS H1: bj 0记 A = (0,0,1,0,0) ej+1, 其中 1 位于第 j +1 个位置, c = 0. 则 jb bb bcA记(GG)-1 = C = (cij)A (GG) -1A -1 = cj+1,j+1-1121,121,21,11,1 HjijjjjjjnpjjjnpjjQQccFcstcsFt或b bb bb bb bb b