基于FDS準(zhǔn)則的響應(yīng)曲面設(shè)計(jì)預(yù)測(cè)能力比較

周 祥，陳雪平，嚴(yán) 靜，林金官

（1.蘇州高等職業(yè)技術(shù)學(xué)校（江蘇聯(lián)合職業(yè)技術(shù)學(xué)院蘇州分院）教務(wù)處，江蘇 蘇州 215000；2.江蘇理工學(xué)院 數(shù)理學(xué)院，江蘇 常州 213001；3.南京審計(jì)大學(xué) 統(tǒng)計(jì)與數(shù)據(jù)科學(xué)學(xué)院，南京 211815）

0 引言

試驗(yàn)設(shè)計(jì)方法（Design of Experiment）已被廣泛應(yīng)用在工學(xué)、農(nóng)學(xué)、醫(yī)學(xué)、管理、經(jīng)濟(jì)等各個(gè)領(lǐng)域。一個(gè)好的設(shè)計(jì)方案可以提高產(chǎn)量、降低成本、增加效益。隨著生產(chǎn)發(fā)展和應(yīng)用的需要，諸多試驗(yàn)設(shè)計(jì)方法不斷被人們發(fā)現(xiàn)和使用。如農(nóng)業(yè)時(shí)期的區(qū)組設(shè)計(jì)、正交設(shè)計(jì)等，工業(yè)時(shí)期的均勻設(shè)計(jì)、混料設(shè)計(jì)、最優(yōu)設(shè)計(jì)、響應(yīng)曲面設(shè)計(jì)（Response Surface Design）[1]等。其中，響應(yīng)曲面設(shè)計(jì)常常被用來探索響應(yīng)變量與多個(gè)預(yù)測(cè)變量之間的關(guān)系，被廣泛應(yīng)用于各類工程實(shí)踐、質(zhì)量控制領(lǐng)域中。響應(yīng)曲面設(shè)計(jì)的分析包含了設(shè)計(jì)、建模、數(shù)據(jù)分析、最優(yōu)化等過程。

隨著研究的不斷深入，學(xué)者們陸續(xù)提出了各類響應(yīng)曲面設(shè)計(jì)方案，比如CCD 設(shè)計(jì)、BBD設(shè)計(jì)、SCD設(shè)計(jì)[1]、APD設(shè)計(jì)[2]等。在研究一類藥物組合試驗(yàn)設(shè)計(jì)時(shí)，Xu等（2014）[3]提出了一類基于正交表的響應(yīng)曲面設(shè)計(jì)——OACD（Orthogonal-array Composite Design），該設(shè)計(jì)比前述響應(yīng)曲面設(shè)計(jì)方法的D效率更高；Chen等（2018）[4]在此基礎(chǔ)上提出了基于極小極大準(zhǔn)則的改進(jìn)設(shè)計(jì)OACM（Orthogonal-array Composite Minmax Design），該設(shè)計(jì)具有在隨機(jī)缺失一個(gè)響應(yīng)數(shù)據(jù)時(shí)，原始設(shè)計(jì)的最大損失達(dá)到最小的特性。

本文重點(diǎn)關(guān)注從預(yù)測(cè)角度出發(fā)的優(yōu)良性準(zhǔn)則，即對(duì)于一類設(shè)計(jì)，希望其在設(shè)計(jì)點(diǎn)上的預(yù)測(cè)方差較小，從而使預(yù)測(cè)具有一定的穩(wěn)健性，如D-最優(yōu)、A-最優(yōu)、E-最優(yōu)、G-最優(yōu)等最優(yōu)設(shè)計(jì)[5—7]。從預(yù)測(cè)方差v(x)出發(fā)，Giovannitti-Jensen 和Myers（1989）[8]提出用一種圖像評(píng)價(jià)方法，即方差散布圖VDG（Variance Dispersion Graph）來比較試驗(yàn)區(qū)域上不同點(diǎn)的預(yù)測(cè)方差性質(zhì)。隨后，不同學(xué)者對(duì)不同設(shè)計(jì)方案的預(yù)測(cè)方差性質(zhì)進(jìn)行了研究。Vining 和Myers（1991）[9]采用均方誤差圖來評(píng)價(jià)響應(yīng)曲面設(shè)計(jì)；Borkowski（1995）[10]研究了普通中心組合設(shè)計(jì)和Box-Behnken 設(shè)計(jì)的預(yù)測(cè)方差性質(zhì)。然而，VDG 準(zhǔn)則僅討論了在給定球面上設(shè)計(jì)點(diǎn)的預(yù)測(cè)方差性質(zhì)，無法綜合考慮在整個(gè)區(qū)域上預(yù)測(cè)方差的大小。為解決這一問題，Zahran等（2003）[11]提出了全區(qū)域預(yù)測(cè)方差占比準(zhǔn)則——FDS 準(zhǔn)則（Fraction of Design Space Criterion），與VDG準(zhǔn)則不同，F(xiàn)DS 準(zhǔn)則是在給定可容許預(yù)測(cè)方差值的情況下，比較預(yù)測(cè)方差小于給定值的區(qū)域占比，因此，F(xiàn)DS值越大，表明小于給定可容許預(yù)測(cè)偏差的設(shè)計(jì)區(qū)域越大，設(shè)計(jì)就越優(yōu)。

本文應(yīng)用FDS準(zhǔn)則，首先，介紹了OACD 和OACM，并比較他們的預(yù)測(cè)方差性質(zhì)；然后，應(yīng)用三種預(yù)測(cè)方差，包括在給定距離半徑上點(diǎn)的最大預(yù)測(cè)方差、最小預(yù)測(cè)方差和平均預(yù)測(cè)方差，采用好格子點(diǎn)法比較了OACD和OACM在不同因子數(shù)下的FDS值；最后，討論了OACD以及OACM的優(yōu)良性。

1 基于正交表的響應(yīng)曲面設(shè)計(jì)

通常，響應(yīng)曲面設(shè)計(jì)用一個(gè)二次模型[1—4]來描述響應(yīng)變量y與x1,…,xk之間的關(guān)系：

其中，β0、βj、βjj和βj1j2分別為截距項(xiàng)及線性項(xiàng)、二次項(xiàng)和交叉乘積項(xiàng)的系數(shù)，?i一般假定為獨(dú)立同分布于均值為0、方差為σ2的正態(tài)分布隨機(jī)誤差。在模型（1）下，對(duì)于設(shè)計(jì)點(diǎn)x，記：

相應(yīng)的預(yù)測(cè)方差函數(shù)為：

從信息矩陣出發(fā)，可以構(gòu)造各類最優(yōu)設(shè)計(jì)，常見的響應(yīng)曲面設(shè)計(jì)主要包括以中心組合設(shè)計(jì)為基礎(chǔ)的各種變換和改進(jìn)。下面簡(jiǎn)要介紹中心組合設(shè)計(jì)的一般構(gòu)造。

定義1：一般地，一個(gè)含有k個(gè)試驗(yàn)因子(x1,…,xk)的中心組合設(shè)計(jì)由三部分試驗(yàn)點(diǎn)組成：

（1）將編碼值-1 與1 看成每個(gè)因子的兩個(gè)水平，采用二水平正交表安排試驗(yàn)，可以是全因子試驗(yàn)，也可以是1/2實(shí)施、1/4實(shí)施等，記其試驗(yàn)次數(shù)為n1。

（2）星號(hào)點(diǎn)：依次取試驗(yàn)點(diǎn)的某一個(gè)因子水平為-α、α，其余因子水平均為0，記其試驗(yàn)次數(shù)為n2。

（3）對(duì)于中心點(diǎn)(0,…,0)，記其試驗(yàn)次數(shù)為n0。

中心組合設(shè)計(jì)具有以下特點(diǎn)：該方案的總試驗(yàn)次數(shù)為n=n1+n2+n0，二水平正交表的分辨度可以達(dá)到V，星號(hào)點(diǎn)的試驗(yàn)次數(shù)n2=2k。于是每個(gè)因子都可以取5個(gè)水平，所以該方案分布的試驗(yàn)點(diǎn)范圍比較廣，有很大的靈活性[1,2]。而且在方案中留有兩個(gè)待定參數(shù)（中心點(diǎn)的試驗(yàn)次數(shù)n0和星號(hào)點(diǎn)的水平值α），以便使其滿足如設(shè)計(jì)正交性或者旋轉(zhuǎn)性的要求。

為了進(jìn)一步改進(jìn)設(shè)計(jì)的效率，Xu等（2014）[3]提出了基于正交表的響應(yīng)曲面設(shè)計(jì)（OACD）。研究表明，OACD 在設(shè)計(jì)D效率上優(yōu)于常見的中心組合設(shè)計(jì)。

定義2：一般地，一個(gè)含有k個(gè)試驗(yàn)因子(x1,…,xk)的OACD由三部分試驗(yàn)點(diǎn)組成：

（1）將編碼值-1 與1 看成每個(gè)因子的兩個(gè)水平，采用二水平正交表安排試驗(yàn)，可以是全因子試驗(yàn)，也可以是部分實(shí)施，記其試驗(yàn)次數(shù)為n1。

（2）將編碼值-1、0、1看成每個(gè)因子的三個(gè)水平，采用三水平正交表安排試驗(yàn)，可以是全因子試驗(yàn)，也可以是部分實(shí)施，記其試驗(yàn)次數(shù)為n2。

（3）對(duì)于中心點(diǎn)(0,…,0)，記其試驗(yàn)次數(shù)為n0。

于是，基于正交表的中心組合設(shè)計(jì)與普通中心組合設(shè)計(jì)相比，區(qū)別在于第二部分設(shè)計(jì)點(diǎn)的選擇。表1為一個(gè)34行的5 個(gè)因子的OACD，其中1 至16 行為二水平正交設(shè)計(jì)，其生成字為E=ABCD，17 至34 行為三水平正交設(shè)計(jì)OA（18），中心點(diǎn)處的試驗(yàn)次數(shù)為0。表2列出了不同因子數(shù)下相應(yīng)OACD的具體構(gòu)造[4]。本文二水平正交設(shè)計(jì)的最優(yōu)選擇主要依賴于最大分辨度準(zhǔn)則，而三水平非正規(guī)正交設(shè)計(jì)的最優(yōu)選擇主要是最小低階混雜準(zhǔn)則。可以發(fā)現(xiàn)，表2中因子數(shù)k=5 的OACD即表1。

表1 一個(gè)34行的OACD（α=1,n0=0）

表2 本文所采用的OACD（k=3,…,8）

在實(shí)際的試驗(yàn)過程中，數(shù)據(jù)缺失往往很難避免，特別是在采用飽和設(shè)計(jì)，或者超飽和設(shè)計(jì)的情形下，缺失數(shù)據(jù)對(duì)模型分析、參數(shù)估計(jì)可能有較大的影響。于是，人們總希望能找到缺失數(shù)據(jù)對(duì)效率損失影響較小的設(shè)計(jì)。很多學(xué)者研究了在面對(duì)缺失數(shù)據(jù)時(shí)，如何保持試驗(yàn)設(shè)計(jì)的有效性，如Akhtar 和Prescott（1986）[12]討論了中心組合設(shè)計(jì)在缺失試驗(yàn)數(shù)據(jù)時(shí)的最小損失設(shè)計(jì)；Ahmad 和Gilmour（2010）[13]討論了擴(kuò)展配對(duì)設(shè)計(jì)（Augmented Pairs Design）在含有缺失數(shù)據(jù)情形時(shí)的最優(yōu)設(shè)計(jì)。針對(duì)OACD，Chen等（2018）[4]基于極小極大準(zhǔn)則，通過TA門限搜索算法，給出了在帶有缺失數(shù)據(jù)時(shí)，最大損失達(dá)到最小的設(shè)計(jì)方案，該設(shè)計(jì)通過確定參數(shù)α來達(dá)到最小損失，具體參數(shù)α的選擇見表3。

表3 OACM的相對(duì)D效率（相比于OACD）

定義3：一般地，一個(gè)含有k個(gè)試驗(yàn)因子x1,…,xk的OACM由三部分試驗(yàn)點(diǎn)組成：

（1）將編碼值-1 與1 看成每個(gè)因子的兩個(gè)水平，采用二水平正交表安排試驗(yàn)，可以是全因子試驗(yàn)，也可以是部分實(shí)施，記其試驗(yàn)次數(shù)為n1。

（2）將編碼值 -α、0、α看成每個(gè)因子的三個(gè)水平，采用三水平正交表安排試驗(yàn)，可以是全因子試驗(yàn)，也可以是部分實(shí)施。α的選取是基于極大極小準(zhǔn)則得到的，其中，li(D)為設(shè)計(jì)缺失第i行后的效率損失，記其試驗(yàn)次數(shù)為n2。

（3）對(duì)于中心點(diǎn)(0,…,0)，記其試驗(yàn)次數(shù)為n0。

比如，對(duì)于表2的OACD，其因子數(shù)k=5，當(dāng)中心點(diǎn)處的試驗(yàn)次數(shù)為1 時(shí)，對(duì)應(yīng)OACM 中參數(shù)的選取可以在表3中查得，α=1.1775 。該設(shè)計(jì)的D 效率優(yōu)于相應(yīng)的OACD ，相對(duì)D 效率為1.2851，即該OACM 是相應(yīng)OACD D 效率的1.2851 倍。當(dāng)k=5,n0=2 時(shí)，OACM對(duì)應(yīng)的參數(shù)α=1.1799，且OACM相對(duì)于OACD的D效率為1.2946。從表3中可以看出，除k=4,n0=1的情形外，OACM的D效率均優(yōu)于相應(yīng)的OACD。本文將進(jìn)一步比較這兩類設(shè)計(jì)在預(yù)測(cè)方差性質(zhì)上的差異。

2 預(yù)測(cè)方差模擬比較

本文通過方差散布圖VDG[8,10]來比較不同試驗(yàn)設(shè)計(jì)方法在同一球面上預(yù)測(cè)方差的差異。記設(shè)計(jì)包含的因子數(shù)為k，離球心距離為r的球面點(diǎn)集為，相應(yīng)設(shè)計(jì)點(diǎn)xi=(xi1,xi2,…,xip)的預(yù)測(cè)方差函數(shù)為即式（2）。記xi表示試驗(yàn)因子，于是，以下三個(gè)預(yù)測(cè)方差分布情況[8—11]的指標(biāo)經(jīng)常用來表示三種VDG。

定義4：平均預(yù)測(cè)方差Vavg、最大預(yù)測(cè)方差Vmax和最小預(yù)測(cè)方差Vmin為：

其中，R=dx表示區(qū)域Sr的面積。

三種VDG都在一定程度上反映了設(shè)計(jì)的預(yù)測(cè)方差性質(zhì)，Vavg反映了在給定球面上所有預(yù)測(cè)方差的平均大小，Vmax反映了在給定球面上所有預(yù)測(cè)方差的最大值，Vmin反映了在給定球面上所有預(yù)測(cè)方差的最小值。

平均值受極值的影響較大，有時(shí)會(huì)導(dǎo)致采用平均值進(jìn)行比較時(shí)不夠穩(wěn)定，而用最大值或者最小值進(jìn)行比較時(shí)，并不清楚最大值或者最小值特征區(qū)域的比例值。為進(jìn)一步解決這些問題，比如，實(shí)際實(shí)驗(yàn)者想要了解小于特定預(yù)測(cè)方差值的設(shè)計(jì)區(qū)域范圍，于是給出了區(qū)域占比指標(biāo)準(zhǔn)則——FDS準(zhǔn)則（Fraction of Design Space Criterion）。

定義5：對(duì)于一個(gè)給定的預(yù)測(cè)方差值v，預(yù)測(cè)方差小于v的試驗(yàn)區(qū)域占比為：

其中，k為設(shè)計(jì)的因子數(shù)，? 為整個(gè)設(shè)計(jì)區(qū)域的體積，Φ={(x1,…,x2):v(x)

顯然，0

其中，g=v(x)-v，于是式（3）等價(jià)于：

其中，D表示整個(gè)設(shè)計(jì)區(qū)域。

在模擬計(jì)算中，可采用好格子點(diǎn)法進(jìn)行近似，即在整個(gè)設(shè)計(jì)區(qū)域均勻抽樣N個(gè)點(diǎn)，記滿足條件集合Φ中的點(diǎn)個(gè)數(shù)為nA，于是FDS近似為nA/N。

例1：由于k=2 是二維平面的形式，可以清晰展示預(yù)測(cè)方差的性質(zhì)，因此本文以k=2 時(shí)的正交表復(fù)合設(shè)計(jì)OACD為例進(jìn)行說明。假如想確定預(yù)測(cè)方差小于5.5的設(shè)計(jì)區(qū)域，從Vmax和Vmin可以知道，在整個(gè)設(shè)計(jì)區(qū)域里，預(yù)測(cè)方差滿足3.2

圖1 OACD在k=2 時(shí)的預(yù)測(cè)方差圖

對(duì)于任意的兩個(gè)設(shè)計(jì)，通常希望FDS值在給定v時(shí)越大越好，下面給出FDS值的模擬算法：

第1步：給定試驗(yàn)因子數(shù)k，根據(jù)表2和表3確定相應(yīng)的OACD和OACM。

第2 步：給定預(yù)測(cè)方差值v，在試驗(yàn)區(qū)域內(nèi)隨機(jī)生成N個(gè)試驗(yàn)點(diǎn)，并計(jì)算他們的預(yù)測(cè)方差v(x)。

第3步：計(jì)算預(yù)測(cè)方差小于v的試驗(yàn)點(diǎn)個(gè)數(shù)，記為nA，于是FDS可以近似為。

采用上述模擬算法，可以給出所有情況下的FDS圖。下頁(yè)圖2分別給出了OACD和OACM在不同因子數(shù)k下的FDS預(yù)測(cè)方差比較結(jié)果，其中，因子數(shù)k=3,4,…,8，相應(yīng)設(shè)計(jì)區(qū)域的球面半徑范圍為，設(shè)計(jì)的構(gòu)造取自表2。比如，當(dāng)設(shè)計(jì)因子數(shù)k=3 時(shí)，在相應(yīng)的OACD中，二水平部分為完全設(shè)計(jì)23，此時(shí)，n1=8；三水平部分為正交表OA（9）的1—3 列，此時(shí)，n2=9。根據(jù)表3，在k=3,n0=1時(shí)，基于極小極大準(zhǔn)則得到的OACM的參數(shù)為α=1.0486，設(shè)計(jì)D 效率相比原始OACD 提高了4.97%，其余設(shè)計(jì)依次類推。

圖2 OACD和OACM在不同因子數(shù)k 下的FDS 預(yù)測(cè)方差比較圖

模擬結(jié)論：從圖2中可以清晰地看出，OACD與OACM的預(yù)測(cè)方差比較結(jié)果可以分為兩個(gè)部分：一方面，在預(yù)測(cè)方差較小的某一范圍里(v

表4 OACD設(shè)計(jì)與OACM設(shè)計(jì)的臨界值

例2：以k=5 為例，由表4可知，在預(yù)測(cè)方差較小的某一范圍里，即v

表5 OACM與OACD在k=5 時(shí)的FDS 值比較

比如，當(dāng)v=11 時(shí)，v

當(dāng)k=4 時(shí)，OACD和OACM的差異本身較小，因此，兩者的FDS值差異也較小。

4 結(jié)論

預(yù)測(cè)能力是評(píng)價(jià)設(shè)計(jì)優(yōu)劣的一個(gè)重要準(zhǔn)則，本文基于預(yù)測(cè)方差占比指標(biāo)FDS，討論了基于正交表的響應(yīng)曲面設(shè)計(jì)（OACD）以及基于極小極大準(zhǔn)則的改進(jìn)設(shè)計(jì)（OACM）的性能。模擬結(jié)果表明，在不同預(yù)測(cè)方差范圍內(nèi)，兩類設(shè)計(jì)各有各的優(yōu)點(diǎn)：在預(yù)測(cè)方差較小的某一范圍（v