基于Matlab矩陣運(yùn)算的多元回歸函數(shù)編程

摘 要：多元線(xiàn)性回歸是計(jì)量分析中最常用的方法之一，在經(jīng)濟(jì)學(xué)、金融學(xué)和管理學(xué)中運(yùn)用廣泛。線(xiàn)性回歸方法依賴(lài)于矩陣和向量的運(yùn)算，因此以矩陣運(yùn)算為內(nèi)核的Matlab是理想的編程平臺(tái)。本文以開(kāi)發(fā)多元線(xiàn)性回歸函數(shù)作為例子，著重分析Matlab自定義函數(shù)的設(shè)計(jì)理念，同時(shí)突出varargin、nargin等特殊函數(shù)的代碼編寫(xiě)。通過(guò)交叉運(yùn)用線(xiàn)性代數(shù)、概率統(tǒng)計(jì)等數(shù)學(xué)工具，本文為計(jì)量研究人員提供了新的編程視角。

關(guān)鍵詞：矩陣運(yùn)算；Matlab；多元線(xiàn)性回歸；函數(shù)編程

中圖分類(lèi)號(hào)：TP311.1 文獻(xiàn)標(biāo)識(shí)碼：A

1 引言（Introduction）

普通最小二乘法（Ordinary Least Square，以下簡(jiǎn)稱(chēng)OLS）是最重要的多元線(xiàn)性回歸法，該方法被廣泛地應(yīng)用于各種定量分析領(lǐng)域[1]。目前，許多主流計(jì)量分析軟件都支持OLS回歸法，比如SAS、R、Eviews等。Matlab也能夠?qū)崿F(xiàn)多元線(xiàn)性回歸，其自帶的Statistics Toolbox和Econometrics Toolbox提供了比較豐富的計(jì)量分析功能。既然有現(xiàn)成的回歸分析軟件，為什么還要自己編寫(xiě)原始代碼？原因有以下兩點(diǎn)。第一，計(jì)量經(jīng)濟(jì)學(xué)的發(fā)展是日新月異的，上述主流回歸分析軟件的更新無(wú)法跟上計(jì)量經(jīng)濟(jì)學(xué)的前進(jìn)步伐。通過(guò)對(duì)原始代碼的詳細(xì)分解，能夠讓廣大研究人員掌握回歸分析的核心編程思路。第二，計(jì)量經(jīng)濟(jì)學(xué)的難點(diǎn)在于線(xiàn)性代數(shù)和概率統(tǒng)計(jì)，在實(shí)踐中學(xué)習(xí)是掌握各種計(jì)量工具的關(guān)鍵。通過(guò)編程，有利于更好地吸收回歸分析知識(shí)。

多元線(xiàn)性回歸的代碼設(shè)計(jì)包括數(shù)據(jù)儲(chǔ)存、矩陣運(yùn)算、結(jié)果輸出等環(huán)節(jié)。本文先從Matlab自定義函數(shù)編程出發(fā)，利用varargin、nargin等特殊函數(shù)處理變量數(shù)據(jù)。之后利用Matlab強(qiáng)大的矩陣運(yùn)算功能得到回歸分析統(tǒng)計(jì)量，并將回歸結(jié)果以表格形式輸出。目前，國(guó)內(nèi)有少量的文獻(xiàn)涉及Matlab計(jì)量編程，這些研究存在一些可改進(jìn)之處。首先，回歸方程的變量數(shù)目應(yīng)該由用戶(hù)自定義。在計(jì)量經(jīng)濟(jì)學(xué)的實(shí)際應(yīng)用中，控制變量的篩選是不確定的。借助Matlab內(nèi)置的varargin、nargin等函數(shù)，本文的程序能夠有效地處理多個(gè)解釋變量，即可以實(shí)現(xiàn)變量數(shù)目的自定義。其次，回歸結(jié)果應(yīng)該給出完整的統(tǒng)計(jì)量，不應(yīng)該僅僅提供系數(shù)大小。本文的Matlab程序給出了回歸系數(shù)、標(biāo)準(zhǔn)誤、P值等統(tǒng)計(jì)量，符合計(jì)量分析的學(xué)術(shù)標(biāo)準(zhǔn)。再次，回歸結(jié)果可以用表格形式輸出。Matlab的內(nèi)置table函數(shù)提供了豐富的輸出選項(xiàng)，在多元回歸的情況下可以用table函數(shù)批量地輸出各種統(tǒng)計(jì)量。

2 線(xiàn)性回歸法的理論回顧（Review of linear regression）

用表示被解釋變量，表示自變量，表示干擾項(xiàng)，那么總體線(xiàn)性方程可以表示為[2]：

用表示被解釋變量的向量形式，表示自變量的矩陣形式，表示樣本容量，我們可以得到回歸系數(shù)：

在同方差條件下，顯著性檢驗(yàn)所需的標(biāo)準(zhǔn)誤可以由下式計(jì)算得到[3]：

3 常用的內(nèi)置函數(shù)（Essential built-in functions）

為了構(gòu)建完整的多元計(jì)量分析程序，我們將代碼分為多個(gè)模塊。一個(gè)完整的計(jì)量分析程序可以被分為三個(gè)模塊，數(shù)理統(tǒng)計(jì)模塊、矩陣運(yùn)算模塊和結(jié)果輸出模塊。數(shù)理統(tǒng)計(jì)模塊是回歸參數(shù)檢驗(yàn)的重要部分，我們需要借助Matlab內(nèi)置的數(shù)學(xué)函數(shù)來(lái)計(jì)算各種統(tǒng)計(jì)量。矩陣運(yùn)算是Matlab的強(qiáng)項(xiàng)[3]，在求回歸參數(shù)時(shí)需要用到diag函數(shù)。輸出模塊也非常關(guān)鍵，內(nèi)置的table函數(shù)能夠使回歸結(jié)果以簡(jiǎn)潔的表格形式輸出。在正式編寫(xiě)回歸分析程序之前，我們先回顧一下將要涉及的Matlab內(nèi)置函數(shù)。

3.1 用于顯著性檢驗(yàn)的tcdf函數(shù)

在計(jì)量分析中t分布具有重要作用，t分布可以用來(lái)檢驗(yàn)?zāi)硞€(gè)解釋變量的顯著性。為了得到P值，我們需要計(jì)算t統(tǒng)計(jì)量。假定某個(gè)自變量的系數(shù)為，那么該變量的t統(tǒng)計(jì)量計(jì)算公式為，其中表示解釋變量的標(biāo)準(zhǔn)誤。利用Matlab內(nèi)置的tcdf函數(shù)，我們可以得到任意自由度的t統(tǒng)計(jì)量。比如，某個(gè)變量的t值為-2.39，自由度為296，那么求相應(yīng)P值的代碼為（這里還用到了求絕對(duì)值的abs函數(shù)）：

pValue=（1-tcdf（abs（-2.39），296））*2

3.2 矩陣運(yùn)算與diag函數(shù)

多元回歸方法涉及大量的矩陣運(yùn)算，因此我們可以借助Matlab進(jìn)行與矩陣相關(guān)的操作。比如，先用ones內(nèi)置函數(shù)構(gòu)造一個(gè)元素全為1的5×5矩陣，將該矩陣命名為mat1。再將該矩陣與其轉(zhuǎn)置矩陣相乘，將得到的新矩陣命名為mat2。利用Matlab內(nèi)置的diag函數(shù)，我們可以很快地得到矩陣mat2主對(duì)角線(xiàn)上的所有元素。上述過(guò)程的代碼如下：

mat1=ones（5，5）

mat2=mat1*mat1'

diag（mat1）

3.3 用于生成表格的table函數(shù)

回顧國(guó)內(nèi)已有的Matlab學(xué)術(shù)文獻(xiàn)，一些學(xué)者采用文本對(duì)話(huà)框的形式輸出計(jì)量回歸結(jié)果。但是，在回歸統(tǒng)計(jì)量數(shù)目極多的情況下，文本對(duì)話(huà)框無(wú)法像表格那樣簡(jiǎn)潔地顯示層次化信息。Matlab提供了強(qiáng)大的表格輸出命令，通過(guò)內(nèi)置的table函數(shù)我們可以很輕松地得到專(zhuān)業(yè)化的表格效果。假設(shè)有兩個(gè)變量var1和var1，其中var1的取值為21和11，var2的取值為35和51。如果用table函數(shù)來(lái)實(shí)現(xiàn)var1和var2數(shù)據(jù)的表格形式輸出，可以編寫(xiě)以下代碼：

table（[21；11]，[35；51]，...

'VariableNames'，{'var1' 'var2'}）

4 Matlab代碼設(shè)計(jì)（Matlab code design）

4.1 自定義函數(shù)的開(kāi)頭

將這個(gè)Matlab多元線(xiàn)性回歸程序命名為mvReg，因此將函數(shù)源代碼以mvReg.m作為文件名保存。為了使該程序能夠靈活處理多元線(xiàn)性回歸的自變量數(shù)目，我們使用varargin函數(shù)記錄輸入的因變量和自變量。代碼如下：

function mvReg（varargin）

4.2 回歸系數(shù)的估計(jì)

由第二部分理論回歸可知，回歸系數(shù)的表達(dá)式可以簡(jiǎn)化為。首先，生成元素全為1的矩陣MatX。其次，用循環(huán)語(yǔ)句將矩陣MatX的元素替換為varargin中的數(shù)據(jù)。求回歸系數(shù)的代碼如下：

NumIV=nargin；

NumCross=length（varargin{1}）；

VecY=varargin{1}；

MatX=ones（NumCross，NumIV）；

for i=2：NumIV

MatX（：，i）=varargin{i}；

end

VecB=inv（MatX'*MatX）*MatX'*VecY；

4.3 顯著性檢驗(yàn)的統(tǒng)計(jì)量

在回歸分析當(dāng)中，比較重要的統(tǒng)計(jì)量有t值、標(biāo)準(zhǔn)誤和P值等。t值等于回歸系數(shù)與標(biāo)準(zhǔn)的比值，表達(dá)式為。P值根據(jù)t分布得到，極小的P值表明不拒絕原假設(shè)“”。用VecSE、VecT和VecP分別表示標(biāo)準(zhǔn)誤向量、t值向量和P值向量，計(jì)算各統(tǒng)計(jì)量的代碼如下：

ehat=VecY-MatX*VecB；

VecSE=[diag（sum（ehat.^2）/…

（NumCross-NumIV）*inv（MatX'…

*MatX））].^0.5；

VecT=VecB./VecSE；

VecP=（1-tcdf（abs（VecT），…

NumCross-NumIV））*2；

4.4 將回歸結(jié)果以表格形式輸出

Matlab的內(nèi)置函數(shù)table提供了強(qiáng)大的表格顯示功能，因此本回歸程序利用table函數(shù)來(lái)實(shí)現(xiàn)回歸結(jié)果的輸出。由于回歸結(jié)果的表格涉及文本內(nèi)容和數(shù)值內(nèi)容，我們將數(shù)據(jù)以cell形式儲(chǔ)存。代碼如下：

VecVN=cell（nargin，1）；

VecVN{1，1}='_cons'；

for i=2：nargin

VecVN{i，1}=inputname（i）；

end

Regression=table（VecVN，VecB，…

VecSE，VecP，'VariableNames'，…

{'Variable'，'Coef'，'Std_Err'，'P'}）

至此代碼編寫(xiě)完畢，接下來(lái)對(duì)多元回歸函數(shù)進(jìn)行測(cè)試。先導(dǎo)入Excel數(shù)據(jù)文件，因變量為wage，自變量為educ和exper。多元回歸的代碼為：

mvReg（wage，educ，exper）

最終的回歸結(jié)果如圖1所示。輸出結(jié)果的第一列為被估計(jì)的回歸系數(shù)，解釋變量educ的系數(shù)大小為0.64427，標(biāo)準(zhǔn)誤為0.053。通過(guò)計(jì)算t統(tǒng)計(jì)量，我們可以得到其P值為0。所以，解釋變量educ是非常統(tǒng)計(jì)性顯著的。同理，解釋變量exper的系數(shù)大小為0.070095，標(biāo)準(zhǔn)誤為0.010978。該變量P值小數(shù)點(diǎn)的后三位都為零，因此變量exper在1%水平上統(tǒng)計(jì)性顯著。參照主流統(tǒng)計(jì)軟件的輸出格式，我們用“_cons”表示截距。從圖1的回歸結(jié)果來(lái)看，截距P值小數(shù)點(diǎn)的后三位都為零，因此也是在1%水平上統(tǒng)計(jì)性顯著。

圖1 Matlab程序運(yùn)行結(jié)果

Fig.1 Matlab' s output

在得到多元線(xiàn)性回歸的結(jié)果之后，我們還可以利用內(nèi)置函數(shù)來(lái)實(shí)現(xiàn)數(shù)據(jù)可視化的功能。在單個(gè)解釋變量的情況下，scatter函數(shù)可以得到關(guān)于因變量與被解釋變量的散點(diǎn)圖。借助plot函數(shù)，還可以畫(huà)出線(xiàn)性回歸的擬合曲線(xiàn)。借助基本的數(shù)理統(tǒng)計(jì)知識(shí)，用戶(hù)可以將這些圖形功能充分應(yīng)用到數(shù)據(jù)分析當(dāng)中。值得注意的是，Matlab提供了大量的高級(jí)統(tǒng)計(jì)函數(shù)?；谶@些內(nèi)置函數(shù)，用戶(hù)可以非常高效地編寫(xiě)屬于自己的新程序。

5 結(jié)論（Conclusion）

Matlab作為以矩陣運(yùn)算為特點(diǎn)的綜合性編程平臺(tái)，其高效快捷的代碼命令能夠處理眾多的計(jì)量分析難題。本文借助Matlab線(xiàn)性代數(shù)模塊、概率分布模塊實(shí)現(xiàn)了多元線(xiàn)性回歸法，并得到了重要的檢驗(yàn)統(tǒng)計(jì)量。在此基礎(chǔ)上，可以進(jìn)一步增加更高級(jí)的計(jì)量分析功能，比如異方差穩(wěn)健回歸、多重共線(xiàn)性檢驗(yàn)等。

參考文獻(xiàn)（References）

[1] Eric Noreen.An Empirical Comparison of Probit and OLS

Regression Hypothesis Tests. Journal Of Accounting Research，

1988，26（1）：119-133.

[2] Sanford Weisberg. Applied Linear Regression [M].John Wiley

& Sons，2013.

[3] David Thomas.Give Meaning to Matrices with MATLAB[J].

Journal of Computers in Mathematics and Science Teaching，

1990，9（3）：73-85.

作者簡(jiǎn)介：

鄔汝威（1988-），男，碩士生.研究領(lǐng)域：算法設(shè)計(jì)，面板數(shù)據(jù)

分析.