1. 簡介
隨機(jī)森林是機(jī)器學(xué)習(xí)中的一種分類算法，這個(gè)算法在之前分享的TCGA生存模型的構(gòu)建以及模型預(yù)測和評(píng)估中使用過孽查。在介紹隨機(jī)森林之前鳄炉，非常有必要了解決策樹這種分類器囱淋。
決策樹是一種分類器嗜闻，通過訓(xùn)練集構(gòu)建一顆決策樹斩芭，從而可以對(duì)新的數(shù)據(jù)預(yù)測其分類乏冀。一顆構(gòu)建好的決策樹如下：
可以看到這顆決策樹的目標(biāo)是將數(shù)據(jù)分成 "接受" 和 "拒絕" 兩類蝶糯，分類的條件有樹中的節(jié)點(diǎn)來決定；
而隨機(jī)森林算法辆沦，可以看作好多顆決策樹構(gòu)成的分類器昼捍，首先通過有放回的抽樣從原始數(shù)據(jù)集中構(gòu)建多個(gè)子數(shù)據(jù)集识虚，然后利用每個(gè)子數(shù)據(jù)集構(gòu)建一顆決策樹，最終的分類效果由多顆決策樹預(yù)測得到的眾數(shù)決定妒茬；
之所以叫做隨機(jī)森林舷礼，是因?yàn)閮蓚€(gè)核心觀點(diǎn)：
1）子數(shù)據(jù)集的構(gòu)建，通過隨機(jī)抽樣得到郊闯，所以有隨機(jī)這個(gè)關(guān)鍵詞
2）在這個(gè)分類器中妻献，有多顆決策樹，所以有森林這個(gè)關(guān)鍵詞
隨機(jī)森林的算法可以用如下幾個(gè)步驟概括：
1）用有抽樣放回的方法（bootstrap）從樣本集中選取n個(gè)樣本作為一個(gè)訓(xùn)練集
2）用抽樣得到的樣本集生成一棵決策樹团赁。在生成的每一個(gè)結(jié)點(diǎn)：
隨機(jī)不重復(fù)地選擇d個(gè)特征
3）利用這d個(gè)特征分別對(duì)樣本集進(jìn)行劃分育拨，找到最佳的劃分特征（可用基尼系數(shù)、增益率或者信息增益判別）
4）重復(fù)步驟1到步驟2共k次欢摄，k即為隨機(jī)森林中決策樹的個(gè)數(shù)熬丧。
用訓(xùn)練得到的隨機(jī)森林對(duì)測試樣本進(jìn)行預(yù)測，并用票選法決定預(yù)測的結(jié)果怀挠。
2. 隨機(jī)森林之特征選擇
隨機(jī)森林具有一個(gè)重要特征：能夠計(jì)算單個(gè)特征變量的重要性析蝴。并且這一特征在很多方面能夠得到應(yīng)用。例如在銀行貸款業(yè)務(wù)中能否正確的評(píng)估一個(gè)企業(yè)的信用度绿淋，關(guān)系到是否能夠有效地回收貸款闷畸。但是信用評(píng)估模型的數(shù)據(jù)特征有很多，其中不乏有很多噪音吞滞，所以需要計(jì)算出每一個(gè)特征的重要性并對(duì)這些特征進(jìn)行一個(gè)排序佑菩，進(jìn)而可以從所有特征中選擇出重要性靠前的特征。
- 一：特征重要性
在隨機(jī)森林中某個(gè)特征X的重要性的計(jì)算方法如下：
a：對(duì)于隨機(jī)森林中的每一顆決策樹,使用相應(yīng)的OOB(袋外數(shù)據(jù))數(shù)據(jù)來計(jì)算它的袋外數(shù)據(jù)誤差,記為errOOB1.
b: 隨機(jī)地對(duì)袋外數(shù)據(jù)OOB所有樣本的特征X加入噪聲干擾(就可以隨機(jī)的改變樣本在特征X處的值),再次計(jì)算它的袋外數(shù)據(jù)誤差,記為errOOB2.
c：假設(shè)隨機(jī)森林中有Ntree棵樹,那么對(duì)于特征X的重要性=∑(errOOB2-errOOB1)/Ntree,之所以可以用這個(gè)表達(dá)式來作為相應(yīng)特征的重要性的度量值是因?yàn)椋喝艚o某個(gè)特征隨機(jī)加入噪聲之后,袋外的準(zhǔn)確率大幅度降低,則說明這個(gè)特征對(duì)于樣本的分類結(jié)果影響很大,也就是說它的重要程度比較高裁赠。
- 二：特征選擇
在論文 Variable Selection using Random Forests中詳細(xì)的論述了基于隨機(jī)森林的特征選擇方法殿漠，這里我們進(jìn)行一些回顧。
首先特征選擇的目標(biāo)有兩個(gè)：
1）找到與應(yīng)變量高度相關(guān)的特征變量佩捞。
2）選擇出數(shù)目較少的特征變量并且能夠充分的預(yù)測應(yīng)變量的結(jié)果绞幌。
其次一般特征選擇的步驟為：
1）初步估計(jì)和排序
a: 對(duì)隨機(jī)森林中的特征變量按照VI（Variable Importance）降序排序。
b: 確定刪除比例,從當(dāng)前的特征變量中剔除相應(yīng)比例不重要的指標(biāo)一忱，從而得到一個(gè)新的特征集莲蜘。
c: 用新的特征集建立新的隨機(jī)森林,并計(jì)算特征集中每個(gè)特征的VI,并排序。
d: 重復(fù)以上步驟,直到剩下m個(gè)特征掀潮。
2）根據(jù)1中得到的每個(gè)特征集和它們建立起來的隨機(jī)森林,計(jì)算對(duì)應(yīng)的袋外誤差率(OOB err),將袋外誤差率最低的特征集作為最后選定的特征集菇夸。
3. randomForest 包的使用
randomForest 包提供了利用隨機(jī)森林算法解決分類和回歸問題的功能；我們這里只關(guān)注隨機(jī)森林算法在分類問題中的應(yīng)用
安裝randomForest包
install.packages("randomForest")
library(randomForest)
載入iris演示數(shù)據(jù)集
iris 數(shù)據(jù)集共有150行仪吧，5列庄新，其中第5列Species為分類變量，共有3種分類情況。
這個(gè)數(shù)據(jù)集可以看做有150個(gè)樣本择诈，將Species這個(gè)分類變量作為因變量械蹋，使用其他4個(gè)指標(biāo)來預(yù)測花屬于哪個(gè)Species。
data(iris)
str(iris)
# 'data.frame': 150 obs. of 5 variables:
# $ Sepal.Length: num 5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
# $ Sepal.Width : num 3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
# $ Petal.Length: num 1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
# $ Petal.Width : num 0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
# $ Species : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...
head(iris)
# Sepal.Length Sepal.Width Petal.Length Petal.Width Species
# 1 5.1 3.5 1.4 0.2 setosa
# 2 4.9 3.0 1.4 0.2 setosa
# 3 4.7 3.2 1.3 0.2 setosa
# 4 4.6 3.1 1.5 0.2 setosa
# 5 5.0 3.6 1.4 0.2 setosa
# 6 5.4 3.9 1.7 0.4 setosa
接下來調(diào)用randomForest 函數(shù)進(jìn)行分類
調(diào)用該函數(shù)時(shí)羞芍，通過一個(gè)表達(dá)式指定分類變量 Species 和對(duì)應(yīng)的數(shù)據(jù)集data 就可以了哗戈。Species ~ .的意思是iris數(shù)據(jù)集中以Species為被解釋變量，其他列為解釋變量荷科。后面的importance 和 proximity 是計(jì)算每個(gè)變量的重要性和樣本之間的距離
set.seed(71)
iris.rf <- randomForest(Species ~ ., data=iris, importance=TRUE,
proximity=TRUE)
分類器構(gòu)建完畢之后唯咬，首先看一下這個(gè)分類器的準(zhǔn)確性
print(iris.rf)
# Call:
# randomForest(formula = Species ~ ., data = iris, importance = TRUE, proximity = TRUE)
# Type of random forest: classification
# Number of trees: 500
# No. of variables tried at each split: 2
# OOB estimate of error rate: 4.67%
# Confusion matrix:
# setosa versicolor virginica class.error
# setosa 50 0 0 0.00
# versicolor 0 47 3 0.06
# virginica 0 4 46 0.08
print 的結(jié)果中，OOB estimate of error rate 表明了分類器的錯(cuò)誤率為4.67%畏浆， Confusion matrix 表明了每個(gè)分類的詳細(xì)的分類情況胆胰；
對(duì)于setosa 這個(gè)group而言，基于隨機(jī)森林算法的分類器刻获，有50個(gè)樣本分類到了setosa 這個(gè)group, 而且這50個(gè)樣本和iris 中屬于setosa 這個(gè)group的樣本完全一致蜀涨，所以對(duì)于setosa 這個(gè)group而言，分類器的錯(cuò)誤率為0蝎毡；
對(duì)于versicolor 這個(gè)group而言厚柳，基于隨機(jī)森林算法的分類器，有47個(gè)樣本分類到了versicolor 這個(gè)group, 3個(gè)樣本分類到了virginica 這個(gè)group沐兵，有3個(gè)樣本分類錯(cuò)誤别垮，在iris 中屬于versicolor 這個(gè)group的樣本有50個(gè)，所以對(duì)于versicolor 這個(gè)group而言痒筒，分類器的錯(cuò)誤率為3/50 = 0.06 宰闰；
對(duì)于virginica 這個(gè)group而言茬贵，基于隨機(jī)森林算法的分類器簿透，有3個(gè)樣本分類到了versicolor 這個(gè)group, 47個(gè)樣本分類到了virginica 這個(gè)group，有3個(gè)樣本分類錯(cuò)誤解藻，在iris 中屬于virginica 這個(gè)group的樣本有50個(gè)老充，所以對(duì)于virginica這個(gè)group而言，分類器的錯(cuò)誤率為3/50 = 0.06 螟左；
然后看一下樣本之間的距離
iris.mds <- cmdscale(1 - iris.rf$proximity, eig=TRUE)
通過調(diào)用cmdscale 函數(shù)進(jìn)行樣本之間的距離啡浊，proximity 是樣本之間的相似度矩陣，所以用1減去之后得到樣本的類似距離矩陣的一個(gè)矩陣
iris.mds 的結(jié)果如下
str(iris.mds)
# List of 5
# $ points: num [1:150, 1:2] -0.565 -0.565 -0.566 -0.566 -0.565 ...
# ..- attr(*, "dimnames")=List of 2
# .. ..$ : chr [1:150] "1" "2" "3" "4" ...
# .. ..$ : NULL
# $ eig : num [1:150] 24 21.15 2.3 1.34 1.18 ...
# $ x : NULL
# $ ac : num 0
# $ GOF : num [1:2] 0.737 0.799
head(iris.mds$points)
# [,1] [,2]
# 1 -0.5649559 -0.04385526
# 2 -0.5651790 -0.04419585
# 3 -0.5656595 -0.04355633
# 4 -0.5657237 -0.04388197
# 5 -0.5651469 -0.04390029
# 6 -0.5655295 -0.04356653
在iris.mds 中points可以看做每個(gè)樣本映射到2維空間中的坐標(biāo)胶背，每一維空間是一個(gè)分類特征巷嚣，但是不是最原始的4個(gè)特征，而是由4個(gè)特征衍生得到的新的分類特征钳吟，根據(jù)這個(gè)坐標(biāo)廷粒，可以畫一張散點(diǎn)圖，得到每個(gè)樣本基于兩個(gè)分類變量的分組情況
plot(iris.mds$points, col = rep(c("red", "blue", "green"), each = 50))
圖中不同分類的樣本用不同的顏色標(biāo)注，可以看到基于兩個(gè)新的分類特征坝茎，樣本的分組效果還是很好的涤姊，不同組的樣本明顯區(qū)分開來。
最后嗤放，在看一下4個(gè)特征思喊，每個(gè)特征的重要性
iris.rf$importance
# setosa versicolor virginica MeanDecreaseAccuracy MeanDecreaseGini
# Sepal.Length 0.02837234 0.0169917984 0.043822619 0.029560910 9.371080
# Sepal.Width 0.01019719 0.0007153689 0.008415053 0.006240133 2.446698
# Petal.Length 0.32043066 0.2938837690 0.284731709 0.297205369 42.134304
# Petal.Width 0.34286020 0.3205494688 0.281664051 0.312822368 45.284763
之前調(diào)用randomForest 函數(shù)時(shí)，通過指定importance = TRUE 來計(jì)算每個(gè)特征的importance , 在 iris.rf$importance 矩陣中次酌，有兩個(gè)值是需要重點(diǎn)關(guān)注的MeanDecreaseAccuracy 和 MeanDecreaseGini
我們還可以利用
varImpPlot(iris.rf, main = "Top 30 - variable importance")
圖中和坐標(biāo)為importance 結(jié)果中的MeanDecreaseAccuracy 和 MeanDecreaseGini 指標(biāo)的值恨课，縱坐標(biāo)為對(duì)應(yīng)的每個(gè)分類特征，該函數(shù)默認(rèn)畫top30個(gè)特征岳服，由于這個(gè)數(shù)據(jù)集只有4個(gè)分類特征庄呈，所以4個(gè)都出現(xiàn)了。
??????使用randomForest對(duì)轉(zhuǎn)錄組的差異基因進(jìn)行進(jìn)一步的篩選
比如轉(zhuǎn)錄組測序中最后篩選出12個(gè)差異基因派阱，想要對(duì)他們進(jìn)行重要度的排序诬留，需要把expr表格整理成如下格式：前面12列是12個(gè)基因的表達(dá)量，必須是數(shù)值型贫母。最后一列是分組文兑，必須是factor型。然后就可以建模
rf <- randomForest(group ~ ., data=expr, importance=TRUE, proximity=TRUE) varImpPlot(rf, main = "Variable importance")
4. 劃分訓(xùn)練腺劣、測試集
前面的演示中绿贞，150個(gè)樣本全部作為了訓(xùn)練集，在這里演示一下劃分了訓(xùn)練橘原、測試集的模型
4.1 載入包和演示數(shù)據(jù)集
library(randomForest)
library(caret)
library(pROC)
data(iris)
4.2 劃分訓(xùn)練籍铁、測試集（一般來說比例是8:2，或者7:3）
??訓(xùn)練集和測試集的變量的數(shù)目必須是一致的趾断，否則會(huì)報(bào)錯(cuò)（這里都是5）拒名。
# 80%是訓(xùn)練集，20%是測試集芋酌。
trainlist <- createDataPartition(iris$Species,p=0.8,list = FALSE)
trainset <- iris[trainlist,]
testset <- iris[-trainlist,]
dim(trainset)
# [1] 120 5
dim(testset)
# [1] 30 5
4.3 使用訓(xùn)練集建模
set.seed(6666)
rf.train <- randomForest(as.factor(Species) ~ ., data=trainset, importance=TRUE,
proximity=TRUE,na.action = na.pass)
rf.train
# Call:
# randomForest(formula = as.factor(Species) ~ ., data = trainset, importance = TRUE, proximity = TRUE, na.action = na.pass)
# Type of random forest: classification
# Number of trees: 500
# No. of variables tried at each split: 2
# OOB estimate of error rate: 6.67%
# Confusion matrix:
# setosa versicolor virginica class.error
# setosa 40 0 0 0.000
# versicolor 0 37 3 0.075
# virginica 0 5 35 0.125
plot(rf.train,main = 'randomforest origin') #繪圖查看訓(xùn)練效果
4.4 預(yù)測測試集
set.seed(6666)
rf.test <- predict(rf.train,newdata = testset,type = 'class')
# 11 19 20 21 31 35 37 44
# setosa setosa setosa setosa setosa setosa setosa setosa
# 47 49 60 63 65 66 69 75
# setosa setosa versicolor versicolor versicolor versicolor versicolor versicolor
# 76 79 93 97 105 106 108 109
# versicolor versicolor versicolor versicolor virginica virginica virginica virginica
# 112 113 116 118 121 142
# virginica virginica virginica virginica virginica virginica
# Levels: setosa versicolor virginica
統(tǒng)計(jì)預(yù)測結(jié)果
rf.cf <- caret::confusionMatrix(as.factor(rf.test),as.factor(testset$Species))
4.5 ROC, AUC
首先要注意的是旱易，ROC需要的是概率，而不是預(yù)測的種類。所以要再預(yù)測一測咒唆，把type改成prob
rf.test2 <- predict(rf.train,newdata = testset,type = 'prob')
head(rf.test2)
# setosa versicolor virginica
# 11 1.00 0.00 0
# 19 0.98 0.02 0
# 20 1.00 0.00 0
# 21 1.00 0.00 0
# 31 1.00 0.00 0
# 35 1.00 0.00 0
繪制roc（因?yàn)槲覀冇?個(gè)class届垫，所以使用的是multiclass.roc，多種類的ROC）
roc.rf <- multiclass.roc(testset$Species,rf.test2)
roc.rf
# Call:
# multiclass.roc.default(response = testset$Species, predictor = rf.test2)
# Data: multivariate predictor rf.test2 with 3 levels of testset$Species: setosa, versicolor, virginica.
# Multi-class area under the curve: 1
可以看到全释，AUC=1装处。
