第十講(1)回歸正交設(shè)計

第七章回歸的正交設(shè)計

古典回歸是被動的，對試驗不做任何要求，回歸方程精度也很少研究。這不僅盲目增加了試驗次數(shù)，且數(shù)據(jù)信息不充分，達不到試驗?zāi)康摹?/p>

尋求最佳農(nóng)藝措施和模式化載培及建立數(shù)學模型，主動地把試驗安排、數(shù)據(jù)處理和回歸方程的精度統(tǒng)一起來考慮。即選擇試驗點使數(shù)據(jù)含有最大信息、減少試驗次數(shù)、數(shù)據(jù)分析具有較好的性質(zhì)。這就是“回歸設(shè)計與分析”研究的問題。

§7.1什么是回歸設(shè)計1

如糧食生產(chǎn)農(nóng)藝最優(yōu)化,在不完全了解生產(chǎn)過程物理、化學和生物原理的情況下,回歸分析是一個有效方法。

類似“黑箱”理論，輸入和輸出已知，但內(nèi)部結(jié)構(gòu)不清楚。輸入是一些因子z1,z2,…，zp(如N肥、P肥、密度等)，輸出就是要最優(yōu)化的農(nóng)藝指標y(如產(chǎn)量)，當然可以有二個或更多的農(nóng)藝指標。

y總可以表示為因子z1,z2,…，zp的一個函數(shù)，即

y=f(z1,z2,…,zp)并且可以用一個多項式來近似表達這個函數(shù)。2

以因子z1,z2,…，zp

為坐標的空間稱為因子空間，尋找最優(yōu)化農(nóng)藝措施，也就是要在因子空間中尋找y為最大值或最小值的點。

如何在因子空間中選擇適當?shù)脑囼烖c，使得能以最快的速度建立一個有效的多項式，從而解決最優(yōu)化問題呢？這就要用到回歸的設(shè)計與分析這一工具。

回歸設(shè)計內(nèi)容相當豐富，有回歸正交、旋轉(zhuǎn)、Ｄ-最優(yōu)設(shè)計等。主要介紹這三種設(shè)計。按照回歸模型的次數(shù)，回歸設(shè)計又分為一次回歸設(shè)計、二次回歸設(shè)計等，首先討論一次回歸正交設(shè)計。

3

一次回歸的正交設(shè)計主要是運用二水平正交表（如L4(23),L8(27),L12(211)等）進行設(shè)計，設(shè)計與分析的主要步驟和性質(zhì)如下：1.確定因子的變化范圍§7.2一次正交回歸設(shè)計為因子zj的零水平，

如要研究p個因子z1,z2,…，zp與某項指標y的數(shù)量關(guān)系，那么首先要確定每個因子zj變化的下界和上界，假如試驗就在水平z1j和z2j上進行，那么分別稱z1j和z2j為因子zj的下水平和上水平，并稱它們的算術(shù)平均41.確定因子的變化范圍為因子zj的變化區(qū)間。

如要研究p個因子z1,z2,…，zp與某項指標y的數(shù)量關(guān)系，那么首先要確定每個因子zj變化的下界和上界，假如試驗就在水平z1j和z2j上進行，那么分別稱z1j和z2j為因子zj的下水平和上水平，并稱它們的算術(shù)平均為因子zj的零水平，它們差的一半5例7.1

為了研究小麥高產(chǎn)栽培技術(shù)，選擇影響小麥產(chǎn)量三個主要因素水分狀況z1（全生育期土壤濕度占田間持水量％,75~95）、肥料z2(公斤/畝，20~40)、密度z3(萬株/畝，45~65)進行回歸正交試驗。

z1:z11=75,z12=95,z01=(75+95)/2=85,1=(95-75)/2=10z2:z12=20,z22=40,z02=(20+40)/2=30,2=(40-20)/2=10z3:z13=45,z23=65,z03=(45+65)/2=55,3=(65-45)/2=1062.對每個因子zj的水平進行編碼編碼就是對因子的取值作如下線性變換：這樣就建立了因子zj與xj取值的一一對應(yīng)關(guān)系:下水平z1j←→-1零水平z0j←→0上水平

z2j←→+1

具體編碼工作可在表上進行。7

從上面知道，當zj在區(qū)間［z1j,z2j］內(nèi)變化時，它的編碼值xj就在區(qū)間［-1,+1］內(nèi)變化。因子水平編碼表因子

z1

z2……zp

下水平(-1)零水平(0)上水平(+1)

z11

z12……z1pz01

z02……z0pz21

z22……z2p變動區(qū)間(j)

1

2

……p

對因子zj編碼后。y對z1,z2,…,zp的回歸問題，就轉(zhuǎn)化為y對x1,x2,…,xp的回歸問題。因此可以在x1,x2,…,xp為坐標軸的編碼空間中選擇試驗點進行回歸設(shè)計。8對例7.1的因子水平編碼如下表。因子水平編碼表因子

z1

z2

z3

下水平(-1)零水平(0)上水平(+1)

752045853055954065變動區(qū)間(j)

10101093.選擇適當?shù)亩秸槐?/p>

運用二水平正交表時，需用“-1”代換通常的二水平正交表中“2”。代換后表中的“+1”既表示因子水平的不同狀態(tài)，也表示因子水平變化的數(shù)量大小。經(jīng)代換后，正交表的交互作用列還可直接由表中相應(yīng)幾列的對應(yīng)元素相乘而得到，因此交互列也就不必要了。如L4(23)試驗號

x1

x2

x1x21234

1111-1-1-11-1-1-1110

究竟用哪一張二水平正交表，要根據(jù)因子的個數(shù)而定。正交表確定以后，把各變量放入正交表的某些列上，把這些列取出就組成了一張試驗計劃。如把三個變量放在正交表L8(27)的x1,x2,x3列上，就得到一張如下的全因子試驗計劃。試驗號

x1

x2

x312345678

11111-11-111-1-1-111-11-1-1-11-1-1-111

用二水平正交表編制的試驗計劃具有正交性。顯然，若以xij表示第i次試驗中第j個變量的編碼值，于是在試驗計劃中有任一列的和、任兩列的內(nèi)積等于0，即試驗號

x1

x2

x312345678

11111-11-111-1-1-111-11-1-1-11-1-1-1

由此稱這種設(shè)計為正交設(shè)計，所得的的試驗計劃為正交計劃。

12對例7.1,選L8(27)正交表,根據(jù)因子水平編碼表,“因素順序上列,水平對號入座”,每行為一個試驗處理,得試驗計劃如下:試驗號

x1(z1)x2(z2)x3(z3)12345678

1(95)1(40)1(65)11-1(45)1-1(20)11-1-1-1(75)11-11-1-1-11-1-1-1試驗計劃134.回歸系數(shù)的計算與統(tǒng)計檢驗

根據(jù)正交設(shè)計進行n次試驗,試驗結(jié)果為y1,y2,…,yn。則一次回歸的數(shù)學模型為

于是參數(shù)的最小二乘估計yi=0+1xi1,+2xi2+…+pxip+ei，i=1,2,…,nb=(XTX)-1XTY=A-1B由正交設(shè)計的性質(zhì)知，A-1=diag(1/n,1/n,…,1/n)B=XTY=(B0,B1,…,Bp)T，B0=Σyi，Bj=Σxijyi，j=1,…,p從而b0=Σyi/n，bj=Σxijyi/n，j=1,…,p(1)回歸系數(shù)的計算

14由于COV(b)=2(XTX)-1=2diag(1/n,1/n,…,1/n),所以回歸系數(shù)bj之間不相關(guān)。若某個回歸系數(shù)不顯著,可將該項直接剔除,不需重新計算其它回歸系數(shù)。

例7.1進行試驗后得結(jié)果如下:試驗號

x1(z1)x2(z2)x3(z3)試驗結(jié)果y12345678

1(95)1(40)1(65)2.111-1(45)2.31-1(20)13.31-1-14.0-1(75)115.0-11-15.6-1-116.9-1-1-17.815b0=Σyi/n=(2.1+…+7.8)/8=37/8=4.625試驗號

x1(z1)x2(z2)x3(z3)試驗結(jié)果y12345678

1(95)1(40)1(65)2.111-1(45)2.31-1(20)13.31-1-14.0-1(75)115.0-11-15.6-1-116.9-1-1-17.8b1=Σxi1yi/n=[12.1+…+(-1)7.8]/8=-13.6/8=-1.7同理得b2=-0.875，b3=-0.3，則回歸方程為16(2)回歸方程的檢驗

平方和當FF(p,n-p-1)時，回歸方程顯著。對例7.1，lyy=(2.12+…+7.82)-(2.1+…+7.8)2/8=30.275u=8(1.72+0.8752+0.32)=29.965，Q=30.275-29.965=0.31F=(29.965/3)/(0.31/4)=128.88>16.69=F0.01(3,4)所以，回歸方程極顯著。17(3)回歸系數(shù)的檢驗

當FjF(1,n-p-1)時，回歸系數(shù)bj顯著。對例7.1：>21.2=F0.01(1,4)>21.2=F0.01(1,4)>7.71=F0.01(1,4)回歸系數(shù)都顯著。185.零水平的重復(fù)試驗

一次回歸方程顯著，也不能保證一次回歸模型是最好的。為此在零水平(z01,z02,…,z0p)安排一些重復(fù)試驗，如安排m0次重復(fù)試驗，其試驗結(jié)果為y01,y02,…,y0m0。

模型檢驗(失擬檢驗)，即檢驗：零水平試驗結(jié)果的均值與回歸方程中常數(shù)項b0是否有顯著差異。記若則在水平下，試驗中心一次回歸與實測值擬合較好。19若則在水平下，試驗中心一次回歸與實測值擬合較好。否則一次回歸擬合不好，需要更高次模型。

對例7.1，在零水平重復(fù)2次試驗，m0=2,y01=4.5,y02=4.3，S0=0.02，Q=0.31，b0=4.625，n=8，p=3

|t|=1.1078<2.571=t/2(5)，所以，一次回歸方程與實測值還是擬合得較好的，此回歸方程是適合的。

20

對例7.1，在零水平重復(fù)2次試驗，m0=2,y01=4.5,y02=4.3，S0=0.02，Q=0.31，b0=4.625，n=8，p=3

|t|=1.1078<2.571=t/2(5)，所以，一次回歸方程與實測值還是擬合得較好的，此回歸方程是適合的。

將編碼公式

x1=(z1-85)/10，x2=(z2-30)/10，x3=(z3-55)/10，代入回歸方程

得原始變量的回歸方程

216.一次回歸正交設(shè)計的旋轉(zhuǎn)性

設(shè)用正交設(shè)計所求得的回歸方程為那么回歸值的方差

而是p維編碼空間內(nèi)的一個球面，球心在原點，半徑為。

位于同一球面上的點的預(yù)測值y的方差是相等的。這個性質(zhì)稱為旋轉(zhuǎn)性。22

預(yù)測值的方差相等，可以直接比較各預(yù)測值的好壞，從而找出預(yù)測值相對較優(yōu)的區(qū)域。

在古典回歸分析中，由于預(yù)測值的方差強烈地依賴于點的位置，所以就不能這樣做。今后我們還要專門研究具有旋轉(zhuǎn)性的設(shè)計——旋轉(zhuǎn)設(shè)計。

顯然用回歸方程作預(yù)報時,預(yù)報值y的方差愈小,預(yù)報就愈準確,方差愈大,預(yù)報值的誤差就愈大。有了旋轉(zhuǎn)性以后，這種預(yù)報的誤差程度可以簡單地用

來表示,

大(即點離球心遠)誤差就大,

小誤差也就小。

23

二次回歸的正交設(shè)計要比一次回歸的正交設(shè)計復(fù)雜一些，p個變量的二次回歸方程§7.3二次正交回歸設(shè)計

為計算二次回歸方程的系數(shù)，每個變量的水平數(shù)應(yīng)大于3，因而做的試驗次數(shù)往往比較多。如當p=4時，三水平全因子試驗次數(shù)是81次，它比4個變量的二次回歸的系數(shù)C24+2=15要多4倍以上，以致剩余自由度過大。

共有回歸系數(shù)q=1+C1p+C2p+C1p=C2p+2個，為了得到二次回歸方程，試驗次數(shù)n當然應(yīng)不小于q。

24

采用“組合設(shè)計”。所謂組合設(shè)計，就是選擇幾類具有不同特點的點組合起來形成試驗計劃。以p=2為例，說明組合設(shè)計中試驗點在因子空間中的分布。

在二個變量x1,x2場合下，組合設(shè)計由n=9個點組成:x1x211-1-1(0,-)(-,0)(,0)(0,)(1,1)(1,-1)(-1,-1)(-1,1)(0,0)25試驗號

x1

x21(11)2(1-1)3(-11)4(-1-1)5(0)6(-0)7(0)8(0-)9(00)同理，可以寫出p=3的試驗計劃。組成兩水平(+1和-1)全因子試驗22分布在x1和x2軸上的星號位置由x1和x2的零水平組成的中心試驗點

26一般p個變量的組合設(shè)計由下列n個點組成：n=mc+2p+m02p—分布在p個坐標軸上的星號位置點，它們與中心的距離稱為星號臂，是待定參數(shù)。根據(jù)一定的要求調(diào)節(jié)

，可得到各種具有很好性質(zhì)的設(shè)計（如正交設(shè)計、旋轉(zhuǎn)設(shè)計）。

m0—在各變量都取零水平的中心點的重復(fù)試驗次數(shù)。它可以只做一次，也可以重復(fù)二次或多次。

其中：mc—二水平(+1和-1)的全因子試驗的試驗點個數(shù)2p

，或部分實施的試驗點個數(shù)2p-1，2p-2等。

27組合設(shè)計的優(yōu)點：

試驗點比全因子試驗要少得多，但卻仍保持足夠的剩余自由度(見表7-12）；它在一次回歸基礎(chǔ)上獲得，如果一次回歸不顯著，只要在一次回歸試驗的基礎(chǔ)上，再在星號點和中心點補充做一些試驗，就可求得二次回歸方程

要使組合設(shè)計成為正交設(shè)計，還要確定適當?shù)男翘柋?。如p=2時，二次回歸組合設(shè)計的設(shè)計矩陣為

28

要使組合設(shè)計成為正交設(shè)計，還要確定適當?shù)男翘柋邸Ｈ鏿=2時，二次回歸組合設(shè)計的設(shè)計矩陣為

試驗號

x0x1

x2x1x2

x12

x221111111211-1-11131-11-11141-1-111151

002061-0020710002810-0029000000

對一次計劃加入星號點后，并不破壞一次變量和交互效應(yīng)的正交性。只是被x0和xj2破壞了。

29

為使組合設(shè)計具有正交性，必須在使得矩陣C=(XX)-1為對角陣的條件下定出的值。可以證明：當mc=2P時（全因子試驗情況），有

4+2P2-2P-1(p+0.5m0)=0當mc=2P-1時（部分實施情況），有

4+2P-12-2P-2(p+0.5m0)=0

對于給定的p和m0，就可計算2的值。常用的2值已計算在表7-14上。然后對平方列進行中心化即可。

30

對于給定的p和m0，就可計算2的值。常用的2值已計算在表7-14上。然后對平方列進行中心化即可。

表7-142值表m0p2345(1/2實施)123456789101.0001.4762.0002.391.1601.6502.1982.581.3171.8312.3902.771.4752.0002.5802.951.6062.1642.7703.141.7422.3252.9503.311.8732.4813.1403.492.0002.6333.3103.662.1232.7823.4903.832.2432.9283.6604.0031

運用組合設(shè)計進行二次回歸正交設(shè)計的步驟和統(tǒng)計分析分述如下：§7.4二次正交回歸設(shè)計的統(tǒng)計分析1.確定因子的變化范圍變化區(qū)間為

如某研究中有p個因子z1,z2,…，zp，其第j個因子zj變化的下界和上界分別為z1j和z2j(j=1,2,…,p)，零水平為32

例2

玉米高產(chǎn)栽培試驗，選擇三個主要因子密度z1：上下界z21=4353株/畝，z11=1647株/畝；施N量z2：z22=16.75公斤/畝，z12=3.25公斤/畝；施P量z3：z23=8.4公斤/畝，z13=1.6公斤/畝。采用二次回歸正交設(shè)計，取m0=3，查表得=1.353，零水平與變動區(qū)間分別為：同理得：z02=10，z03=5同理得：2=2，3=2.5332.編制因子水平編碼表與一次回歸設(shè)計類似，對因子作如下線性變換：則有因子水平編碼表xj因子z1z2zp10-1-z21

z22……z2pz01+1

z02+2……z0p+pz01

z02……z0pz01-1

z02-2……z0p-pz11

z12……z1p變化區(qū)間j1

2……p34xj因子密度z1(株/畝)施N量z2(公斤/畝)施P量z3(公斤/畝)1.35310-1-1.3534353400030002000164716.75151053.258.47.552.51.6變化區(qū)間j100052.5對例2

z21=4353，z11=1647；z22=16.75，z12=3.25；z23=8.4，z13=1.6；=1.353，z01=3000，z02=10，

z03=51=1000，2=2，3=2.5因子水平編碼表35對例7.2，試驗方案與試驗結(jié)果如下：3.選擇相應(yīng)的組合設(shè)計

進行n=mc+2p+m0次試驗。如對于p=3,mc=23,m0=3的情況，由x1、x2、x3所占的列組成試驗計劃。xj因子密度z1(株/畝)施N量z2(公斤/畝)施P量z3(公斤/畝)1.35310-1-1.3534353400030002000164716.75151053.258.47.552.51.6變化區(qū)間j100052.5因子水平編碼表36x1(z1)x2(z2)x3(z3)y-1(2000)-1(5)-1(2.5)1136-1-11(7.5)1088-11(15)-11128-11111971(4000)-1-112811-11134911-110751111275xj密度z1(株/畝)施N量z2(公斤/畝)施P量z3(公斤/畝)1.35310-1-1.3534353400030002000164716.75151053.258.47.552.51.6mc=23=8個試驗點?！皩μ柸胱?7x1(z1)x2(z2)x3(z3)y-1.353(1647)0(10)0(5)11801.353(4353)0013440(3000)-1.353(3.25)0120101.353(16.75)0116500-1.353(1.6)1102001.353(8.4)1218000123600012060001221xj密度z1(株/畝)施N量z2(公斤/畝)施P量z3(公斤/畝)1.35310-1-1.3534353400030002000164716.75151053.258.47.552.51.62p+m0=23+3

=9個試驗點。38x1(z1)x2(z2)x3(z3)y-1(2000)-1(5)-1(2.5)1136-1-11(7.5)1088-11(15)-11128-11111971(4000)-1-112811-11134911-110751111275-1.353(1647)0(10)0(5)11801.353(4353)0013440(3000)-1.353(3.25)0120101.353(16.75)0116500-1.353(1.6)1102001.353(8.4)