胡純嚴 ,胡良平 ,2*
(1.軍事科學院研究生院,北京 100850;2.世界中醫(yī)藥學會聯(lián)合會臨床科研統(tǒng)計學專業(yè)委員會,北京 100029*通信作者:胡良平,E-mail:lphu927@163.com)
正交設計是安排多因素的一種高效方法,在滿足某些特定條件時,正確地運用正交設計,不僅可以獲得比較精準的試驗結果,而且所需的樣本含量較少。本文將介紹正交設計的基本概念、設計方法、計算公式以及基于SAS軟件實現(xiàn)正交設計定量資料方差分析的方法。
在幾何學上,正交的意思是兩條直線互相垂直。在代數(shù)學上,若兩個向量的內(nèi)積(或稱數(shù)量積)等于0,就定義這兩個向量是互相正交的[1]。假設有兩個向量如下:a=(-1,-1,-1,-1,1,1,1,1),b=(-1,1,-1,1,-1,1,-1,1)。將它們對應位置上的數(shù)字相乘后求和,可記為(a,b)=-1×(-1)+(-1)×1+……+1×1=0,它就是向量a與向量b的“內(nèi)積”。于是,在數(shù)學上,稱上面兩串數(shù)字排成的序列或向量是互相正交的。
在統(tǒng)計學上,統(tǒng)計學家構造出一系列類似于上面給出的“向量”,將它們排列在一張表格中,表格中任何兩列數(shù)據(jù)的“內(nèi)積”都等于0,這種表格就叫做正交表[2]。以L8(27)正交表為例,其中“L”代表正交表,數(shù)字“8”代表此表有8行,數(shù)字“2”代表每列有2個不同的水平,數(shù)字“7”代表此表有7列。L8(27)正交表見表1。表1共有7列,每列可以代表一個試驗因素的2個水平各重復出現(xiàn)4次的一種排列;實際使用時,每列可以安排一個二水平因素,每行代表全部因素的一種水平組合結果,例如,第1行所有因素都取“1水平”。若用“-1”取代表中的“2”,則表1中的任何兩列的“內(nèi)積”都等于0。
表1 L8(27)正交表Table 1L8(27)Orthogonal table
正交設計是利用一系列規(guī)格化的正交表來安排各試驗因素及其水平組合的方法。正交表是經(jīng)過嚴格的數(shù)學推導編制而成,正交表中的每一行代表全部試驗因素的一種水平組合,稱為一個試驗點;正交表中的每一列代表一個試驗因素(或交互作用項)及其水平的排列和重復出現(xiàn)的情況,將全部因素及其擬考察的交互作用項恰當?shù)匕才旁谡槐淼谋眍^上,就叫做正交設計。
正交設計有4個突出特點:①由正交表挑選出來的試驗點在空間上具有“均勻分散性”;②由正交表挑選出來的試驗點在統(tǒng)計分析時具有“整齊可比性”;③某些未包括在正交表中的試驗點,可以通過統(tǒng)計分析將其發(fā)現(xiàn);④在某些假設(例如某些因素之間的交互作用不存在或可以忽略不計,特定試驗條件下試驗結果的穩(wěn)定性很好)成立的條件下,正交設計比其他多因素設計(排除均勻設計)所需的樣本含量更少。
在采用正交設計安排試驗時,應預先考慮好如何估計試驗誤差:在正交表中留有空白列(可用于估計不同試驗條件之間所產(chǎn)生的試驗誤差,被稱為第一類誤差),或者在挑選出來的試驗點上進行m次獨立重復試驗(m≥2,可用于估計相同試驗條件下不同重復試驗次數(shù)之間所產(chǎn)生的試驗誤差,被稱為第二類誤差)。一般來說,第一類誤差大于第二類誤差。在對資料進行方差分析時,若只有第一類誤差,就以其均方為分母,構造F檢驗統(tǒng)計量;若有兩類誤差,就需要將它們的離均差平方和、自由度分別進行合并,用合并后的離均差平方和除以合并后的自由度,構造出誤差的均方。誤差項的自由度越大(意味著所需要的樣本含量越大),所得結論的可信度就越高。
正交表可分為兩大類:第一類是飽和正交表,第二類是非飽和正交表,也稱不完備正交表[3]。在不做重復試驗的飽和正交表中,總自由度=各列自由度之和=表的行數(shù)-1,各列自由度=各列中不同水平數(shù)-1??傋杂啥炔坏扔诟髁凶杂啥戎偷恼槐恚褪欠秋柡驼槐?。例如,在L8(27)正交表中,df總=8-1=7,每一列的自由度=2-1=1,7列自由度之和=7,滿足“總自由度=各列自由度之和”的要求,故L8(27)正交表是飽和正交表。
根據(jù)正交表中各列的水平數(shù)是否相等,還可將正交表分為同水平正交表與混合水平正交表[4-5]。常見的同水平正交表有:二水平的正交表[L4(23)、L8(27)、L16(215)、L32(231)、L64(263)],三水平的正交表[L9(34)、L27(313)、L81(340)],四水平的正交表[L16(45)、L64(421)],五水平的正交表[L25(56)、L125(531)]。常見混合水平的正交表有:L8(4×24)、L12(3×24)、L16(42×29)、L16(43×26)、L18(2×37)、L18(6×36)、L27(9×39)、L32(45×216)、L36(6×312)、L32(16×216)、L24(3×4×24)、L36(62×35×2)、L32(8×46×26)。
有些混合水平正交表是非飽和正交表,例如,在L12(3×24)正交表中,總自由度=12-1=11,而各列自由度之和=(3-1)+4×(2-1)=6。
正交表的種類很多,且因素的個數(shù)、水平數(shù)以及需要考察的交互作用項數(shù)都會隨著具體問題發(fā)生變化,還涉及是否進行重復試驗。因此,無法給出一個統(tǒng)一的正交設計定量資料一元方差分析計算公式。現(xiàn)假定在L8(4×24)的前3列上分別安排因素A(4水平)、因素B(2水平)和因素C(2水平),扼要介紹正交設計定量資料一元方差分析的計算思路。計算正交表中每一列上的離均差平方和的方法,與單因素兩水平或多水平設計定量資料一元方差分析計算組間離均差平方和的方法相同,可參見文獻[6]。
第一步,計算最后一列試驗結果Y的總離均差平方和。
第二步,計算A、B、C列試驗結果Y的總離均差平方和。
第三步,計算誤差E的離均差平方和。
第四步,計算總自由度和因素A、B、C的自由度以及誤差E的自由度。
第五步,計算誤差的均方和因素A、B、C的均方。
第六步,計算與因素A、B、C對應的檢驗統(tǒng)計量F的值。
若基于手工計算,需根據(jù)計算檢驗統(tǒng)計量F值的分子和分母的自由度去查方差分析用的F臨界值表,再比較檢驗統(tǒng)計量F值與其對應的F臨界值的大小,做出統(tǒng)計推斷。因篇幅所限,此處從略。本文將基于SAS軟件實現(xiàn)方差分析[7]。
【例1】某化工廠為了提高產(chǎn)品的收率,根據(jù)具體情況和經(jīng)驗確定試驗條件如下:因素A為反應溫度(A1=30℃,A2=40℃);因素B為反應時間(B1=60 min,B2=90 min);因素C為兩種原料比(C1=1∶1,C2=1.5∶1);因素D為攪拌速度(D1=慢,D2=快)。選用L8(27)正交表安排上述4個因素,試驗結果見表2[4]。試分析哪些因素對收率Y的影響有統(tǒng)計學意義。
表2 4個試驗因素對產(chǎn)品收率影響的試驗結果Table 2 Experimental results of the influence of four experimental factors on the product yield
【例2】某電解腐蝕試驗考察4個三水平因素(各因素的具體含義從略),交互作用均不考慮。試驗指標為產(chǎn)品質(zhì)量,按照規(guī)定標準進行綜合評分。選用L9(34)正交表安排試驗,在每個試驗條件下進行3次重復試驗。試驗安排和試驗結果見表3[3]。試分析4個因素對試驗結果的影響是否有統(tǒng)計學意義。
表3 4個三水平因素對電解腐蝕作用的正交設計及試驗結果Table 3 Orthogonal design and the experimental results of four three-level factors on the role of the electrolytic corrosion
【例3】某農(nóng)業(yè)科研所為了提高小麥產(chǎn)量,擬考察A、B、C三個因素對小麥畝產(chǎn)的影響。各因素各水平的含義如下:因素A為品種,有4個水平,1~4分別代表泰山4號、濟6矮、矮稈早、139號;因素B為施肥量,有2個水平,1和2分別代表高肥和低肥;因素C為播種期,有2個水平,1和2分別代表早播和晚播。試驗安排和試驗結果見表4[5]。試分析三個因素對小麥畝產(chǎn)Y的影響是否有統(tǒng)計學意義?
表4 3個因素對小麥畝產(chǎn)影響的正交設計及試驗結果Table 4 Orthogonal design and the experimental results of the influence of three factors on the wheat yield per mu
3.2.1 對例1的分析與解答
【分析與解答】所需要的SAS程序如下:
【SAS程序說明】在SAS過程步中,有四個“model語句”,第一個是可執(zhí)行的,后三個放在注釋語句“/**/”中是不可執(zhí)行的。每次運行SAS程序,只能讓一個“model語句”處于可執(zhí)行狀態(tài)。
【SAS輸出結果及解釋】
由第一個“model語句”的輸出結果可知,4個因素均無統(tǒng)計學意義,其中,因素D的作用最小。因篇幅所限,具體輸出結果從略,下同。
由第二個“model語句”的輸出結果可知,因素C(F=11.54,P=0.042 6)和交互作用項A*B(F=11.54,P=0.04 26)具有統(tǒng)計學意義,而因素A(F=2.88,P=0.188 0)和因素B(F=2.88,P=0.188 0)均無統(tǒng)計學意義。
由第三個“model語句”的輸出結果可知,因素C(F=11.54,P=0.042 6)具有統(tǒng)計學意義,而交互作用項A*B(F=5.77,P=0.092 0)無統(tǒng)計學意義。
由第四個“model語句”的輸出結果可知,當模型中僅有因素C時,它對定量觀測結果的影響變得無統(tǒng)計學意義(F=3.41,P=0.114 4)。
由以上分析過程和結果可知:多因素方差分析結果是相對的,不是一成不變的。因為隨著模型中待分析的項的改變(在本質(zhì)上,是模型誤差的均方在改變),各因素及交互作用項的統(tǒng)計顯著性也會發(fā)生變化。
3.2.2 對例2的分析與解答
【分析與解答】所需要的SAS程序如下:
【SAS輸出結果及解釋】
由第一個過程步的輸出結果可知,因素C(F=0.54,P=0.591 0)對試驗結果的影響無統(tǒng)計學意義,因素A(F=11.20,P=0.000 7)、因素B(F=5.17,P=0.016 9)、因素D(F=5.01,P=0.018 6)對試驗結果的影響均有統(tǒng)計學意義。因篇幅所限,具體輸出結果從略,下同。
由第二個過程步的輸出結果可知,因素A(F=11.74,P=0.000 4)、因素B(F=5.41,P=0.013 2)、因素D(F=5.25,P=0.014 7)對試驗結果的影響均有統(tǒng)計學意義。
因為結果的平均值越大越好(為節(jié)省篇幅,4個因素各水平下的均值計算結果從略),故因素A應取第三個水平,因素B應取第二個水平,因素C應取第三個水平,因素D應取第三個水平。
3.2.3 對例3的分析與解答
【分析與解答】所需要的SAS程序如下:
【SAS程序說明】在SAS過程步中,有三個“model語句”,第一個是可執(zhí)行的,后兩個放在注釋語句“/**/”中是不可執(zhí)行的。每次運行SAS程序,只能讓一個“model語句”處在可執(zhí)行狀態(tài)。
【SAS輸出結果及解釋】
由第一個過程步的輸出結果可知,僅因素C(F=154.27,P=0.006 4)對產(chǎn)量的影響具有統(tǒng)計學意義。
由第二個過程步的輸出結果可知,在淘汰掉P值最大的因素B之后,仍然只有因素C(F=61.39,P=0.004 3)對產(chǎn)量的影響具有統(tǒng)計學意義。
由第三個過程步的輸出結果可知,在淘汰掉因素A和因素B之后,因素C(F=23.98,P=0.002 7)對產(chǎn)量的影響仍有統(tǒng)計學意義。
標準正交表對因素的水平數(shù)有嚴格要求,但在某些實際問題中,很難完全與之符合,此時,就需要對原先的正交表進行改造。例如,多個因素有3個水平,少數(shù)因素有2個水平,又沒有相應的混合水平正交表可以使用。此時,可采用“擬水平法”將二水平因素改造成三水平因素[5]。在有的試驗中,根據(jù)研究者經(jīng)驗,在安排方案時已知某些因素之間存在依賴關系,一個因素用量或水平的選取需要隨另一個因素用量或水平而定,這時就需要采用“活動水平法”[5];在有的試驗中,按正交表做了一些試驗后,用趨勢圖分析發(fā)現(xiàn)某一因素對試驗指標的影響存在某種變化趨勢,研究者認為有必要對該因素添加新的水平并追加幾次試驗,以便對它的影響有更全面的了解,此時,就需要采用“部分追加法”[5]。
值得注意的是,對標準正交表進行改動后,對定量資料進行方差分析的方法也需相應改變,否則,會影響計算結果的精確度;當正交表中的評價指標或結果變量有多個時,若它們之間存在一定的關聯(lián),就需要同時對它們進行分析。常規(guī)的統(tǒng)計分析方法是正交設計定量資料多元方差分析[8],還有一種統(tǒng)計分析方法是“廣義方差分析”[5]。因篇幅所限,本文不予贅述。
本文介紹了正交設計的基本概念、具體實施方法以及定量資料一元方差分析的計算公式,通過三個實例,展示了二水平、三水平和混合水平正交設計的具體實施方法;還給出了無重復試驗和有重復試驗定量資料方差分析的SAS實現(xiàn)方法。