SAS學(xué)習(xí)系列25. 非線性回歸

25.非線性回歸現(xiàn)實(shí)世界中嚴(yán)格的線性模型并不多見，它們或多或少都帶有某種程度的近似；在不少情況下，非線性模型可能更加符合實(shí)際。對(duì)變量間非線性相關(guān)問題的曲線擬合，處理的方法主要有：（1）首先確定非線性模型的函數(shù)類型，對(duì)于其中可線性化問題則通過變量變換將其線性化，從而歸結(jié)為前面的多元線性回歸問題來解決；（2）若實(shí)際問題的曲線類型不易確定時(shí)，由于任意曲線皆可由多項(xiàng)式來逼近，故?？捎枚囗?xiàng)式回歸來擬合曲線；（3）若變量間非線性關(guān)系式已知（多數(shù)未知），且難以用變量變換法將其線性化，則進(jìn)行數(shù)值迭代的非線性回歸分析。一）可變換為線性的非線性回歸表1典型的函數(shù)及線性化方法"窗數(shù)窯稱.函數(shù)盍達(dá)式中線性化方法.取曲線函數(shù)J1b—=d——.■y工V——u——V工V=111JU-111X指數(shù)函數(shù).y=ae^?v=Inyu=Xry=◎嚴(yán)-X對(duì)數(shù)函數(shù).y=戊+&Ln兀"v=yu=InXh，W型函數(shù)十v=—u=e~x4

在很多場(chǎng)合，可以對(duì)非線性模型進(jìn)行線性化處理，尤其是可變換為線性的非線性回歸，運(yùn)用最小二乘法進(jìn)行推斷，對(duì)線性化后的線性模型，可以應(yīng)用REG過程步進(jìn)行計(jì)算。例1有實(shí)驗(yàn)數(shù)據(jù)如下：X1.11.21.31.41.51.61.71.81.922.12.22.32.4Y109.9540.4520.0924.5311.027.394.952.721.821.490.820.30.20.22試分別采用指數(shù)回歸（y=aebx）方法進(jìn)行回歸分析。代碼：dataexam25_1;inputxy;cards?;1.1109.951.240.451.320.091.424.531.511.021.67.391.74.951.82.721.91.8221.492.10.822.20.32.30.22.4?0.22;run;procsgplotdata=exam251;scatterx=xy=y;run;proccorrdata=exam251;varxy;run;datanewl;setexam25_1;v=log(y);run;procsgplotdata=newl;scatterx=xy=v;title'變量代換后數(shù)據(jù)'；run;procregdata=newl;varxv;modelv=x;printcli;title'殘差圖';plotresidual.*predicted.;run;datanew2;setexam25_1;y1=14530?28*exp(-4.73895*x);run;procgplotdata=new2;ploty*x=1y1*x=2/overlay;symbolv=doti=nonecv=red;symbol2i=smcolor=blue;title'指數(shù)回歸圖'；run;運(yùn)行結(jié)果：

CORR過程2變量：kV筒單統(tǒng)計(jì)量N均值標(biāo)準(zhǔn)差總和最小值最大值XI41.750000.4133324500001100002.40000y16.13929294SD3722595DOOo.?noon109.95000Pearson相關(guān)系數(shù)，N_14

Prob>|r|underHO:Hho=0Xy1.00000■0.711040.0044-[)./11IMrooooo0.0()14變咼忙膜后歆據(jù)讀取的觀測(cè)數(shù)14使用的觀測(cè)數(shù)E方差分析源自由度平方和均方F值Pr>F模型151.0910651.09106750.31<.0001誤差1?0.810630.血麻校正合計(jì)1351.90170均方根逞差0.25991R方09844國(guó)塑量均值1.29084調(diào)整R方0.9631變異系數(shù)20,13498參數(shù)怙計(jì)值自由度參數(shù)估計(jì)值準(zhǔn)差標(biāo)誤L值Pr>|t|Intercept19.583990.309-1530.97<.oaoiX1-473895().1m?-兒b()<.()001輸出統(tǒng)計(jì)星觀測(cè)因變量標(biāo)準(zhǔn)誤差均值預(yù)測(cè)95%GL預(yù)測(cè)殘差1-1.70001S7I20.I31S3.了砂5.00810.32892U./U013.89730.11753.27584.5137-0.1^/233.000230.KM12.0133-1一（）畑】-0,4231d3.19992.9+950.03202.34S73.55020.250452.39972.-1/bG0.08171.88193.0692-0.U/bLI62.0001?.00170.07411.4l?fi2.5305-0.00154971.59941.527B0.U/000.㈣132.11430.0/1681.00061.05390.(17000.-IfiM1.fi40d-0.0533甘0.59300.58000.0741-0.0088761.ICS90.OIOS100.3988[].10N10.0811-0.48750.6997(1.292711-0.1995-0.367(50.OJ1PO-096S50.23290.IC9312-1.2040-0.841f0.1041-1.4517-0.2317-0.362313-1.6094-1.3156o.urn-1.-o.cfi-n-0.加豹1d-1.5141-I.73950.1310-2.4244-1.15450.2754□.Hi-aa.H44-□.H31TOC\o"1-5"\h\z測(cè)U境掘2唉芒口由屢12N5EDDB75F?方0嗆44出車R冉DWJ11.82.12.2程序說明：（1）調(diào)整后的R2=0.9831，說明擬合程度很好；F檢驗(yàn)的P值=0.0001<a=0.05，拒絕原假設(shè)，故直線回歸的斜率不為0；（2）將線性回歸系數(shù)代入，得到原回歸方程y=14530.28*e-4.73895x（3）殘差圖趨勢(shì)，符合殘差隨機(jī)正態(tài)分布的假設(shè)（不帶其它明顯趨勢(shì)）。二、多項(xiàng)式回歸一般函數(shù)都可用多項(xiàng)式來逼近，故多項(xiàng)式回歸分析可用來處理相當(dāng)廣泛的非線性問題。對(duì)觀測(cè)數(shù)據(jù)（xt,yt）,t=l,...,N.多項(xiàng)式回歸模型為：yt—+力］兀+玄X；+…++str匸12…N"令則模型可寫為：Y=XB+s當(dāng)X列滿秩時(shí)，用最小二乘估計(jì)B=（XX）-1XY可求得其多項(xiàng)式回歸方程。但由于（XX）-1的計(jì)算既復(fù)雜又不穩(wěn)定，故一般采用正交多項(xiàng)式法來進(jìn)行多項(xiàng)式回歸。多項(xiàng)式模型可以直接應(yīng)用GLM（廣義線性模型）求解。例2重慶市種畜場(chǎng)奶牛群1—12月份(xl)，產(chǎn)犢母牛平均產(chǎn)奶量(y)的資料如程序數(shù)據(jù)步中，試對(duì)該資料配置一個(gè)合適的回歸方程。代碼：dataexam25_2;inputxly@@;x2=x1*x1;datalines;3833.433476.763811.583466.223769.473395.423565.743807.083481.993817.033372.823884.52?run;procsgplotdata=exam25_2;scatterx=x1y=y;title'原始數(shù)據(jù)散點(diǎn)圖'；run;procregdata=exam25_2;modely=x1x2;run;運(yùn)行結(jié)果：

原蠟妁埠敬點(diǎn)圖3WD3TOD2.57.510.0原蠟妁埠敬點(diǎn)圖3WD3TOD2.57.510.0REfl過程

核型:MODEI1

因?qū)D量：/

讀取的觀測(cè)數(shù)吃

使用的觀測(cè)裁壯方差分析源自由奩平方和均方F值Pr>F模型233259416629/1G.890.0009誤差988S019844.56254校正合計(jì)II42II95均方根溟差99.21977R方0【伽因鑾量均值3640.17167調(diào)整R方0.7429蠻異系數(shù)2.72509參數(shù)估計(jì)值翌量自由度參數(shù)怙計(jì)値標(biāo)準(zhǔn)溟差土值Pr>|t|Intercept14117.20136102.54G4940.15<.0001xl1-2049366836.2S837-5.650.0003115.785702,715895.810.om

GLM過程源自由度平方和均方1值Pr>1模型332594.3259166297,1630n.nans*1=1主恢差88601.06299844,5625樓正合計(jì)11421195.3388R方變異系數(shù)根MSt丫均值0.7896442.72550999.219773040.172源自由度I型SS均方F值Pr>FX1I11.005011.00500.000.9741x1*x1I332583.3210332583.321033.780.0003源自由庚I11型均方卜值Pr>1-xl1314325.5414314325.541431.930.0003x1*x11332583.3210332593.32100.0003參數(shù)怙計(jì)值標(biāo)準(zhǔn)誤差tfiPr>|t|IntercEipt4117.201364102.546490540.15<.0001xl-204.93607836.2683725-A650.CKXJ3xl+x115.7856992.7158977L).81U.ZUS程序說明：（1）觀察數(shù)據(jù)的散點(diǎn)圖，更適合二次多項(xiàng)式擬合，也可以測(cè)試幾種不同次數(shù)的多項(xiàng)式擬合選擇其中最優(yōu)的；（2）將回歸系數(shù)代入多項(xiàng)式方程得到：y=4117.20136-204.93668x1+15.78570x12三、不能變換為線性的非線性回歸該類非線性回歸分析就是利用最小二乘準(zhǔn)則來估計(jì)回歸系數(shù)卩，使得殘差平方和最小。一般來用數(shù)值迭代法來進(jìn)行，先選定回歸系數(shù)的初值卩0,按照給定的步長(zhǎng)和搜索方向逐步迭代，直到殘差平方和達(dá)到最小。有5種常用的非線性回歸迭代方法:高斯-牛頓法（Gauss-Newton）、最速下降法（梯度法）、牛頓法（Newton）、麥夸特法（Marquardt）、正割法（DUD）。高斯-牛頓法在初值選取適當(dāng)，且可逆時(shí)非常有效，但在其他情形，其求解較為困難，對(duì)此，Marguardt對(duì)其中的正則系數(shù)陣作適當(dāng)修正，得到了改進(jìn)算法。(二)PROCNLIN過程步對(duì)于不能線性化的非線性模型。其估計(jì)不能直接運(yùn)用經(jīng)典的最小二乘法，而需要運(yùn)用其他估計(jì)方法，如加權(quán)最小二乘法、直接搜索法、直接最優(yōu)法與Taylor級(jí)數(shù)展開法進(jìn)行線性逼近。此時(shí)，可以利用NLIN過程步實(shí)現(xiàn)相應(yīng)的計(jì)算，它是采用最小誤差平方法及迭代推測(cè)法來建立一個(gè)非線性模型，估計(jì)參數(shù)默認(rèn)采用高斯-牛頓迭代法。NLIN過程不保證一定可以算出符合最小誤差平方法之標(biāo)準(zhǔn)的參數(shù)估計(jì)值。基本語法：PROCNLINdata=數(shù)據(jù)集</可選項(xiàng)>;PARMS參數(shù)名=數(shù)值;MODEL因變量=表達(dá)式</可選項(xiàng)>;BOUNDS不等式;>CONTROL變量;>der.參數(shù)名=表達(dá)式;><OUTPUTout=輸出數(shù)據(jù)集</可選統(tǒng)計(jì)量>;>說明：NLIN的可選項(xiàng)包括outest=輸出數(shù)據(jù)集輸出每步迭代的結(jié)果；best=n只輸出最好的n組殘差平方和；method=gauss|marquardt|newton|gradient]dud|設(shè)定參數(shù)估計(jì)的迭代方法，默認(rèn)為gauss(沒有der.語句)；(2)PARMS語句指定參數(shù)并賦值，一般包括參數(shù)名、初始值(GridSearch可以幫助選擇合適的初始值)、迭代準(zhǔn)則；例如：parmsb0=0b1=1to10b2=1to10by2b3=1,10,100;bounds語句用于設(shè)定參數(shù)的約束，主要是不等式約束，約束間用逗號(hào)分隔。例如，boundsa<=20,b>30,1<=c<=10;der.語句用于計(jì)算模型關(guān)于各參數(shù)的偏導(dǎo)數(shù)，相應(yīng)格式為:一階偏導(dǎo)數(shù)：der.參數(shù)名=表達(dá)式；二階偏導(dǎo)數(shù)：der.參數(shù)名.參數(shù)名=表達(dá)式；例如，對(duì)于模型modely=b0*(1—exp(-b1*x));二階偏導(dǎo)數(shù)表達(dá)式：der?b0?b1=x*exp(-b1*x);例3根據(jù)對(duì)已有數(shù)據(jù)的XY散點(diǎn)圖的觀察和分析，發(fā)現(xiàn)Y隨X增長(zhǎng)趨勢(shì)是減緩的，并且Y趨向一個(gè)極限值，我們認(rèn)為用負(fù)指數(shù)增長(zhǎng)曲線來擬合模型較為合適。代碼：dataexpd;inputxy@@;datalines;0200.570300.720400.810500.870600.910700.940800.950900.971000.981100.991201.001300.991400.991501.001601.001700.991801.001901.002000.992101.00procnlindata=expdbest=10method=gauss;parmsb0=0to2by0.5b1=0.01to0.09by0.01;modely=b0*(1-exp(-b1*x));der.b0=1-exp(-b1*x);der.b1=b0*x*exp(-b1*x);outputout=expoutp=ygs;run;goptionsreset=globalgunit=pctcback=whiteborderhtitle=6htext=3ftext=swissbcolors=(back);procgplotdata=expout;ploty*xygs*x/haxis=axis1vaxis=axis2overlay;symbol1i=nonev=pluscv=redh=2.5w=2;symbol2i=joinv=nonel=1h=2.5w=2;axis1order=20to210by10;axis2order=0.5to1?1by0.05;title1'y=b0*(1-exp(-b1*x)';title2'procnlinmethod=gauss';run;運(yùn)行結(jié)果：NLIN過程因姿量網(wǎng)格搜索hOhi平方和1.【)1燦00400I).[)01d()1.0000o.oroo0.1.0000n.0.05521.00000.03000.06661.00000.07000.09731.00000.08000.13651.00000.09000.17081.￡)0(10t).矽E】o.lias1.MXN)0.01000.1.()0(J()001002.If曲

NLIN過程因變量7

方法：Gauss-Newton迭代階段迭代bOb1平方和0I00000.04000.0014010.99(510.0-1190.000^8020.S9620.00057730.S9620.04200.000577斗0,9962NOIL;Convergencecriterionmet.估計(jì)匯總方法Uauss-NciMton進(jìn)代1R2.E52E1PPG(b1)1.O2ftE-0RPC(b1)394E7Object2.b(5F-1()□標(biāo)0.000577讀取的觀測(cè)2014用的觀測(cè)2()缺失的觀測(cè)UNoteNote;AninterueplwtisriuLspeci1iedlotthismodeI.源自由度平方和均方FVa1ue近似Pr>F模型217.67178.8359275733<.0001誤差180.0005770.000032未校正合計(jì)2017.6723參數(shù)怙計(jì)值近似標(biāo)準(zhǔn)淒差近似95%査信限bO0.馭為0.001610.99280