基于改進的距離函數(shù)法的多響應穩(wěn)健參數(shù)設計

何 楨，王 晶,2 ， 李湧 范

(1. 天津大學管理學院，天津 300072；2. 南開大學泰達學院，天津 300457；3. 中國質(zhì)量認證中心，北京 100070)

響應曲面方法(response surface methodology，RSM)最早由Box和Wilson提出，由于其設計思想先于田口穩(wěn)健性設計方法，所以大多數(shù) RSM 都沒有考慮穩(wěn)健參數(shù)設計問題．Vining和 Myers[1]提出了一種將響應曲面方法與穩(wěn)健性思想結(jié)合起來的雙響應方法，通過對一個質(zhì)量特性的均值和方差進行擬合，根據(jù)擬合的響應曲面模型來優(yōu)化均值并減少方差，從而得到穩(wěn)健的過程參數(shù)．然而，該方法主要解決的是單響應穩(wěn)健參數(shù)設計問題．在現(xiàn)實生產(chǎn)過程中，經(jīng)常會遇到多響應參數(shù)設計問題，如手機、半導體、鋰電池生產(chǎn)、冶金、化工等產(chǎn)品設計和過程設計中，往往存在多個質(zhì)量特性值．多響應問題通常不存在一組特定的輸入變量，可使所有響應變量同時達到最優(yōu)．解決多響應問題的一個簡單有效的方法就是重疊等值線圖法，該方法在一個圖形中為多個響應繪制等值線圖，然后通過觀察找出最優(yōu)解．但這種方法通常只適于輸入變量和響應變量很少的情況．另一類解決多響應問題的思路是運用數(shù)學的方法將多個響應轉(zhuǎn)化為單個響應，然后再進行優(yōu)化，如滿意度函數(shù)法、多元損失函數(shù)法和廣義距離函數(shù)法等．Derringer與Suich[2]提出的滿意度函數(shù)法由于未考慮響應之間的相關(guān)性，該方法難以應對在諸如化工、電子等行業(yè)中經(jīng)常遇到的多響應相關(guān)的情況．Pignatiello[3]、Vining[4]和Ko等[5]提出的質(zhì)量損失函數(shù)法雖然考慮了過程的經(jīng)濟性，但是其成本矩陣的確定過多依賴于定性分析．廣義距離函數(shù)法是由Khuri和Conlon等[6]提出的，使用的是能夠消除方差和交互作用負面影響的馬氏距離，因此也稱為馬氏距離函數(shù)法．該方法使用響應間的方差-協(xié)方差矩陣，考慮了響應間的相關(guān)性，是一種應用廣泛的多響應優(yōu)化方法，但仍未考慮產(chǎn)品或過程的穩(wěn)健性問題．近年來，雖有很多學者對該問題進行研究[7-8]，但是，有關(guān)如何充分利用 RSM 的響應曲面信息，對多響應穩(wěn)健性參數(shù)設計進行的研究還非常有限[9]．

如何解決多響應問題中的穩(wěn)健參數(shù)設計問題，是當今國內(nèi)外質(zhì)量工程領域研究的熱點問題．多響應穩(wěn)健參數(shù)設計方法發(fā)展至今，大多忽略了可控因子的波動對響應的影響．然而在實際生產(chǎn)過程中，可控因子也存在一定變異，當可控因子發(fā)生微小變動時，若響應輸出值的變動較大，則認為該產(chǎn)品或過程不穩(wěn)健．筆者針對目前研究中存在的不足，考慮響應對可控因子波動的穩(wěn)健性，對廣義距離函數(shù)進行了改進，考慮了多響應問題的最優(yōu)性及穩(wěn)健性，提出了穩(wěn)健優(yōu)化總體廣義距離函數(shù)方法，并用該方法對文獻[10]中的實例進行了應用驗證．

1 最優(yōu)點及穩(wěn)健點的確定

試驗區(qū)域中各個響應的最優(yōu)值通常利用Derringer與Suich[2]的滿意度函數(shù)法通過對各個響應單獨優(yōu)化求得．令 Ti表示目標值，Li表示規(guī)格下限，Ui表示規(guī)格上限，則對于望大特性的響應變量 yi，其滿意度函數(shù)為

對于望小特性的響應變量yi，其滿意度函數(shù)為

對于望目特性的響應變量yi，其滿意度函數(shù)為

通過最大化各個響應的滿意度函數(shù)，可以分別得到各響應的最優(yōu)點和最優(yōu)值．

在多響應穩(wěn)健參數(shù)設計中，結(jié)果的最優(yōu)性和穩(wěn)健性同等重要．本文利用穩(wěn)健滿意度函數(shù)法[11]尋找試驗區(qū)域中的穩(wěn)健點．對于試驗區(qū)域中的任意一點，以該點為中心定義一個半徑為 r的圓形區(qū)域，用以描述該點處可控因子的波動范圍．在此區(qū)域中可分別計算各響應的極差

則對于響應變量yi，其穩(wěn)健滿意度函數(shù)為

對于各個響應變量，使其穩(wěn)健滿意度函數(shù)最大的點可認為是試驗區(qū)域中該響應變量的穩(wěn)健點．在各響應的穩(wěn)健點處，響應受可控因子波動的影響最?。?/p>

2 改進的穩(wěn)健優(yōu)化總體廣義距離函數(shù)方法

傳統(tǒng)的廣義距離函數(shù)法使用廣義距離來衡量各響應與其最優(yōu)值之間的偏差，用響應的方差-協(xié)方差結(jié)構(gòu)來估計各響應距其最優(yōu)值距離的權(quán)重，將多響應優(yōu)化問題轉(zhuǎn)化為使廣義距離最小化問題，來確定最佳的操作條件．但是，該方法僅考慮了響應的最優(yōu)性問題．考慮到可控因子的波動對響應穩(wěn)健性的影響，利用穩(wěn)健滿意度函數(shù)法計算出各響應的穩(wěn)健解，定義一個穩(wěn)健廣義距離函數(shù)來衡量各個響應與其穩(wěn)健值之間的偏差，然后使其最小化，從而使得多個響應同時達到或接近其穩(wěn)健值，考慮了響應的穩(wěn)健性要求．然而從統(tǒng)計分析及質(zhì)量管理的角度出發(fā)，響應的最優(yōu)性以及穩(wěn)健性同等重要．穩(wěn)健優(yōu)化總體廣義距離函數(shù)法通過對優(yōu)化廣義距離和穩(wěn)健廣義距離進行權(quán)衡，求二者的幾何均值，然后使其最小化，得到穩(wěn)健最優(yōu)的操作條件，使優(yōu)化的響應對可控因子的波動具有穩(wěn)健性．

假設某一生產(chǎn)過程中所有m個響應都依賴于同一組 k個輸入變量{x1，x2，…，xk}，共計 n次試驗，且均可以用同階次的多項式回歸模型表示．則第i個響應的二階模型可以用向量形式表示為

式中：yi為一個n×1響應矩陣，表示第i個響應的n個觀測值；X為n×p設計矩陣，p =(k + 1)(k +2)/2；βi為表示回歸系數(shù)的p×1列向量；εi為一個n×1列向量，表示第 i個響應的隨機誤差，服從均值為 0、方差為σii的正態(tài)分布．進一步假設第i個響應的隨機誤差與第 j個響應的隨機誤差的協(xié)方差cov(εi,εj)=σijIn，i ≠ j = 1,2,… ,m，其中In為n×n單位矩陣．βi的最小二乘估計值可推導為=(X ′X )?1，則第 i個響應的二階擬合模型為

var[y?(x ) ] = z ′(x ) (X ′X )?1z ( x)σ i = 1,2,…,m

i ii y?i(x)和y?j(x)的協(xié)方差為

其中ijσ為方差-協(xié)方差矩陣Σ的第(i，j)個元素．Σ的無偏估計可表示為

用?Σ替代式(4)中的Σ可以得到? var[()]yx 的無偏估計值

用 θi表示響應yi在試驗區(qū)域內(nèi)單獨優(yōu)化時得到的最優(yōu)值，θ ＝ [θ1,θ2,… ,θm]′ 表示由這些最優(yōu)值構(gòu)造的m×1向量，則衡量各個響應與其最優(yōu)值之間偏差的優(yōu)化廣義距離函數(shù)為

因為式(6)中 v ar[y?(x ) ]是正定矩陣，故而 { var[y?(x)]}?1也是正定矩陣，因此當y?(x)≠ θ 時，距離ρ[y?(x),θ]嚴格為正．將響應的方差-協(xié)方差矩陣的估計形式(5)代入式(6)中，可得

通過最小化廣義距離 ρ?[y?(x),θ]可確定最佳的操作條件，使各響應同時達到或接近最優(yōu)值．

用φi表示響應yi在試驗區(qū)域內(nèi)得到的各自穩(wěn)健點的響應值，φ ＝ [φ1,φ2,… ,φm]′ 表示由這些穩(wěn)健值構(gòu)造的 m×1向量，則衡量各個響應與其穩(wěn)健值之間偏差的穩(wěn)健廣義距離函數(shù)為

將式(5)代入上式，可得其估計值為

通過最小化 ρ?′[y?(x),φ]可確定最穩(wěn)健的操作條件．

考慮多響應的最優(yōu)性和穩(wěn)健性的權(quán)衡問題，穩(wěn)健優(yōu)化總體廣義距離表達式為

式中 w1和 w2分別為優(yōu)化廣義距離和穩(wěn)健廣義距離的權(quán)重，均為正實數(shù)且滿足 w1+ w2= 1．使穩(wěn)健優(yōu)化總體廣義距離最小的解可認為是過程的穩(wěn)健優(yōu)化解．根據(jù)多響應問題的最優(yōu)性及穩(wěn)健性的重要程度，權(quán)重w1和w2可選擇不同的取值．當w1＝1、w2=0時，ρoverall=ρ[y?(x),θ]，只考慮多響應的最優(yōu)性問題；當w1＝0、w2=1 時，ρoverall= ρ ′[y?(x),φ]，只考慮多響應的穩(wěn)健性問題；當響應的最優(yōu)性及穩(wěn)健性具有相同的權(quán)重，即 w1＝w2=0.5時，ρoverall= ( ρ [ y? (x ) ,θ ] )12(ρ′[y?(x ),φ] )12，在實際生產(chǎn)中這種情況應用最為廣泛．

3 實例分析

本文以 Myers與 Montgomery[10]介紹過的一個經(jīng)典的多響應試驗設計問題為例，利用本文提出的改進的穩(wěn)健優(yōu)化總體廣義距離函數(shù)方法對其進行計算，并與傳統(tǒng)方法進行比較分析．該試驗通過調(diào)整一個化學過程的反應時間(x1)和反應溫度(x2)，從而使其產(chǎn)量(y1)最大化，黏性(y2)達到目標值 65，同時使相對分子質(zhì)量(y3)最小化，假定各響應變量可接受的范圍分別是 y1≥78.5，62≤y2≤68，y3≤3,300．該試驗是一個含有 5個中心點的中心復合設計，試驗數(shù)據(jù)見表 1．

首先分別擬合出各響應的二階響應曲面模型為

表1 試驗數(shù)據(jù)Tab.1 Test data

方差-協(xié)方差矩陣的無偏估計值可通過式(3)計算得到，即

因為該試驗設計是一個球形中心復合設計，因此試驗區(qū)域是一個以原點為圓心、半徑為 2的圓，根據(jù)擬合得到的二階響應曲面模型，對各響應單獨優(yōu)化得到的最優(yōu)值θi以及取得最優(yōu)值的可控因子的參數(shù)組合如表 2所示．各響應在試驗區(qū)域內(nèi)得到的穩(wěn)健值iφ及其對應的穩(wěn)健點坐標見表3．

表2 各響應的最優(yōu)值θi與其對應點坐標Tab.2 Optimal value θi and its corresponding location for each response

表3 各響應的穩(wěn)健值iφ及其對應點坐標Tab.3 Robust value iφ and its corresponding location for each response

首先利用Minitab軟件中多響應優(yōu)化的滿意度函數(shù)方法解決該問題．選擇T1＝80作為響應產(chǎn)量y1的目標值，L1=78.5作為規(guī)格下限，為望大特性．令 T2＝65作為響應黏性 y2的目標值，L2=62作為規(guī)格下限，U2＝68作為規(guī)格上限，為望目特性．選擇響應相對分子質(zhì)量y3的目標值T3=3.100，規(guī)格上限為U3＝3.300，為望小特性．各響應單獨滿意度函數(shù)中所有權(quán)重均取 1．若考慮對3個響應同時優(yōu)化，Minitab給出一組可控變量的最佳參數(shù)組合：x1=-0.157.1，x2=-0.785,7．在該水平下，產(chǎn)量的預測值為78.76，黏性的預測值為 66.47，相對分子質(zhì)量的預測值為3.229.5．但是滿意度函數(shù)法沒有考慮最優(yōu)點處是否穩(wěn)?。谇蟮玫淖顑?yōu)點(－0.157.1，－0.785.7)處，定義可控因子的變動幅度為 0.01，則在此范圍內(nèi)各響應的極差分別為：R1＝0.065.1，R2＝0.215.2，R3＝7.526.2．

然后利用廣義距離函數(shù)法對該問題進行分析．前面已經(jīng)求得方差-協(xié)方差矩陣的無偏估計值和各響應單獨優(yōu)化時得到的最優(yōu)值θi．利用廣義距離函數(shù)法求得一組可控因子的最優(yōu)參數(shù)組合：x1＝0.048,2，x2=-0.626,2．在該水平下，產(chǎn)量的預測值為79.26，黏性的預測值為 68，相對分子質(zhì)量的預測值為3,300．同樣，若考慮響應的穩(wěn)健性，在求得的最優(yōu)點(0.048,2，-0.626,2)處，定義可控因子的變動幅度為0.01，在該區(qū)域內(nèi)各個響應的極差分別為：R1＝0.049,7，R2＝0.158,5，R3＝7.034,3．與滿意度函數(shù)法相比，利用廣義距離函數(shù)法計算得到的最優(yōu)點處各響應的極差均較小，因此可以認為各響應在(0.048,2，-0.626,2)處，較之在(-0.157,1，-0.785,7)處更穩(wěn)?。?/p>

在實際生產(chǎn)過程中，找到既穩(wěn)健又可使響應達到最優(yōu)的參數(shù)組合十分重要．下面利用本文提出的穩(wěn)健優(yōu)化總體廣義距離函數(shù)法求解．當 w1=1、w2＝0時，得到的結(jié)果與利用傳統(tǒng)廣義距離函數(shù)法的結(jié)果相同．當 w1＝0、w2＝1時，求得的可控因子的穩(wěn)健參數(shù)組合為 x1＝0.058,3，x2＝-0.641,6．在該水平下，產(chǎn)量的預測值為 79.24，黏性的預測值為 67.89，相對分子質(zhì)量的預測值為 3,301，各響應的極差值為 R1＝0.048,9，R2＝0.161,0，R3＝6.886,3．雖然在該參數(shù)組合下響應的極差值最小，但是由于相對分子質(zhì)量的預測值超出其規(guī)格上限，因此是不可接受的．通常情況下，需要同時考慮響應的最優(yōu)性和穩(wěn)健性問題，因此賦予二者相同的權(quán)重 w1＝w2＝0.5，利用該方法求得一組可控因子的穩(wěn)健最優(yōu)參數(shù)組合：x1＝0.051,2，x2＝-0.649．在該水平下，產(chǎn)量的預測值為79.22，黏性的預測值為 67.83，相對分子質(zhì)量的預測值為3,298.5．當可控因子的變動幅度為0.01時，各響應的極差值為 R1＝0.049,5，R2＝0.163,5，R3＝6.904,0．比較利用3種方法計算得到的結(jié)果，如表4所示．

表4 改進的廣義距離函數(shù)方法與傳統(tǒng)方法的結(jié)果對比Tab.4 Comparison of results between improved generalized distance function method and traditional methods

若僅比較預測的響應值，3種方法中利用滿意度函數(shù)法可以得到最好的結(jié)果，但是預測的響應極差較大，說明最優(yōu)點處的穩(wěn)健性較差．雖然另 2種方法得到的響應的預測值稍差，但是響應的穩(wěn)健性都要好于滿意度函數(shù)法計算的結(jié)果．另一方面，利用廣義距離函數(shù)法預測的響應 y2和 y3都在響應的規(guī)格限上，可控因子的微小波動就可能會導致響應變量超出規(guī)格限，造成大量廢品的產(chǎn)生．利用穩(wěn)健優(yōu)化廣義距離函數(shù)進行分析：當 w1=1、w2=0時，可簡化為傳統(tǒng)廣義距離函數(shù)法；當w1=w2=0.5時，預測的響應y1和y3的極差都要小于傳統(tǒng)方法，表明點(0.051.2，－0.649)處的響應y1和y3更穩(wěn)健．雖然響應y2的極差R2相對廣義距離函數(shù)法更大，但是預測的響應值 67.83與廣義距離函數(shù)法相比更接近目標值，不再是響應的規(guī)格限，可大大減少廢品的產(chǎn)生．這說明利用改進的距離函數(shù)法對多響應問題進行優(yōu)化時，權(quán)衡了響應的穩(wěn)健性和最優(yōu)性．隨著穩(wěn)健性的權(quán)重w2變大，預測的響應值的極差相應變小，但是同時響應值也逐漸遠離目標值，當 w1=0、w2=1時，響應 y3的預測值 3,301已超出其規(guī)格限 3,300，不能滿足顧客要求．因此在實際生產(chǎn)中通常選取 w1=w2=0.5，即式(10)取幾何均值作為穩(wěn)健優(yōu)化廣義距離函數(shù)．該方法與傳統(tǒng)方法相比，考慮了響應的方差-協(xié)方差結(jié)構(gòu)和可控因子的波動，使優(yōu)化的響應對可控因子的波動具有穩(wěn)健性．

4 結(jié) 語

多響應穩(wěn)健參數(shù)設計是改進產(chǎn)品或過程質(zhì)量的一個重要方法．本文基于廣義距離函數(shù)法對多響應穩(wěn)健參數(shù)設計進行了研究，考慮到響應的方差-協(xié)方差結(jié)構(gòu)和可控因子的波動，提出了改進的穩(wěn)健優(yōu)化總體廣義距離函數(shù)法，通過對優(yōu)化廣義距離和穩(wěn)健廣義距離進行權(quán)衡，將多響應穩(wěn)健優(yōu)化問題轉(zhuǎn)化為使穩(wěn)健優(yōu)化總體廣義距離最小化問題，從而得到穩(wěn)健最優(yōu)的操作條件，使優(yōu)化的響應對可控因子的波動具有穩(wěn)健性．

［1］ Vining G G，Myers R H. Combining Taguchi and response surface philosophies：A dual response approach[J]. Journal of Quality Technology，1990，22(1)：38-45.

［2］ Derringer G C，Suich R. Simultaneous optimization of several response variables[J]. Journal of Quality Technology，1980，12(4)：214-219.

［3］ Pignatiello J J. Strategies for robust multiresponse quality engineering [J]. IIE Transactions，1993，25(3)：5-15.

［4］ Vining G G. A compromise approach to multiresponse optimization[J]. Journal of Quality Technology，1998，30(4)：309-313.

［5］ Ko Y H，Kim K J，Jun C H. A new loss function-based method for multiresponse optimization[J]. Journal of Quality Technology，2005，37(1)：50-59.

［6］ Khuri A I，Conlon M. Simultaneous optimization of multiple responses represented by polynomial regression functions[J]. Technometrics，1981，23(4)：363-375.

［7］ Park K S，Kim K J. Optimizing multi-response surface problems：How to use multi-objective optimization techniques[J]. IIE Transactions，2005，37(6)：523-532.

［8］ Onur K?ksoy. A nonlinear programming solution to robust multi-response quality problem[J]. Applied Mathematics and Computation，2008，196(2)：603-612.

［9］ Beyer H G，Sendhoff B. Robust optimization—A comprehensive survey[J]. Computer Methods in Applied Mechanics and Engineering，2007，196(33/34)：3190-3218.［10］ Myers R H，Montgomery D C. Response Surface Methodology[M]. 2nd ed. New York：John Wiley & Sons，2002.

［11］ Wang Jing，He Zhen，Oh J H，et al. Multi-response robust optimization using desirability function[C] // IEEE Symposium on Advanced Management of Information for Globalized Enterprises. Tianjin，China，2008：313-315.