Matlab多元線性回歸預(yù)測(cè)城市建成區(qū)面積

曹應(yīng)舉，張永彬，溫繼滿，汲 姣，劉佳麗

(1.華北理工大學(xué) 礦業(yè)工程學(xué)院，河北 唐山 063210；2.開灤安全技術(shù)培訓(xùn)中心，河北 唐山 063000)

0 引言

隨著經(jīng)濟(jì)的發(fā)展，城市化速度加快，城市的面積在不斷的擴(kuò)張。為了了解城市的擴(kuò)張動(dòng)態(tài)，獲取城市擴(kuò)張的歷史軌跡，許多學(xué)者對(duì)城市建成區(qū)面積的預(yù)測(cè)問題進(jìn)行了研究：李愛民利用遙感技術(shù)對(duì)不同年份城市建成區(qū)多期遙感影像進(jìn)行提取分類，分析了城市擴(kuò)張的時(shí)空規(guī)律[1]；劉柯運(yùn)用BP神經(jīng)網(wǎng)絡(luò)方法建立預(yù)測(cè)模型，使用多期數(shù)據(jù)作為學(xué)習(xí)樣本和檢驗(yàn)樣本，對(duì)2005年北京市城市建成區(qū)面積進(jìn)行了模擬預(yù)測(cè)[2]；雷波通過建立多元回歸模型，選用13個(gè)社會(huì)和經(jīng)濟(jì)驅(qū)動(dòng)因子，對(duì)福州城市建成區(qū)面積的擴(kuò)張驅(qū)動(dòng)進(jìn)行了回歸分析[3]。在模擬預(yù)測(cè)方面，多元回歸模型是一種常用的預(yù)測(cè)模型，其預(yù)測(cè)精度高、應(yīng)用領(lǐng)域廣泛，例如，周永生等通過綜合影響糧食產(chǎn)量的多種因素，將多元回歸分析應(yīng)用于糧食產(chǎn)量的預(yù)測(cè)中[4]；付倩嬈通過對(duì)不同大氣成分濃度數(shù)據(jù)進(jìn)行分析，將多元回歸分析應(yīng)用在霧霾的預(yù)測(cè)中[5]；葉鋒通過綜合考慮經(jīng)濟(jì)技術(shù)等多種影響因子，將多元回歸分析應(yīng)用于油田產(chǎn)量的預(yù)測(cè)中[6]；韋浩采用不同方法對(duì)滑坡距離進(jìn)行預(yù)測(cè)，得出結(jié)論是“與傳統(tǒng)預(yù)測(cè)方法相比，多元回歸分析法所建立的預(yù)測(cè)模型精度較高”[7]。但傳統(tǒng)的用純數(shù)學(xué)方法來實(shí)現(xiàn)多元線性回歸方程求解的過程比較繁瑣，為了簡(jiǎn)化計(jì)算難度，提高計(jì)算的速度，本文以Matlab為語(yǔ)言平臺(tái)、以回歸分析為數(shù)學(xué)統(tǒng)計(jì)方法，建立關(guān)于建成區(qū)面積與其影響因子的多元回歸模型，并應(yīng)用其進(jìn)行城市建成區(qū)面積的預(yù)測(cè)。

1 線性回歸模型

變量Y的值主要受影響因子X(由X1，X2，…，Xk確定，可以表示為X1，X2，…，Xk的某個(gè)函數(shù)關(guān)系式：Y=f(X1，X2，…，Xk))和隨機(jī)誤差ε的影響，本研究將自變量寫成如下形式：

Y=f(X1，X2，…，Xk)ε。

(1)

其中，隨機(jī)變量Y稱為被解釋變量或因變量；X1，X2，…，Xk稱為解釋變量或自變量。f(X1，X2，…，Xk)為一般變量X1，X2，…，Xk的確定性關(guān)系，ε為隨機(jī)誤差。

當(dāng)概率模型式(1)中回歸函數(shù)為線性函數(shù)時(shí)，即有

Y=β0+β2X1+…βkXk+ε。

(2)

其中β0，β1，β2，…，βk是k+1個(gè)未知參數(shù)，β0稱為回歸常數(shù)，β1，β2，…，βk稱為回歸系數(shù)，Y稱為被解釋變量(因變量)，而X1，X2，…，Xk是k個(gè)可以控制的一般變量，稱為解釋變量(自變量)。k=1時(shí)，為一元線性回歸模型；k≥2時(shí)，稱為多元線性回歸模型。

線性回歸模型的“線性”是針對(duì)未知參數(shù)βi(i=0，1，2，…，k)而言的。對(duì)于回歸解釋變量的線性是非本質(zhì)的，因?yàn)榻忉屪兞渴欠蔷€性的，?？梢酝ㄟ^變量的替換把它轉(zhuǎn)化成線性的。

若(Xi1，Xi2，…，Xik；Yi)，i=1，2，…，n是式(2)中變量(X1，X2，…，Xk；Y)的一組觀測(cè)值，則線性回歸模型可表示為

Yi=β0+β1Xi1+…βkXik+εi。

(3)

其中E(Yi)=β0+β1Xi1+…βkXik為多元回歸方程組，其相應(yīng)的矩陣表達(dá)式為

Y=xβ。

(4)

其中回歸系數(shù)的最小二乘估計(jì)為

(5)

2 建成區(qū)面積案例分析

通過搜集《中國(guó)統(tǒng)計(jì)年鑒》(2016)和地方統(tǒng)計(jì)年鑒整理獲得25個(gè)省會(huì)城市2015年的建成區(qū)面積。鑒于影響建成區(qū)面積有諸多因素，本文在參照已有研究成果的基礎(chǔ)上，結(jié)合實(shí)際情況，構(gòu)建了城市建成區(qū)面積的經(jīng)濟(jì)社會(huì)影響因子體系(如圖1所示)，并通過建立多元回歸模型進(jìn)行分析和檢驗(yàn)，從而對(duì)未來的城市建成區(qū)面積進(jìn)行預(yù)測(cè)。

圖1 建成區(qū)面積影響因子

3 數(shù)據(jù)處理

Matlab多元線性回歸數(shù)據(jù)處理過程包括數(shù)據(jù)預(yù)處理、建模、分析及預(yù)測(cè)四大部分，其具體步驟如圖2所示。

圖2 線性回歸流程圖

3.1 利用regress函數(shù)實(shí)現(xiàn)多元線性回歸

在Matlab中使用命令regress實(shí)現(xiàn)多元線性回歸，調(diào)用格式為：b=regress(Y，X)，[b，bint，r，rint，stats]=regress(Y，X，alpha)[8]，其中：Y表示一個(gè)n-1的因變量數(shù)據(jù)矩陣；X是n-p矩陣，自變量X是一列具有相同行數(shù)，值是1的矩陣的組合；alpha為顯著性水平(缺省時(shí)設(shè)定為0.05)；輸出向量b為回歸系數(shù)最小二乘估計(jì)值；bint為b的置信區(qū)間；r，rint為殘差及其置信區(qū)間。

stats是用于檢驗(yàn)回歸模型的統(tǒng)計(jì)量，第一個(gè)是R2，其中R是相關(guān)系數(shù)；第二個(gè)是F統(tǒng)計(jì)量值；第三個(gè)是與統(tǒng)計(jì)量F對(duì)應(yīng)的概率P；第四個(gè)是S2，為誤差方差估計(jì)值。R2越接近1，說明回歸方程越顯著；F>F1-alpha(P，n-p-1)時(shí)拒絕H0，F(xiàn)越大，回歸方程越顯著；與F對(duì)應(yīng)的概率P

3.2 原始數(shù)據(jù)歸一化

為了方便數(shù)據(jù)處理，需要對(duì)樣本中的各影響因子進(jìn)行歸一化處理。歸一化是一種消除指標(biāo)之間量綱影響的簡(jiǎn)化計(jì)算方法，常用的數(shù)據(jù)歸一化方法有“最小-最大標(biāo)準(zhǔn)化”“Z-score標(biāo)準(zhǔn)化”和“按小數(shù)定標(biāo)標(biāo)準(zhǔn)化”等，本文采用的是“最小-最大標(biāo)準(zhǔn)化”方法對(duì)原數(shù)據(jù)進(jìn)行線性變換。將A中的一個(gè)原始值X通過“最小-最大標(biāo)準(zhǔn)化”映射成在區(qū)間[0，1]中的值X′，其形式如式(6)所示。

(6)

其中：X′為歸一化后的數(shù)據(jù)，X為原始數(shù)據(jù)，Xmin和Xmax分別是X的最小值和最大值。

3.3 整體最小二乘(TLS)多元回歸

通常認(rèn)為數(shù)據(jù)矩陣X是給定的，不存在誤差，但是如果數(shù)據(jù)矩陣X也存在誤差或者擾動(dòng)，那么最小二乘估計(jì)從統(tǒng)計(jì)觀點(diǎn)看就不再是最優(yōu)的，它將是有偏的，而且偏差的協(xié)方差將由于X的噪聲誤差的作用而增加。因此，當(dāng)X也存在誤差時(shí)，應(yīng)該使用整體最小二乘進(jìn)行回歸[9]。

(7)

或

(8)

(9)

或等價(jià)為

(B+D)Z=0。

(10)

‖D‖F(xiàn)=min。

(11)

其中‖D‖F(xiàn)是D的F(Frobenius)范數(shù)。

3.4 相關(guān)性檢驗(yàn)

若想使得所擬直線有實(shí)際意義，必須保證建成區(qū)面積(變量Y)與其影響因子(自變量X)存在線性相關(guān)性[10]，描述它們之間相關(guān)性系數(shù)的定義為：

(12)

其估值為：

(13)

當(dāng)ρ越接近±1時(shí)，表明隨機(jī)變量Y與X的相關(guān)性越密切，即所建立的線性模型和實(shí)際的試驗(yàn)情況越接近。

在Matlab中使用corrcoef函數(shù)可以求兩個(gè)序列的相關(guān)度，corrcoef(X，Y)表示序列X和序列Y的相關(guān)系數(shù)，得到的結(jié)果是一個(gè)2*2矩陣。

相關(guān)系數(shù)的大小所表示的意義通常如表1所示。

表1 相關(guān)系數(shù)表示的意義

依據(jù)相關(guān)系數(shù)表，編寫相應(yīng)程序代碼，將相關(guān)程度在微相關(guān)范圍內(nèi)的影響因子(戶籍人口數(shù)量X1，第一產(chǎn)業(yè)X3，人均生產(chǎn)總值X8，衛(wèi)生機(jī)構(gòu)數(shù)量X9，普通高校數(shù)量X10，公共圖書館數(shù)量X11)予以剔除，保留相關(guān)系數(shù)大于0.3的影響因子，剩余5個(gè)影響因子(生產(chǎn)總值X2，第二產(chǎn)業(yè)X4，工業(yè)X5，建筑業(yè)X6，第三產(chǎn)業(yè)X7)的相關(guān)程度顯著，對(duì)因變量建成區(qū)面積Y的解釋程度較高。

創(chuàng)建城市建成區(qū)面積Y和5個(gè)影響因子(生產(chǎn)總值X2，第二產(chǎn)業(yè)X4，工業(yè)X5，建筑業(yè)X6，第三產(chǎn)業(yè)X7)的一元線性回歸預(yù)測(cè)方程，如表2所示。

表2 一元線性回歸預(yù)測(cè)方程

根據(jù)表2的擬合方程結(jié)果來看，建成區(qū)面積Y與各個(gè)影響因子X的擬合效果良好。對(duì)于R2來說數(shù)值越大擬合效果越好，各個(gè)方程的判定系數(shù)最低為0.527，總體來說判定系數(shù)比較高，所以建成區(qū)面積與各個(gè)影響因子擬合效果總體較好。為了進(jìn)一步探討建成區(qū)面積與各個(gè)影響因子之間的關(guān)系，以建成區(qū)面積Y為因變量，以影響因子X2-X7為自變量進(jìn)行多元線性回歸分析。

繪制因變量Y與自變量X2-X7之間的散點(diǎn)圖，如圖3所示。

圖3 Y與X2，X4，X5，X6，X7的散點(diǎn)圖

模型建立過程中繪制的殘差圖會(huì)有異常點(diǎn)而影響模型的正確性，需要對(duì)異常點(diǎn)進(jìn)行剔除，利用Matlab編寫簡(jiǎn)單的循環(huán)語(yǔ)句可以實(shí)現(xiàn)以上操作，剔除完成后繪制殘差圖如圖4所示。

圖4 異常點(diǎn)剔除前后的殘差圖

剔除殘差后利用regress再次求解參數(shù)，結(jié)果如表3所示。

表3 參數(shù)求解

t檢驗(yàn)是逐一對(duì)參數(shù)的顯著性進(jìn)行檢驗(yàn)。其原理是|t(f)|>tα/2時(shí)，接受H0，查表可得，tα/2=2.262，經(jīng)檢驗(yàn)，t檢驗(yàn)統(tǒng)計(jì)量小于tα/2的常數(shù)項(xiàng)以及自變量X6顯著性不明顯，應(yīng)予以剔除。stats中第1個(gè)參數(shù)R2是回歸平方和與總離差平方和的比值，其值越大越好，該模型達(dá)到0.979 9；第2個(gè)參數(shù)f統(tǒng)計(jì)量，越大越好，本模型為87.546 5；第3個(gè)參數(shù)為P的顯著性概率，應(yīng)該小于0.05，越接近0越好，本模型基本為0；第4個(gè)參數(shù)為估計(jì)誤差方差，本模型估計(jì)誤差方差為0.002，綜上，stats中4個(gè)參數(shù)充分說明回歸方程顯著，該回歸模型成立。

因此，最終建立的多元線性回歸方程式為：

Y=1.746X2+0.519X4+0.388X5+1.856X7。

3.5 結(jié)果預(yù)測(cè)

對(duì)城市建成區(qū)面積進(jìn)行回歸分析的目的主要是進(jìn)行預(yù)測(cè)和控制[11]。未來城市建成區(qū)的面積通過單一的往年建成區(qū)增減面積來預(yù)測(cè)難以定性和定量，所以在實(shí)際生活中，都是通過近幾年的數(shù)據(jù)構(gòu)建主要影響因子和建成區(qū)面積的方程，通過觀測(cè)影響因子，來預(yù)測(cè)未來的城市建成區(qū)面積。

Matlab自帶程序rstool可以實(shí)現(xiàn)回歸分析的控制預(yù)測(cè)功能，在本例中選擇5個(gè)城市的建成區(qū)面積留做預(yù)測(cè)，選出未參與回歸處理的5組數(shù)據(jù)，將其自變量影響因子X輸入相應(yīng)位置，Matlab將自動(dòng)計(jì)算預(yù)測(cè)結(jié)果，如圖5所示。

圖5 rstool預(yù)測(cè)處理

將Matlab計(jì)算結(jié)果與經(jīng)過歸一化后的真實(shí)值進(jìn)行比對(duì)，預(yù)測(cè)值在限差范圍內(nèi)，說明了回歸方程的正確性，結(jié)果比對(duì)如表4所示。

表4 實(shí)際值與預(yù)測(cè)值比對(duì)

結(jié)果顯示，第一、三、四、五組數(shù)據(jù)的實(shí)際值都在預(yù)測(cè)區(qū)間內(nèi)，預(yù)測(cè)的可靠性高，第二組數(shù)據(jù)預(yù)測(cè)結(jié)果不理想，導(dǎo)致該結(jié)果的原因可能是由于各種因素之間產(chǎn)生了更為復(fù)雜的相互作用，使得城市建成區(qū)面積不再表現(xiàn)為線性關(guān)系中那種按比例的規(guī)則變化，而代之以不按比例、不規(guī)則的變化或突變。例如，城市區(qū)劃調(diào)整使城市建成區(qū)面積和城市經(jīng)濟(jì)社會(huì)指標(biāo)突然變化，這必然使傳統(tǒng)的基于線性假設(shè)的回歸預(yù)測(cè)模型產(chǎn)生很大的誤差。

4 小結(jié)

城市建成區(qū)面積是人口、經(jīng)濟(jì)、社會(huì)、環(huán)境等多因素綜合影響的結(jié)果，從本研究中發(fā)現(xiàn)，2015年所選25個(gè)省會(huì)城市的建成區(qū)面積受第三產(chǎn)業(yè)和國(guó)民生產(chǎn)總值的影響最大，第二產(chǎn)業(yè)和工業(yè)的影響次之。說明要想對(duì)城市進(jìn)行合理建設(shè)使之健康發(fā)展，應(yīng)該以大力調(diào)控第三產(chǎn)業(yè)、國(guó)民生產(chǎn)總值為出發(fā)點(diǎn)。

本文利用建成區(qū)面積作為研究變量，運(yùn)用Matlab軟件建立多元線性回歸模型，利用最小二乘原理求解數(shù)學(xué)模型中的最優(yōu)解，不僅簡(jiǎn)化和優(yōu)化了繁瑣的計(jì)算過程，而且通過檢驗(yàn)證明所建立的模型計(jì)算結(jié)果精度較高，對(duì)城市建成區(qū)面積變化的預(yù)測(cè)有一定的參考價(jià)值。然而，由于在實(shí)際情況中城市建成區(qū)面積受復(fù)雜的多種因素影響，在建模過程中對(duì)因素考慮或選擇不同，則會(huì)造成計(jì)算結(jié)果的多樣性，因此，在研究中應(yīng)力求完善。

參考文獻(xiàn)：

[1] 李愛民.基于遙感影像的城市建成區(qū)擴(kuò)張與用地規(guī)模研究[D].鄭州：解放軍信息工程大學(xué)，2009.

[2] 劉柯.基于主成分分析的BP神經(jīng)網(wǎng)絡(luò)在城市建成區(qū)面積預(yù)測(cè)中的應(yīng)用——以北京市為例[J].地理科學(xué)進(jìn)展，2007(6)：129-137.

[3] 雷波.BP神經(jīng)網(wǎng)絡(luò)和多元回歸模型在城市建成區(qū)面積預(yù)測(cè)中的應(yīng)用比較——以福州市為例[J].城市發(fā)展研究，2008(1)：153-155.

[4] 周永生，肖玉歡，黃潤(rùn)生.基于多元線性回歸的廣西糧食產(chǎn)量預(yù)測(cè)[J].南方農(nóng)業(yè)學(xué)報(bào)，2011，42(9)：1165-1167.

[5] 付倩嬈.基于多元線性回歸的霧霾預(yù)測(cè)方法研究[J].計(jì)算機(jī)科學(xué)，2016，43(S1)：526-528.

[6] 葉鋒.多元線性回歸在經(jīng)濟(jì)技術(shù)產(chǎn)量預(yù)測(cè)中的應(yīng)用[J].中外能源，2015，20(2)：45-48.

[7] 韋浩.多元回歸分析法在滑坡空間預(yù)測(cè)中的應(yīng)用[D].西安：長(zhǎng)安大學(xué)，2011.

[8] 張智星.MATLAB程序設(shè)計(jì)與應(yīng)用[M].北京：清華大學(xué)出版社，1993：56-65.

[9] 邱衛(wèi)寧，陶本藻.測(cè)量數(shù)據(jù)處理理論與方法[M].武漢：武漢大學(xué)出版社，2008：109-113.

[10] 劉大杰，陶本藻.實(shí)用測(cè)量數(shù)據(jù)處理方法[M].北京：測(cè)繪出版社，2000：7-14.

[11] 王樂洋，朱建軍.回歸分析、測(cè)量平差與大地測(cè)量反演[J].測(cè)繪通報(bào)，2007(2)：27-30.