R語言與機器學習中的迴歸方法
機器學習中的一些方法如決策樹，隨機森林，SVM，神經網路由於對資料沒有分佈的假定等普通線性迴歸模型的一些約束，預測效果也比較不錯，交叉驗證結果也能被接受。下面以R中lars包包含資料集diabetes為例說明機器學習中的迴歸方法。
一、資料集及交叉驗證辦法描述
Diabetes資料集包含在R的lars包中，資料分為x，y，x2三個部分，因變數為y，資料是關於糖尿病的血液化驗等指標。這個資料集最早被用在偏最小二乘迴歸的文章裡。
交叉驗證採用指標NMSE來評價模型好壞。這一統計量是計算模型預測效能和基準模型的預測效能之間的比率。通常採用目標變數的平均值來作為基準模型。其取值範圍通常為0～1。如果模型表現優於這個非常簡單的基準模型預測，那麼NMSE應明顯小於1。NMSE的值越小，模型的效能就越好。NMSE的值大於1，意味著模型預測還不如簡單地把所有個案的平均值作為預測值！
交叉驗證辦法為將資料集分為5份，取4份作為訓練集，1份作為測試集，共作5次，把誤差平均起來作為衡量標準。選取程式碼如下：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
n<-length（dataset）
index1<-1：n
index2<-rep（1：5，ceiling（n/5））［1：n］
index2<-sample（index2，n）
二、迴歸樹
決策樹是透過一系列規則對資料進行分類的過程。它提供一種在什麼條件下會得到什麼值的類似規則的方法。決策樹分為分類樹和迴歸樹兩種，分類樹對離散變數做決策樹，迴歸樹對連續變數做決策樹。
基本演算法：
1。選擇一個屬性放置在根節點，為每個可能的屬性值產生一個分支
2。將樣本劃分成多個子集，一個子集對應於一個分支
3。在每個分支上遞迴地重複這個過程，僅使用真正到達這個分支的樣本
4。如果在一個節點上的所有樣本擁有相同的類別，即停止該部分樹的擴充套件
構造決策樹（集合劃分）時選擇屬性：
1。ID3：Information Gain
2。C4。5：Gain Ratio
3。CART：Gini Index
在R中我們使用rpart包中的rpart（）函式實現樹迴歸。我們先把rpart包中的兩個十分重要的函式介紹如下：
構建迴歸樹的函式：rpart（）用法如下：
rpart（formula， data， weights， subset，na。action = na。rpart， method，
model = FALSE， x = FALSE， y = TRUE， parms， control， cost， 。。。）
主要引數說明：
fomula迴歸方程形式：例如 y～x1+x2+x3。
data資料：包含前面方程中變數的資料框（dataframe）。
na。action缺失資料的處理辦法：預設辦法是刪除因變數缺失的觀測而保留自變數缺失的觀測。
method根據樹末端的資料型別選擇相應變數分割方法，本引數有四種取值：連續型“anova”；離散型“class”；計數型（泊松過程）“poisson”；生存分析型“exp”。程式會根據因變數的型別自動選擇方法，但一般情況下最好還是指明本引數，以便讓程式清楚做哪一種樹模型。
parms用來設定三個引數：先驗機率、損失矩陣、分類純度的度量方法。
control控制每個節點上的最小樣本量、交叉驗證的次數、複雜性參量：即cp：complexitypamemeter，這個引數意味著對每一步拆分，模型的擬合優度必須提高的程度，等等。
進行剪枝的函式：prune（）用法如下：
prune（tree， cp， 。。。）
主要引數說明：
tree一個迴歸樹物件，常是rpart（）的結果物件。
cp複雜性參量，指定剪枝採用的閾值。cp全稱為complexity parameter，指某個點的複雜度，對每一步拆分，模型的擬合優度必須提高的程度，用來節省剪枝浪費的不必要的時間，R內部是怎麼計算的我不知道，希望讀者能夠補充。
進一步補充一點關於CP值的東西：建立樹模型要權衡兩方面問題，一個是要擬合得使分組後的變異較小，另一個是要防止過度擬合，而使模型的誤差過大，前者的引數是CP，後者的引數是Xerror。所以要在Xerror最小的情況下，也使CP儘量小。如果認為樹模型過於複雜，我們需要對其進行修剪 。（摘自推酷上的《分類-迴歸樹模型（CART）在R語言中的實現》）
執行程式碼：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
library（rpart。plot）
reg<-rpart（y~。，w）
rpart。plot（reg，type=2，faclen=T）
得到迴歸樹：
我們可以透過print（reg）來看到樹的各個節點的細節。
我們來進行交叉驗證，執行程式碼如下：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
w<-diabetes［，2：3］
n<-length（w$y）
index1<-1：n
index2<-rep（1：5，ceiling（n/5））［1：n］
index2<-sample（index2，n）
NMSE<-rep（0，5）
NMSE0<-NMSE
for（i in 1：5）{
m<-index1［index2==i］
reg<-rpart（y~。，w［-m，］）
y0<-predict（reg，w［-m，］）
y1<-predict（reg，w［m，］）
NMSE0［i］<-mean（（w$y［-m］-y0）^2）/mean（（w$y［-m］-mean（w$y［-m］））^2）
NMSE［i］<-mean（（w$y［m］-y1）^2）/mean（（w$y［m］-mean（w$y［m］））^2）
}
R中輸出結果：
> NMSE
［1］ 0。7892592 0。8857756 0。8619379 1。0072968 0。7238316
> NMSE0
［1］ 0。3503767 0。3417909 0。3400387 0。3192845 0。3467186
明顯出現了過擬合現象，應該使用剪枝函式，對模型進行修正。
> reg$cptable
CP n split relerror xerror xstd
1 0。29154165 0 1。0000000 1。0040015 0。05033316
2 0。08785891 1 0。7084583 0。8040962 0。04896667
3 0。05660089 2 0。6205994 0。7227529 0。04657623
4 0。02169615 3 0。5639986 0。6424826 0。04302580
5 0。02093950 4 0。5423024 0。6591446 0。04376777
6 0。01723601 50。5213629 0。6749867 0。04578783
7 0。01678503 6 0。5041269 0。6841483 0。04554068
8 0。01349365 8 0。4705568 0。6845580 0。04429950
9 0。01166564 9 0。4570632 0。7370893 0。04829371
10 0。01089168 11 0。43373190。7452419 0。05041336
11 0。01070564 12 0。42284020。7417955 0。05054043
12 0。01042308 14 0。40142890。7399988 0。05088835
13 0。01003363 15 0。39100580。7566972 0。05143535
14 0。01000000 17 0。37093860。7650728 0。05110011
參照上述結果，選擇合適的cp值。故修正為：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
reg2<-prune（reg，cp=0。025）
rpart。plot（reg2，type=2，faclen=T）
結果為：
再次進行交叉驗證（程式碼略）可以看到：
> NMSE
［1］ 0。5982049 0。6995054 0。6826815 0。8970573 0。6407927
> NMSE0
［1］ 0。5559462 0。5177565 0。4953384 0。5019682 0。5233709
過擬合現象基本消除。
三、boosting迴歸
Boosting方法是一種用來提高弱分類演算法準確度的方法，這種方法透過構造一個預測函式系列，然後以一定的方式將他們組合成一個預測函式。他是一種框架演算法，主要是透過對樣本集的操作獲得樣本子集，然後用弱分類演算法在樣本子集上訓練生成一系列的基分類器。他可以用來提高其他弱分類演算法的識別率，也就是將其他的弱分類演算法作為基分類演算法放於Boosting 框架中，透過Boosting框架對訓練樣本集的操作，得到不同的訓練樣本子集，用該樣本子集去訓練生成基分類器；每得到一個樣本集就用該基分類演算法在該樣本集上產生一個基分類器，這樣在給定訓練輪數 n 後，就可產生 n 個基分類器，然後Boosting框架演算法將這 n個基分類器進行加權融合，產生一個最後的結果分類器，在這 n個基分類器中，每個單個的分類器的識別率不一定很高，但他們聯合後的結果有很高的識別率，這樣便提高了該弱分類演算法的識別率。
Boosting方法簡而言之就是採取少數服從多數的原理，他的合理性在於如果每個迴歸樹的出錯率為40%，那麼1000棵決策樹犯錯誤的機率就降到了4。40753e-11，這就比較小了。
對diabetes資料做boosting迴歸，使用到的函式包為mboost，使用函式為mboost。用法如下：
mboost（formula， data = list（），
baselearner = c（“bbs”， “bols”， “btree”， “bss”， “bns”）， 。。。）
其中formular需要使用到函式btree（）：
btree（。。。， tree_controls = ctree_control（stump = TRUE，
mincriterion = 0，
savesplitstats = FALSE））
試執行下列程式碼：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
library（mboost）
reg<-mboost（y~btree（x2。ltg）+btree（x2。bmi）+btree（x2。bmi。ltg），
data=w，control=boost_control（nu=0。1））
我們可以得到boosting迴歸的資訊，透過plot函式還有視覺化結果。
這裡值得一提的是我在做boosting迴歸時為了減少工作量直接做了變數選擇（這個選擇與迴歸樹的最終選擇是一致的，所以也有一定的道理）。最後交叉驗證的結果為：測試集NMSE為0。5513152，訓練集NMSE為0。4656569。比起上次計算的迴歸樹模型測試集NMSE為0。7036484，訓練集NMSE為0。518876好了不少。
四、bagging迴歸
與boosting迴歸想法類似，bagging迴歸的做法就是不斷放回地對訓練樣本進行再抽樣，對每個自助樣本都建立一棵迴歸樹，對於每一個觀測，每棵樹給一個預測，最後將其平均。
對diabetes資料做bagging迴歸，使用到的函式包為ipred，使用函式為bagging（），用法如下：
bagging（formula， data， subset， na。action=na。rpart， 。。。）
主要引數介紹：
Formula：迴歸方程形式
Data：資料集
Control：對樹枝的控制，使用函式rpart。control（），可以控制諸如cp值，xval等參量。
輸入程式碼：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
library（ipred）
reg<-bagging（y~。，data=w，coob=TRUE，control=rpart。control（cp=0。025））
結果為：
Baggingregression trees with 25 bootstrap replications
Call：bagging。data。frame（formula = y ~ 。， data = w， coob = TRUE， control =rpart。control（cp = 0。025））
Out-of-bagestimate of root mean squared error： 58。3648
使用交叉驗證（程式碼略），得到結果：測試集NMSE為0。5705753，訓練集NMSE為0。3906232。比起上次計算的迴歸樹模型測試集NMSE為0。7036484，訓練集NMSE為0。518876好了不少。
五、隨機森林迴歸
與bagging迴歸相比，隨機森林則更進一步，不僅對樣本進行抽樣，還對變數進行抽樣。
對diabetes資料做隨機森林迴歸，使用到的函式包為randomForest，使用函式為randomForest（），用法如下：
randomForest（formula， data=NULL， 。。。， subset， na。action=na。fail）
這裡值得一提的是，隨機森林有個十分牛逼的性質，不會出現過擬合現象，這也省去了我們確定樹的分類程度等一系列麻煩的事情。得到結果：測試集NMSE為0。08992529，訓練集NMSE為0。08835731，效果顯著提升。隨機森林還可以輸出自變數重要性度量，試執行程式碼：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
library（randomForest）
reg<-randomForest（y~。，data=w，importance=TRUE）
reg$importance
得到結果：
%IncMSE IncNodePurity
x。glu 68。8034199 42207。351
x2。age 22。6784331 18569。370
x2。sex 6。2735713 2808。346
x2。bmi 1379。0675134 371372。494
x2。map 331。3925059 113411。547
x2。tc 18。6080948 14990。179
x2。ldl 24。3690847 17457。214
x2。hdl 216。2741620 64627。209
x2。tch 419。0451857 93688。855
x2。ltg 1514。0912885 379235。430
x2。glu 81。7568020 51984。121
x2。age。2 1。5242836 19364。582
x2。bmi。2 75。6345112 53635。024
x2。map。2 5。9156799 23049。475
x2。tc。2 1。6792683 15631。426
（省略部分輸出）
其中第二列為均方誤差遞減意義下的重要性，第三列為精確度遞減意義下的重要性。
六、其他的迴歸方法
除去我們所說的以上4種方法外，還有人工神經網路迴歸，SVM迴歸，他們可以透過nnet，rminer包中有關函式實現，這裡我們從略。
在結束本文之前，我想我們可以做一件更有意義的事。在介紹diabetes資料集時，我們提到了這個資料最早是一個關於偏最小二乘的例子。那麼想必他與這個模型的契合度還是不錯的，我們可以嘗試去算算利用偏最小二乘得到的訓練集與測試集的NMSE。
程式碼如下：
［plain］ view plain copy print？在CODE上檢視程式碼片派生到我的程式碼片
library（lars）
library（pls）
NMSE<-rep（0，5）
NMSE0<-NMSE
for（i in 1：5）{
m<-index1［index2==i］
reg。pls <- plsr（y ~ x2， 64， data = diabetes［-m，］， validation =“CV”）
y0<-predict（reg。pls，diabetes［-m，］）
y1<-predict（reg。pls，diabetes［m，］）
NMSE0［i］<-mean（（diabetes$y［-m］-y0）^2）/mean（（diabetes$y［-m］-mean（diabetes$y［-m］））^2）
NMSE［i］<-mean（（diabetes$y［m］-y1）^2）/mean（（diabetes$y［m］-mean（diabetes$y［m］））^2）
}
執行結果：測試集NMSE為0。6094071，訓練集NMSE為0。4031211。這裡也僅僅是得到了一個優於迴歸樹，與bagging迴歸相當，不如隨機森林的一個結果。也說明了機器學習方法在實際問題中還是有一定優勢的。但這絕不意味著這樣的迴歸可以代替OLS，畢竟這些資料探勘的方法在表示式的簡潔程度，計算量等諸多方面還是存在著或多或少的不足。