以數據為基礎而建立數學模型的方法稱為數據建模方法， 包括回歸、統計、機器學習、深度學習、灰色預測、主成分分析、神經網絡、時間序列分析等方法， 其中最常用的方法還是回歸方法。 本講主要介紹在數學建模中常用幾種回歸方法的 MATLAB 實現過程。

根據回歸方法中因變量的個數和回歸函數的類型（線性或非線性）可將回歸方法分為:一元線性、一元非線性、多元回歸。另外還有兩種特殊的回歸方式，一種在回歸過程中可以調整變量數的回歸方法，稱為逐步回歸，另一種是以指數結構函數作為回歸模型的回歸方法，稱為 Logistic 回歸。本講將逐一介紹這幾個回歸方法。

1.一元回歸

1.1一元線性回歸

[ 例1 ]近 10 年來，某市社會商品零售總額與職工工資總額（單位：億元）的數據見表3-1，請建立社會商品零售總額與職工工資總額數據的回歸模型。

表1 商品零售總額與職工工資總額

該問題是典型的一元回歸問題，但先要確定是線性還是非線性，然后就可以利用對應的回歸方法建立他們之間的回歸模型了，具體實現的 MATLAB 代碼如下：

（1）輸入數據

clc, clear all, closeall

x=[23.80,27.60,31.60,32.40,33.70,34.90,43.20,52.80,63.80,73.40];

y=[41.4,51.8,61.70,67.90,68.70,77.50,95.90,137.40,155.0,175.0];

（2）采用最小二乘回歸

Figure

plot(x,y,'r*')%作散點圖

xlabel('x（職工工資總額）','fontsize', 12)%橫坐標名

ylabel('y（商品零售總額）', 'fontsize',12)%縱坐標名

set(gca,'linewidth',2);

% 采用最小二乘擬合

Lxx=sum((x-mean(x)).^2);

Lxy=sum((x-mean(x)).*(y-mean(y)));

b1=Lxy/Lxx;

b0=mean(y)-b1*mean(x);

y1=b1*x+b0;

hold on

plot(x, y1,'linewidth',2);

運行本節程序，會得到如圖 1 所示的回歸圖形。在用最小二乘回歸之前，先繪制了數據的散點圖，這樣就可以從圖形上判斷這些數據是否近似成線性關系。當發現它們的確近似在一條線上后，再用線性回歸的方法進行回歸，這樣也更符合我們分析數據的一般思路。

圖1職工工資總額和商品零售總額關系趨勢圖

（3）采用 LinearModel.fit 函數進行線性回歸

m2 = LinearModel.fit(x,y)

運行結果如下：

m2 =

Linear regression model:

y ~ 1 + x1
Estimated Coefficients:

Estimate SE tStat pValue

(Intercept) -23.549 5.1028 -4.615 0.0017215

x1 2.7991 0.11456 24.435 8.4014e-09

R-squared: 0.987, Adjusted R-Squared 0.985

F-statistic vs. constant model: 597, p-value = 8.4e-09

（4）采用 regress 函數進行回歸

Y=y';

X=[ones(size(x,2),1),x'];

[b, bint, r, rint, s] = regress(Y, X)

運行結果如下：

b =

-23.5493

2.7991

在以上回歸程序中，使用了兩個回歸函數LinearModel.fit 和 regress。在實際使用中，只要根據自己的需要選用一種就可以了。函數 LinearModel.fit 輸出的內容為典型的線性回歸的參數。關于 regress，其用法多樣，MATLAB 幫助中關于 regress 的用法，有以下幾種：

b = regress(y,X)

[b,bint] = regress(y,X)

[b,bint,r] = regress(y,X)

[b,bint,r,rint] = regress(y,X)

[b,bint,r,rint,stats] = regress(y,X)

[...] = regress(y,X,alpha)

輸入 y（因變量，列向量），X（1與自變量組成的矩陣）和（alpha，是顯著性水平, 缺省時默認0.05）。

輸出

bint 是 β0，β1的置信區間，r 是殘差（列向量），rint是殘差的置信區間，s包含4個統計量：決定系數 R^2（相關系數為R），F 值，F(1,n-2) 分布大于 F 值的概率 p，剩余方差 s^2 的值。也可由程序 sum(r^2)/(n-2) 計算。其意義和用法如下：R^2 的值越接近 1，變量的線性相關性越強，說明模型有效；如果滿足

則認為變量y與x顯著地有線性關系，其中 F1-α(1，n-2)的值可查F分布表，或直接用 MATLAB 命令 finv(1-α,1, n-2) 計算得到；如果 p<α 表示線性模型可用。這三個值可以相互印證。s^2 的值主要用來比較模型是否有改進，其值越小說明模型精度越高。

1.2一元非線性回歸

在一些實際問題中，變量間的關系并不都是線性的，此時就應該用非線性回歸。用非線性回歸首先要解決的問題是回歸方程中的參數如何估計。下面通過一個實例來說明如何利用非線性回歸技術解決實例的問題。

[ 例2 ]為了解百貨商店銷售額 x 與流通率（這是反映商業活動的一個質量指標，指每元商品流轉額所分攤的流通費用）y 之間的關系，收集了九個商店的有關數據（見表2）。請建立它們關系的數學模型。

表2銷售額與流通費率數據

圖2銷售額與流通費率之間的關系圖

為了得到 x 與 y 之間的關系，先繪制出它們之間的散點圖，如圖 2 所示的“雪花”點圖。由該圖可以判斷它們之間的關系近似為對數關系或指數關系，為此可以利用這兩種函數形式進行非線性擬合，具體實現步驟及每個步驟的結果如下：

（1）輸入數據

clc, clear all, closeall

x=[1.5, 4.5, 7.5,10.5,13.5,16.5,19.5,22.5,25.5];

y=[7.0,4.8,3.6,3.1,2.7,2.5,2.4,2.3,2.2];

plot(x,y,'*','linewidth',2);

set(gca,'linewidth',2);

xlabel('銷售額x/萬元','fontsize', 12)

ylabel('流通費率y/%', 'fontsize',12)

（2）對數形式非線性回歸

m1 = @(b,x) b(1) + b(2)*log(x);
nonlinfit1 = fitnlm(x,y,m1,[0.01;0.01])
b=nonlinfit1.Coefficients.Estimate;
Y1=b(1,1)+b(2,1)*log(x);
hold on
plot(x,Y1,'--k','linewidth',2)

運行結果如下：

nonlinfit1 =

Nonlinear regression model:

y ~ b1 + b2*log(x)

Estimated Coefficients:

Estimate SE tStat pValue

b1 7.3979 0.26667 27.742 2.0303e-08

b2 -1.713 0.10724 -15.974 9.1465e-07

R-Squared: 0.973, Adjusted R-Squared 0.969

F-statistic vs. constant model: 255, p-value = 9.15e-07

（3）指數形式非線性回歸

m2 = 'y ~ b1*x^b2';
nonlinfit2 = fitnlm(x,y,m2,[1;1])
b1=nonlinfit2.Coefficients.Estimate(1,1);
b2=nonlinfit2.Coefficients.Estimate(2,1);
Y2=b1*x.^b2;
hold on
plot(x,Y2,'r','linewidth',2)
legend('原始數據','a+b*lnx','a*x^b')

運行結果如下：

nonlinfit2 =

Nonlinear regression model:

y ~ b1*x^b2

Estimated Coefficients:

Estimate SE tStat pValue

b1 8.4112 0.19176 43.862 8.3606e-10

b2 -0.41893 0.012382 -33.834 5.1061e-09

R-Squared: 0.993, Adjusted R-Squared 0.992

F-statistic vs. zero model: 3.05e+03, p-value = 5.1e-11

在該案例中，選擇兩種函數形式進行非線性回歸，從回歸結果來看，對數形式的決定系數為 0.973 ，而指數形式的為 0.993 ，優于前者，所以可以認為指數形式的函數形式更符合 y 與 x 之間的關系，這樣就可以確定他們之間的函數關系形式了。

2.多元回歸

[ 例3 ]某科學基金會希望估計從事某研究的學者的年薪 Y 與他們的研究成果（論文、著作等）的質量指標 X1、從事研究工作的時間 X2、能成功獲得資助的指標 X3之間的關系，為此按一定的實驗設計方法調查了 24 位研究學者，得到如表3 所示的數據（ i 為學者序號），試建立 Y 與 X1, X2, X3之間關系的數學模型，并得出有關結論和作統計分析。

表3從事某種研究的學者的相關指標數據

該問題是典型的多元回歸問題，但能否應用多元線性回歸，最好先通過數據可視化判斷他們之間的變化趨勢，如果近似滿足線性關系，則可以執行利用多元線性回歸方法對該問題進行回歸。具體步驟如下：

（1）作出因變量 Y 與各自變量的樣本散點圖

作散點圖的目的主要是觀察因變量 Y 與各自變量間是否有比較好的線性關系，以便選擇恰當的數學模型形式。圖3 分別為年薪 Y 與成果質量指標 X1、研究工作時間 X2、獲得資助的指標 X3之間的散點圖。從圖中可以看出這些點大致分布在一條直線旁邊，因此，有比較好的線性關系，可以采用線性回歸。繪制圖3的代碼如下：

subplot(1,3,1),plot(x1,Y,'g*'),

subplot(1,3,2),plot(x2,Y,'k+'),

subplot(1,3,3),plot(x3,Y,'ro'),

圖3因變量Y與各自變量的樣本散點圖

（2）進行多元線性回歸

這里可以直接使用 regress 函數執行多元線性回歸，具體代碼如下：

x1=[3.5 5.3 5.1 5.8 4.2 6.0 6.8 5.5 3.1 7.2 4.5 4.9 8.0 6.5 6.5 3.7 6.2 7.0 4.0 4.5 5.9 5.6 4.8 3.9];

x2=[9 20 18 33 31 13 25 30 5 47 25 11 23 35 39 21 7 40 35 23 33 27 34 15];

x3=[6.1 6.4 7.4 6.7 7.5 5.9 6.0 4.0 5.8 8.3 5.0 6.4 7.6 7.0 5.0 4.0 5.5 7.0 6.0 3.5 4.9 4.3 8.0 5.0];

Y=[33.2 40.3 38.7 46.8 41.4 37.5 39.0 40.7 30.1 52.9 38.2 31.8 43.3 44.1 42.5 33.6 34.2 48.0 38.0 35.9 40.4 36.8 45.2 35.1];

n=24; m=3;

X=[ones(n,1),x1',x2',x3'];

[b,bint,r,rint,s]=regress(Y',X,0.05);

運行后即得到結果如表4所示。

表4對初步回歸模型的計算結果

計算結果包括回歸系數 b = (β0,β1,β2,β3) = (18.0157, 1.0817, 0.3212, 1.2835)、回歸系數的置信區間，以及統計變量 stats（它包含四個檢驗統計量：相關系數的平方R^2，假設檢驗統計量 Ｆ，與 F 對應的概率 p，s^2的值）。因此我們得到初步的回歸方程為：

由結果對模型的判斷：

回歸系數置信區間不包含零點表示模型較好，殘差在零點附近也表示模型較好，接著就是利用檢驗統計量 Ｒ，Ｆ，p 的值判斷該模型是否可用。

1）相關系數 Ｒ 的評價：本例 Ｒ 的絕對值為 0.9542 ，表明線性相關性較強。

2）F 檢驗法：當 F > F1-α(m,n-m-1) ，即認為因變量 y 與自變量x1,x2,...,xm之間有顯著的線性相關關系；否則認為因變量 y 與自變量x1,x2,...,xm之間線性相關關系不顯著。本例 Ｆ＝67.919 > F1-0.05( 3,20 ) = 3.10。

3）p 值檢驗：若 p < α（α 為預定顯著水平），則說明因變量 y 與自變量?x1,x2,...,xm之間顯著地有線性相關關系。本例輸出結果，p<0.0001,顯然滿足 p<α=0.05。

以上三種統計推斷方法推斷的結果是一致的，說明因變量 y與自變量之間顯著地有線性相關關系，所得線性回歸模型可用。s^2當然越小越好，這主要在模型改進時作為參考。

3.逐步歸回

[ 例4 ]（Hald,1960）Hald 數據是關于水泥生產的數據。某種水泥在凝固時放出的熱量 Y（單位：卡/克）與水泥中 4 種化學成品所占的百分比有關：

在生產中測得 12 組數據，見表5，試建立 Y關于這些因子的“最優”回歸方程。

表5水泥生產的數據

對于例 4 中的問題，可以使用多元線性回歸、多元多項式回歸，但也可以考慮使用逐步回歸。從逐步回歸的原理來看，逐步回歸是以上兩種回歸方法的結合，可以自動使得方程的因子設置最合理。對于該問題，逐步回歸的代碼如下：

X=[7,26,6,60;1,29,15,52;11,56,8,20;11,31,8,47;7,52,6,33;11,55,9,22;3,71,17,6;1,31,22,44;2,54,18,22;21,47,4,26;1,40,23,34;11,66,9,12];%自變量數據

Y=[78.5,74.3,104.3,87.6,95.9,109.2,102.7,72.5,93.1,115.9,83.8,113.3];%因變量數據

Stepwise（X,Y,[1,2,3,4],0.05,0.10）% in=[1,2,3,4]表示X1、X2、X3、X4均保留在模型中

程序執行后得到下列逐步回歸的窗口，如圖 4 所示。

圖4逐步回歸操作界面

在圖 4 中，用藍色行顯示變量 X1、X2、X3、X4 均保留在模型中，窗口的右側按鈕上方提示：將變量X3剔除回歸方程（Move X3 out），單擊 Next Step 按鈕，即進行下一步運算，將第 3 列數據對應的變量 X3 剔除回歸方程。單擊 Next Step 按鈕后，剔除的變量 X3 所對應的行用紅色表示，同時又得到提示：將變量 X4 剔除回歸方程（Move X4 out），單擊 Next Step 按鈕，這樣一直重復操作，直到 “Next Step” 按鈕變灰，表明逐步回歸結束，此時得到的模型即為逐步回歸最終的結果。

4.Logistic 回歸

[ 例5 ]企業到金融商業機構貸款，金融商業機構需要對企業進行評估。評估結果為 0 , 1 兩種形式，0 表示企業兩年后破產，將拒絕貸款，而 1 表示企業 2 年后具備還款能力，可以貸款。在表 6 中，已知前 20 家企業的三項評價指標值和評估結果，試建立模型對其他 5 家企業（企業 21-25）進行評估。

表6企業還款能力評價表

對于該問題，很明顯可以用 Logistic 模型來回歸，具體求解程序如下：

% logistic回歸Matlab實現程序

%% 數據準備

clc, clear, closeall

X0=xlsread('logistic_ex1.xlsx', 'A2:C21');% 回歸模型的輸入

Y0=xlsread('logistic_ex1.xlsx', 'D2:D21');% 回歸模型的輸出

X1=xlsread('logistic_ex1.xlsx', 'A2:C26');% 預測數據輸入

%% logistics函數

GM = fitglm(X0,Y0,'Distribution','binomial');

Y1 = predict(GM,X1);

%% 模型的評估

N0 =1:size(Y0,1); N1= 1:size(Y1,1);

plot(N0', Y0, '-kd');

hold on; scatter(N1', Y1, 'b')

xlabel('數據點編號'); ylabel('輸出值');

得到的回歸結果與原始數據的比較如圖5所示。

圖5 回歸結果與原始數據的比較圖

5.小結

本講主要介紹數學建模中常用的幾種回歸方法。在使用回歸方法的時候，首先可以判斷自變量的個數，如果超過 2 個，則需要用到多元回歸的方法，否則考慮用一元回歸。然后判斷是線性還是非線性，這對于一元回歸是比較容易的，而對于多元，往往是將其他變量保持不變，將多元轉化為一元再去判斷是線性還是非線性。如果變量很多，而且復雜，則可以首先考慮多元線性回歸，檢驗回歸效果，也可以用逐步回歸。總之，用回歸方法比較靈活，根據具體情景還是比較容易找到合適的方法的。