乾貨--C5.0與CART演算法實戰

來自專欄 R語言中文社區7 人贊了文章

作者：劉順祥

公眾號：數據分析1480 (微信ID：lsxxx2011)

配套教程：手把手教你做文本挖掘 https://edu.hellobi.com/course/181

在上一期的《淺談C5.0與CART演算法的比較--理論理解》我們詳細講解了有關C5.0決策樹和CART決策樹的理論知識，包括構造樹過程中如何選擇節點變數、節點變數的分割點、如何完成剪枝，避免模型的過擬合，從而增強樹模型的泛化能力。接下來我們將從實際的案例中來比較兩個演算法的實現，希望讀者在閱讀本文時能夠再次查看上面提到的那篇理論理解，這將有助於解讀後文的落地。

本文實戰部分的數據來自於UCI機器學習網站（http://archive.ics.uci.edu/ml/），後文會給出腳本及數據下載的鏈接。該數據的因變數是反映某銀行顧客是否會繳納一項保證金，總共有45211條記錄和16個自變數。各個自變數的含義如下：

age：年齡；

job：工作類型（如行政、管理者、失業等）；

marital：婚姻狀態（已婚、未婚、離異）；

education：教育程度（初等教育、中等教育、高等教育）；

default：是否擁有信用卡；

balance：平均年餘額；

housing：是否有房貸；

loan：是否有個人貸款；

contact：聯繫方式（固定電話、手機）；

day：最後一次聯繫的日；

month；最後一次聯繫的月；

duration；最後一次聯繫的時長（秒為單位）；

campaign：在本次市場活動中聯繫的次數；

pdays：最後一次聯繫的時間距離上一次市場活動的間隔天數（-1表示該用戶在上一次活動中沒有聯繫）；

previous：上一次活動中，聯繫的次數；

poutcome：上一次市場活動的結果（成功、失敗、其他）；

接下來我們就利用上面所說的數據集進行建模：

# 載入所需的第三方包if(!suppressWarnings(require(Hmisc))){ install.packages(Hmisc) require(Hmisc)}if(!suppressWarnings(require(C50))){ install.packages(C50) require(C50)}if(!suppressWarnings(require(rpart))){ install.packages(rpart) require(rpart)}# 讀取數據集mydata <- read.csv(file = file.choose(), sep = ;)head(mydata)

# 篩選出因子型變數，並對這些變數作統計factors <- names(mydata)[sapply(mydata,class) == factor]sapply(mydata[,factors], table)

上圖為截取的一部分。

# 數據清洗# 將job中未知職業的記錄刪除（僅佔0.64%）,並刪除缺失嚴重的變數poutcome（佔82%）。clear.mydata <- subset(mydata, job != unknown, select = -poutcome)dim(mydata)dim(clear.mydata)

# 受教育程度中有1857個未知，我們不妨用眾數（secondary）替補education.impute <- ifelse(clear.mydata$education != unknown, as.character(clear.mydata$education), secondary)# 至少有28%的觀測在contact變數上是缺失的，不妨我們按工作種類分組填補。# 首先將unknown設置為R中的缺失標誌NAclear.mydata$contact <- as.character(clear.mydata$contact)table(clear.mydata$contact)

clear.mydata$contact[clear.mydata$contact == unknown] <- NAtable(clear.mydata$contact, useNA = ifany)

由於R中沒有自帶的眾數函數，這裡我們自定義一個眾數函數。

# 自定義眾數函數stat.mode <- function(x, rm.na = TRUE){ if (rm.na == TRUE){ y = x[!is.na(x)] } res = names(table(y))[which.max(table(y))] return(res)}# 自定義函數，實現分組替補my.impute <- function(data, category.col = NULL, miss.col = NULL, method = stat.mode){ impute.data = NULL for(i in as.character(unique(data[,category.col]))){ sub.data = subset(data, data[,category.col] == i) sub.data[,miss.col] = impute(sub.data[,miss.col], method) impute.data = c(impute.data, sub.data[,miss.col]) } data[,miss.col] = impute.data return(data)}其中，category.col 指定所要分組的變數，miss.col指定需要填補的缺失值變數，默認的方法為眾數填補。clear.mydata <- my.impute(clear.mydata, category.col = job, miss.col = contact)table(clear.mydata$contact, useNA = ifany)

很顯然，那些缺失的觀測全被替補成了"cellular"溝通方式，說明在各組中聯繫方式的眾數為"cellular"。

# 再將字元串變數轉換為因子型變數clear.mydata$contact <- factor(clear.mydata$contact)# 數據合併final.data <- cbind(clear.mydata,education.impute)final.data <- final.data[,-4]# 簡單的了解一下數據str(final.data)

summary(final.data)

上面的過程全都是數據預處理的過程，接下來我們要對處理好的數據進行建模和預測：

# 抽樣，並將總體分為訓練集和測試集set.seed(1)index <- sample(1:nrow(final.data), size = 0.75*nrow(final.data))train <- final.data[index,]test <- final.data[-index,]# 大致查看抽樣與總體之間是否吻合prop.table(table(final.data$y))prop.table(table(train$y))prop.table(table(test$y))

# 構建C5.0決策樹，並對重要變數進行篩選fit <- C5.0(x = train[,-15], y = train[,15], control = C5.0Control(winnow = TRUE))summary(fit)

從結果中看，模型選擇的重要變數為duration,housing,month,campaign,previous,day,pdays,marital,loan,age，接下來我們就利用這些變數，多模型進行修正。

# 建模並預測vars <- c(y,duration,housing,month, campaign,previous,day,pdays, marital,loan,age)train2 <- train[,vars]test2 <- test[,vars]# 建模fit1 <- C5.0(x = train2[,-1], y = train2[,1])# 預測pred1 <- predict(fit1, newdata = test2[,-1])# 混淆矩陣freq1 <- table(pred1, test2[,1])freq1

# 準確率accuracy1 <- sum(diag(freq1))/sum(freq1)accuracy1# 正例的覆蓋率recall1 <- freq1[2,2]/sum(freq1[,2])recall1

雖然模型的準確率達到90%以上，但預測正確的yes在實際的yes中只佔了51.8%，即正例的覆蓋率並不高，模型的準確性值得懷疑。

C5.0演算法可通過錯誤率和損失矩陣進行剪枝，之前的文章提過，默認的alpha（置信水平）為0.25，當alpha設置低於0.25時，將會進行剪枝。為了確定最佳的alpha值，我們自定一個函數，通過遍歷的方式確定alpha。

# 剪枝--基於錯誤率的剪枝法err.rate <- function(train, test, y.index = NULL, y.name = NULL){ alpha <- NULL res <- NULL if(is.null(y.index)){ y.index = which(names(train) == y.name) } for (i in seq(0.25,0.1,-0.01)){ fit <- C5.0(x = train[,-y.index], y = train[,y.index], control = C5.0Control(CF = i)) pred <- predict(fit, newdata = test[,-y.index]) freq <- table(pred, test[,y.index]) accuracy <- sum(diag(freq))/sum(freq) alpha <- c(alpha,i) res <- c(res,accuracy) } return(data.frame(alpha,res))}err.rate(train2, test2, y.name = y)

根據上面的結果，我們確定alpha值為0.23，此時模型的準確率提高了一點點，於是我們基於這個值，再結合損失矩陣再做一次模型的構建和預測。

# 構建損失矩陣（注意必須為矩陣設置行名稱和列名稱）costs <- matrix(c(0,4,1,0), ncol = 2, byrow = TRUE, dimnames = list(unique(train2$y),unique(train2$y)))# 同過control參數設置alpha值fit3 <- C5.0(x = train2[,-1], y = train2[,1], control = C5.0Control(CF = 0.23), costs = costs)# 預測pred3 <- predict(fit3, newdata = test2[,-1])freq3 <- table(pred3, test2[,1])freq3

accuracy3 <- sum(diag(freq3))/sum(freq3)accuracy3recall3 <- freq3[2,2]/sum(freq3[,2])recall3

模型經過改善後，大大提高了正例的覆蓋率，從原來的51.8%提升到目前的80.8%，雖然模型的整體準確率降低了3.4個百分點，但這樣的損失在一定程度上是有助於業務市場的活動，因為能夠預測到更多的yes對象，就可以對這些群體進行營銷，改善業務。

接下來我們再試試CART演算法在該數據集上應用：

# 構建CART演算法fit4 <- rpart(y ~ ., data = train2)# 預測pred4 <- predict(fit4, newdata = test2[,-1], type = class)# 構建混淆矩陣freq4 <- table(pred4, test2[,1])freq4

# 模型準確率accuracy4 <- sum(diag(freq4))/sum(freq4)accuracy4

發現模型的準確率也挺高的，也在90%，但誇張的是正例的覆蓋率只有32.7%，比C5.0模型什麼都不做還差很多，這就需要我們對CART演算法進行剪枝操作。

首先來看一下模型的cp表，可以通過cp值進行「最小代價複雜度」剪枝：

從結果中看，cp值為0.01時，誤差率最低，切xerror+xstd也是達到最小，而模型構造的時候默認就是cp=0.01，故暫不需要通過cp值進行剪枝。那看看是否可以通過損失矩陣進一步優化模型：

# 剪枝--基於損失矩陣costs <- matrix(c(0,1.25,1,0), ncol = 2, byrow = TRUE)fit5 <- rpart(y ~ ., data = train2, parms = list(loss = costs))fit5# 預測pred5 <- predict(fit5, newdata = test2[,-1], type = class)freq5 <- table(pred5, test2[,1])freq5

accuracy5 <- sum(diag(freq5))/sum(freq5)accuracy5

經過不停的嘗試，模型效果一直不夠滿意，如果損失矩陣設置的大一點就會大失所望。。。如下方所示：

決策樹就剩下根節點了，並沒有進行樹的構造。我懷疑該數據集並不適合使用CART演算法進行樹的構造。如果有其他觀點的歡迎朋友們聯繫我。

OK，我們這期的實戰部分就到這裡。歡迎大家多多溝通和交流，通過互相學習，達到取長補短的效果。快要過年了，提前祝福各位網友和朋友2017年新年快樂，心想事成，萬事如意！

數據集合腳本鏈接：

鏈接：http://pan.baidu.com/s/1ge98tG3 密碼：was9