試驗設計與數據處理第六章續(第七,八節)

第七節二次響應面分析法當試驗中考察的指標宜于用多元二次回歸方程來擬合因素與指標的函數關系，就可以分析回歸方程所反應的曲面形狀，如果得到的曲面是凸面(像山丘)或凹面(像山谷)這類簡單曲面，那么預測的最佳指標值(極大值或極小值)可以從所估計的曲面上獲得；如果曲面很復雜，或者預測的最佳點遠離所考察因素的試驗范圍，那么可以通過嶺嵴分析來確定重新進行試驗的方向.這就是應用較廣，頗有實用價值的響應面分析法(ResponseSurfaceAnalysis).響應面分析的三個任務:任務1.模型擬合和方差分析;任務2.典型分析(尋找最優值點);任務3.嶺嵴分析(若最優值點在試驗范圍內達到,則不必進行嶺嵴分析;否則就要進行).例6.7.2求面包酵母高產和高發酵活力的最優工藝參數.指標:y1---面包酵母產量[指最終發酵液的菌體濃度(g/L)]y2---發酵活力[指面團從水中浮起的時間(min)]因素:z1---還原糖濃度,2%~6%;z2---氮源濃度,0.4%~1.2%;z3---磷源濃度,0.2%~0.4%.取三水平,z1:2%,4%,6%;z2:0.4%,0.8%,1.2%;z3:0.2%,0.3%,0.4%.令編碼變量xi=[zi-(zi1+zi2)/2]/(zi2-zi1)/2故可作變換:X1=(z1-4)/2X2=(z2-0.8)/0.4X3=(z3-0.3)/0.1響應面分析的SAS程序如下:DataE62;Inputz1-z3x1-x3y1y2@@;Cards;

40.40.20-1-124.513.9

40.40.40-1135.314.1………40.80.300037.58.1

;Proc

RsReg

data=E62;/*響應面分析*/Modely1y2=z1-z3;Run;

TheRSREGProcedureCodingCoefficientsfortheIndependentVariablesFactorSubtractedoffDividedbyz14.0000002.000000z20.8000000.400000z30.3000000.100000ResponseSurfaceforVariabley1ResponseMean32.400000RootMSE1.329474R-Square0.9775CoefficientofVariation4.1033下面給出程序運行得到的部分結果:

TypeISumRegressionDFofSquaresR-SquareFValuePr>FLinear314.5575000.03712.750.1524Quadratic3116.5325000.297221.980.0026Crossproduct3252.1725000.643147.560.0004TotalModel9383.2625000.977524.090.0013SumofResidualDFSquaresMeanSquareTotalError58.8375001.767500

TheRSREGProcedureCanonicalAnalysisofResponseSurfaceBasedonCodedDataCriticalValueFactorCodedUncodedz1-0.0988443.802311z20.2648400.905936z3-0.3158000.268420Predictedvalueatstationarypoint:38.040882EigenvectorsEigenvaluesz1z2z3-0.5185070.184252-0.7614870.621441-1.7733340.894929-0.131446-0.426408-8.0206590.4063900.6347120.657258Stationarypointisamaximum.如果穩定點不是理想點就要進一步作嶺嵴分析,請看示意圖和例子演示(例見datarubber或datapectic)例6.7.3橡膠樹幼苗對硫35的吸收實驗(2×3析因試驗)因素:溫度t(℃)---7℃,20℃

取樣時間d---2d,7d,12d響應面分析SAS程序如下:datarubber;inputftd@@;cards;數(略);procsort;bytd;/*對自變量td進行sort由小到大排序*/procrsreg;modelf=td/lackfit;/*選項lackfit要求對有重復試驗自由度的回歸模型執行不適合度檢定(lack-of-fittest)(誤差參照為重復試驗的均方誤差),預先應先對自變量進行sort由小到大排序*/run;可以作二個因素的響應面圖(固定其它因素),E62的響應面圖如下:回歸方程作響應面,在對數據集E62進行響應面回歸(Rsreg)后得回歸方程(固定水平x3=0,本例自變量均已經編碼處理)y1=37.9-0.5625x1-0.6x2-2.875x12-3.575x1x2-3.45x22SAS程序如下:

DataE62G;Dox1=-1to1by0.1;/*當自變量是原始變量時取類似的二元網格點,即do變量=左端點

to右端點

by步長*/Dox2=-1to1by0.1;/*當自變量是原始變量時取類似的二元網格點,即do變量=左端點

to右端點

by步長*/Y1=37.9-0.5625*x1-0.6*x2-2.875*x1**2-3.575*x1*x2-3.45*x2**2;Output;End;End;Procprint;Run;/*以上為由回歸方程建立作圖數據集*/等值線圖(contourplot)37.95/*以下利用程序作三維圖*/procg3ddata=E62G;plotX2*X1=Y1;/*作x1(縱)-x2(橫)-y1(豎)三維圖/run;第八節非線性回歸非線性回歸是指模型從參數角度為非線性的且無法線性化的模型,例如y=b0+b1x1+b2x2+b11x12+b12x1x2+b22x22是諸參數bi及bij的線性模型(從自變量x1,x2角度模型是非線性的),y=aebx是可以線性化的非線性模型(線性化后lny=bx+lna),而y=aexp(-exp(b-cx))是非線性的且無法線性化的模型.此時SAS采用nlin非線性回歸程序程序形式:procnlindata=…;model因變量

=非線性函數表達式;

parms

參數1=初值參數2=初值

…;

der.參數1=偏導表達式;der.參數2=偏導表達式;……………;runs;方法:采用迭代求參數,要給出參數的一階偏導表達式der.參數=…;,要賦參數初值parms參數=…;請看例子演示例6.8.1已知牧草重量y與生長天數x的關系是y=aexp(-exp(b-cx))+ε.9次觀察的數據如表6.8.1,求a,b,c的值.datahw;inputxy;cards;數據(略)

；procnlin

data=hw;/*調用nlin進行非線性回歸*/zz=exp(b-c*x);/*簡化記號*/term=exp(-zz);/*簡化記號*/model

y=a*term;

/*因變量=非線性函數表達式,即y=aexp(-exp(b-cx))*/parmsa=70b=1.1c=0.2;/*賦參數初值*/der.a=term;/*給出參數的一階偏導表達式*/der.b=-a*zz*term;der.c=a*x*zz*term;run;

DependentVariableyMethod:Gauss-NewtonSumofIterbcaSquares01.10000.200070.00009161.21-0.1190