韓 冰，唐 碩

（西北工業(yè)大學(xué) 航天學(xué)院，西安 710072）

一種基于正交表的工程無混雜試驗設(shè)計法

韓 冰，唐 碩

（西北工業(yè)大學(xué) 航天學(xué)院，西安 710072）

文章提出了一種靈巧的試驗設(shè)計法，它能根據(jù)工程實際中因素間的交互作用，在正交表中選擇恰當(dāng)?shù)牧凶鳛橹饕蛩亓?，使所有交互作用列與主因素列無混雜。此外，在設(shè)計過程中分別引入動態(tài)規(guī)劃法和微分進(jìn)化法對各主因素列進(jìn)行了均勻性優(yōu)化，并根據(jù)數(shù)值實驗結(jié)果，得出兩種方法的優(yōu)化效果隨因素數(shù)和因素水平數(shù)的變化規(guī)律。

正交設(shè)計；交互作用；均勻性；動態(tài)規(guī)劃；微分進(jìn)化

0 引言

近年來，試驗設(shè)計理論在工程實際中得到了廣泛應(yīng)用。為評價設(shè)計表的統(tǒng)計推斷能力，提出了許多優(yōu)良性準(zhǔn)則。其中，以正交性[1]為設(shè)計準(zhǔn)則的正交設(shè)計與以均勻性[2]為準(zhǔn)則的均勻設(shè)計都獲得了成功。Ma等[3]進(jìn)一步研究了正交性與均勻性的關(guān)系，為正交表的均勻性優(yōu)化提供了理論依據(jù)。

目前，試驗設(shè)計技術(shù)逐漸應(yīng)用于飛行仿真中，如選擇多種飛行器方案，分析干擾因素的影響，優(yōu)化設(shè)計參數(shù)等。飛行仿真試驗的因素多，因素水平高，因素間交互作用強(qiáng)。于是，設(shè)計的分辨力變得十分重要。Box和Hunter[4]對正規(guī)設(shè)計提出了最大分辨力準(zhǔn)則，Ma等[5]進(jìn)一步提出廣義最小低階混雜（GMA），將該準(zhǔn)則推廣到非正規(guī)設(shè)計情形。Fang等[6]描述了2水平設(shè)計的均勻性與混雜的解析關(guān)系。最近，Sun等[7]提出了基于GMA準(zhǔn)則的一種析因設(shè)計。以上設(shè)計法均假設(shè)主因素間的交互作用未知，而在工程實際中，因素間存在的所有顯著的交互作用往往是已知的。這就迫切的需要一種已知交互作用與主因素?zé)o混雜的設(shè)計，能夠排除不存在的交互作用對試驗表設(shè)計過程的干擾，從而保證在實際無混雜的同時最大限度地減小試驗次數(shù)。在正規(guī)正交表的構(gòu)造算法[8]的基礎(chǔ)上，本文提出了一種靈巧的設(shè)計程序，能根據(jù)用戶輸入的交互因素表，自行選取相應(yīng)正交表的恰當(dāng)列作為主因素列，使設(shè)計表能分辨出所有的主因素效應(yīng)與交互效應(yīng)，并對主因素列進(jìn)行均勻性優(yōu)化。

1 交互作用列表設(shè)計法

交互作用列表設(shè)計法，是在正交表中選取主因素列并計算交互作用列的方法，它必須保證所有主因素列不與任一交互作用列混雜。

1.1 正交表的構(gòu)造及交互列的定義

本文涉及的正交表均為Ln(Sm)型正規(guī)正交表。馬希文[8]詳細(xì)描述了它的構(gòu)造方法，并給出了交互作用列的定義，均為本文采用。此種正交表構(gòu)造法為：設(shè)域Z中的元素為0,1,…s-1，Zk為定義在域Z中的所有k維向量組成的線性空間。按字典順序?qū)k中的全體向量推成一個Sk×k的矩陣C，從C中抽取第一個非零元素為1的行，可得到矩陣P0,將P0的列的順序顛倒過來，即得生成矩陣P。設(shè)計矩陣F=CPT。將矩陣F的每個元素加1，即得到Ln(Sm)型正規(guī)正交表。

生成矩陣P的行向量為F中對應(yīng)列的生成向量。于是給出交互作用列的定義：設(shè)D1,D2是F中兩個不同列，則D1,D2的生成向量形成的線性空間中，所包含的矩陣P中的其它行向量在F中的對應(yīng)列叫做D1,D2的交互作用列。推廣至n列情形：設(shè) D1,D2,…,Dn是 F中 n個不同列，則 D1,D2,…,Dn的生成向量形成的線性空間中，所包含的P中的所有行向量在F中的對應(yīng)列，除去所有比D1,D2,…,Dn低階的交互列，得交互作用列。其中，低階交互列包括：D1,D2,…,Dn中各因素兩兩組合的交互效應(yīng)列（稱2階交互列），3階交互列,…,n-1階交互列。

1.2 交互作用的錄入約定

為表述方便，本文將主因素分為兩類：已經(jīng)確定列號的主因素稱為舊因素，待選列號的主因素稱為新因素?？蓪⑺薪换プ饔镁幹瞥杀?下文稱錄入表)，錄入到計算程序中。錄入表的一行稱為一個錄入行，它是一個需計算的交互作用，參與交互作用的所有主因素編號記錄在該行的每一列上。計算程序按行依次讀取，選擇本行的新因素列并計算其交互作用列，同時進(jìn)行必要的校核，這一過程稱為一個計算循環(huán)。為方便計算，對交互表的因素排列順序作如下約定：編號較小的因素的首次錄入先于編號較大的因素；舊因素間的新交互作用，排在含新因素的交互作用前；同一行中舊因素在前，新因素在后。

例：若需計算因素 1×2,1×3,1×4,1×2×3,1×2×4,1×3×4 的交互效應(yīng)，約定錄入表格式形如表1：

1.3 選擇新因素列的原則

新因素列的選擇需要滿足三個原則：

（1）新因素列不能為舊因素列；

（2）新因素列不能為任一交互列；

（3）新因素列不能為任一交互作用的任一低階交互列。

明顯的，違反了第1，2個原則，將發(fā)生混雜現(xiàn)象。而違反原則3的設(shè)計沒有意義。事實上，我們考慮高階交互作用的影響，是基于一切較之低階交互作用都有顯著影響的情形，或者說，在一般情況下，設(shè)計高階交互作用列，必須先設(shè)計一切較之低階的交互作用列。

1.4 新因素列的校核

設(shè)計表中所有未指定因素和交互作用的列稱為賦閑列。賦閑列是滿足原則的必要非充分條件。必要性是明顯的，非充分性是因為在大多數(shù)情況下，并不需要考慮舊因素的所有組合的交互作用，這就使得賦閑列中存在這種 “潛在組合”形成的交互列（下文稱潛在列）。若選取賦閑列中的潛在列作為新因素所在列，將會有違反原則的可能；而將其從可選列中一概除去，會造成很大的浪費，因為并不是所有潛在列均違反原則，而且，在某次計算循環(huán)中違反原則的潛在列，并不是在之后的所有計算循環(huán)中必定違反原則。校核的實質(zhì)是將違反原則的潛在列從本次計算循環(huán)的可選列中除去，從而達(dá)到保留可行的潛在列的目的。潛在列違反原則的情形有兩類，第一類情形是：在錄入行帶有新因素的計算循環(huán)中，選定新因素列后，計算所得的交互作用列與某舊因素列混雜。若錄入行中同時含有新因素及某“潛在組合”中的舊因素，或同時含有兩個新因素，均可引發(fā)這類情形。此時計算程序需重新選擇賦閑列中的其余列并校核是否發(fā)生混雜。值得注意的是，計算程序一次只選定一個新因素列并校核。這就避免了一次校核兩個新因素列的情況。第二類情形是：在錄入行帶有舊因素間的某個新組合的計算循環(huán)中，計算所得的交互列與某舊因素列混雜。潛在列違反原則的錯誤事實上發(fā)生在之前最后一次選取新因素列的計算循環(huán)中。此時計算程序應(yīng)當(dāng)退回至發(fā)生錯誤的計算循環(huán)，并重新選擇其它賦閑列為主因素列。

1.5 交互作用列表設(shè)計算法與算例

完整的交互列表設(shè)計算法如下：

步驟1.初始化程序：主因素1、2所在列自動選為正交表的列 1、列 2。

步驟2.從錄入表錄入一行，本次計算循環(huán)開始。新因素備選列為所有賦閑列。

步驟3.讀入全部舊因素，若舊因素數(shù)大于1，則計算它們的交互作用列并校核。若校核失敗，程序轉(zhuǎn)至記錄的錄入表錄入狀態(tài)處并從備選列中剔除記錄的主因素號；校核成功則繼續(xù)；

步驟4.讀入新因素后記錄當(dāng)前錄入表錄入狀態(tài)（即：錄入表的行列號，新因素編號，舊因素編號，舊因素列號，舊交互作用列號）；

步驟5.將備選列中列號最小的選定為新因素列，并進(jìn)行校核。校核成功則記錄選定的主因素列號；校核失敗則將此列剔出備選列并重新選擇。若直至無備選列時仍校核失敗，則中斷程序，并要求增大正交表規(guī)模；

步驟6.若本行仍有新因素未錄入，則繼續(xù)讀入一個新因素并轉(zhuǎn)步驟4；否則進(jìn)行下一步；

步驟7.本次計算循環(huán)結(jié)束。若錄入表仍未讀完，則轉(zhuǎn)步驟-2；否則算法結(jié)束。

現(xiàn)在給出L81340型正規(guī)正交表的交互作用列表設(shè)計實例。取表1為錄入表，將其輸入程序后得各主因素列與交互列如表2所示。

依交互列的定義，可以很方便的檢驗計算結(jié)果的正確性。

圖1給出了程序設(shè)計流程示意圖。

2 設(shè)計的均勻性優(yōu)化法

2.1 均勻性優(yōu)化與均勻性準(zhǔn)則

方開泰[9]對正規(guī)正交表進(jìn)行了如下均勻性優(yōu)化：將正規(guī)正交表的列的水平作置換，得到 (S！)m個同構(gòu)的正交表，選取均勻性最好的表作為使用表。本文選取中心化L2-偏差（CD2）作為均勻性準(zhǔn)則。對于大規(guī)模的飛行仿真試驗，對同構(gòu)正交表逐一搜索的方法計算量巨大，以至無法完成。本文主要采用動態(tài)規(guī)劃法和微分進(jìn)化算法對正交表進(jìn)行優(yōu)化。

2.2 融合交互作用列表設(shè)計的一體化優(yōu)化方案

對正交表的均勻性優(yōu)化，本質(zhì)是對采樣點在試驗區(qū)域內(nèi)的均勻性優(yōu)化?；蛘哒f，是優(yōu)化主因素所在列，而對交互作用列的優(yōu)化沒有任何意義。故需將交互作用列表設(shè)計與均勻性尋優(yōu)相結(jié)合。

2.2.1 基于動態(tài)規(guī)劃的一體化優(yōu)化方案

表1 一個錄入表示例

表2 算例結(jié)果

表 3 L27313優(yōu)化效果對照表

表 4 L81340優(yōu)化效果對照表

表 5 L64421優(yōu)化效果對照表

表 6 L125531優(yōu)化效果對照表

本方案需對交互作用列表設(shè)計做以下改進(jìn)：首先，原設(shè)計中因素1,2所在列默認(rèn)為原正交表的1,2列。改進(jìn)方案為，因素1所在列仍為原表第1列，而在其余n-1列中順序取出每一列與第1列搭配，并對所選列作水平置換。選取使搭配的均勻性最佳的列的置換作為因素2所在列；其次，原設(shè)計中，當(dāng)讀入新因素時，選取列號最小的備選列作為新因素列。改進(jìn)方案為，在備選列中順序取出每一列與舊因素列搭配，并對所選列作水平置換。選取使搭配的均勻性最優(yōu)的列的置換作為新因素所在列；最后，在設(shè)計結(jié)束后對賦閑列進(jìn)行優(yōu)化：在賦閑列中順序取出每一列作水平置換，并與已選主因素列搭配，選取均勻性最好的，而后再在剩余賦閑列中按相同方法挑選出最優(yōu)列，直至賦閑列全部優(yōu)化完畢。

2.2.2 基于微分進(jìn)化的一體化優(yōu)化方案

本方案需對交互作用列表設(shè)計做以下改進(jìn)：首先，在不考慮交互作用的情況下，對原正交表進(jìn)行編碼并優(yōu)化。正交表的編碼方案為：對表的每列的所有S!種水平的排列（即水平置換方案），按從小到大的順序依次編碼為0,1…,S!。排列的大小定義如下：一個排列的“值”對應(yīng)于一個s位十進(jìn)制數(shù)，它的第k位恰為在排列中位于第k個的數(shù)，排列的大小以該值的大小表征。依上述編碼方式，m列的正交表可編碼為一個m位S!進(jìn)制數(shù)。本文選擇最大進(jìn)化代數(shù)為1000，種群規(guī)模為水平數(shù)的5倍，交叉因子為1，交叉概率為0.1，依微分進(jìn)化算法[10]尋優(yōu)。尋優(yōu)結(jié)束后開始交互作用列表設(shè)計，當(dāng)讀入新因素時，令各備選列（不作水平置換）與舊因素列搭配，選取均勻性最優(yōu)的備選列為新因素列。

2.3 優(yōu)化效果對比實驗

現(xiàn)在分別使用兩種優(yōu)化算法對L27313,L81340,L64421,L125531進(jìn)行優(yōu)化設(shè)計，要求考慮主因素1與2，1與3的交互作用。表5—8給出了不同主因素總數(shù)下的優(yōu)化結(jié)果的均勻性以及程序耗時（本程序在AMD（主頻2.01G）微機(jī)上運行，運行環(huán)境為32位Microsoft Windows操作系統(tǒng)）

通過對比以上實驗結(jié)果，得出如下結(jié)論：首先，通過對照各表結(jié)果，可發(fā)現(xiàn)當(dāng)主因素只有參與交互的三因素時，原設(shè)計表的均勻性偏差很小，此時優(yōu)化效果微弱(CD2值在小數(shù)點后15位才有所降低)。當(dāng)因素總數(shù)增多時，原表均勻性降低（這印證了設(shè)計的均勻性與分辨力的關(guān)系），而兩種算法的優(yōu)化效果逐漸增強(qiáng)。微分進(jìn)化所得結(jié)果的均勻性更好，計算時間也更多。隨著正交表水平與設(shè)計維數(shù)（生成矩陣P列數(shù)）的增加，微分進(jìn)化與動態(tài)規(guī)劃的耗時比先增加，后減少；其次，對比表3與表4，可發(fā)現(xiàn)當(dāng)正交表水平數(shù)一定，而設(shè)計維數(shù)變大時，動態(tài)規(guī)劃的效果有所下降，而微分進(jìn)化幾乎不受影響；而對比表3與表6，可發(fā)現(xiàn)當(dāng)正交表設(shè)計維數(shù)不變，而水平數(shù)增大時，兩種優(yōu)化算法的效果均明顯增強(qiáng)。表5為一個特例，此時動態(tài)規(guī)劃所得某些結(jié)果的均勻性偏差稍稍高于原設(shè)計表，而微分進(jìn)化所得結(jié)果則優(yōu)于原表。這突出反映了基于動態(tài)規(guī)劃法的優(yōu)化方案無法收斂至全局最優(yōu)，這是因為它在優(yōu)化的每一個階段，只保留了最優(yōu)值，而失掉了有潛力的次優(yōu)值。值得注意的是，在大多數(shù)情況下，動態(tài)規(guī)劃法能快速穩(wěn)定的優(yōu)化原表，而微分進(jìn)化優(yōu)化結(jié)果有隨機(jī)性。以上各表中“微分進(jìn)化”的優(yōu)化數(shù)據(jù)均為多次獨立執(zhí)行優(yōu)化程序所得的最好結(jié)果。

3 結(jié)束語

數(shù)值實驗的結(jié)論表明，本文提出的靈巧設(shè)計，能分辨出工程中的各主效應(yīng)和交互效應(yīng)。大樣本強(qiáng)耦合的正交設(shè)計的均勻性快速全局尋優(yōu)較為困難，本文給出的優(yōu)化方案均實現(xiàn)了不同程度的優(yōu)化。對程序代碼段的測試表明，均勻性指標(biāo)值的計算耗時隨正交表規(guī)模的增大而顯著增加，這是制約優(yōu)化算法快速性的主要因素。

[1]田口玄一著，魏錫祿等譯.實驗設(shè)計法[M].北京：機(jī)械工業(yè)出版社，1987.

[2]方開泰.均勻設(shè)計與均勻設(shè)計表[M].北京：科學(xué)出版社，1994.

[3]C.X.Ma,K.T.Fang,Dennis K.J.Lin.A Note on Uniformity and Orthogonality[J].Journal of Statistical Planning and Inference,2003,（113）.

[4]G..E.P.Box,J.S.Hunter.The 2k-p Fractional Factorial Designs[J].Technometrics,1961,（3）.

[5]C.X.Ma,K.T.Fang.A Note on Generalized Aberration in Factorial Designs[J].Metrika,2001,（53）.

[6]K.T.Fang,R.Mukerjee.A Connection between Uniformity and Aberration in RegularFractions ofTwo-levelFactorials[J].Biometrika,2000,（87）.

[7]Fasheng Sun,Min-Qian Liu,Wenrui Hao.An Algorithmic Approach to Finding Factorial Designs with Generalized Minimum Aberration[J].Journal of Complexity,2009,（25）.

[8]馬希文.正交設(shè)計的數(shù)學(xué)理論[M].北京：人民教育出版社，1981.

[9]方開泰，馬長興.正交與均勻試驗設(shè)計[M].北京：科學(xué)出版社，2001.

[10]Rainer Storn,Kenneth Prince.Differential Evolution-A Simple and Efficient Heuristic for Global Optimization over Continuous Spaces[J].Journal of Global Optimization,1997,（11）.

（責(zé)任編輯/易永生）

TP311.1

A

1002－6487（2010）17-0164-03

韓 冰（1983-），男，江蘇南京人，博士研究生，研究方向：飛行仿真。

唐 碩（1963-），男，四川達(dá)州人，教授，博士生導(dǎo)師，研究方向：飛行仿真。