R語言：SMOTE - Supersampling Rare Events in R:用R對非平衡數據的處理方法

SMOTE - Supersampling Rare Events in R:用R對稀有事件進行超級采樣

在這個例子中將用到以下三個包
{DMwR} - Functions and data for the book “Data Mining with R” and SMOTE algorithm:SMOTE算法
{caret} - modeling wrapper, functions, commands：模型封裝、函數、命令

{pROC} - Area Under the Curve (AUC) functions：曲線下面積(ACU)函數

SMOTE算法是為了解決不平衡的分類問題。也就是說，它可以產生一個新的“SMOTEd”數據，解決類不平衡問題集。或者，它也可以運行在這個新的數據集的分類算法，并返回所得到的模型。

我們利用 Thyroid Disease 數據來進行研究。
讓我們清洗一些數據
# 加載數據，刪除冒號和句號，并追加列名
hyper <-read.csv('http://archive.ics.uci.edu/ml/machine-learning-databases/thyroid-disease/hypothyroid.data', header=F)
names <- read.csv('http://archive.ics.uci.edu/ml/machine-learning-databases/thyroid-disease/hypothyroid.names', header=F, sep='\t')`1`
names <- gsub(pattern =":|[.]", replacement="", x = names)
colnames(hyper)<-names
# 我們將第一列的列名從 hypothyroid, negative改成target，并將negative變成0，其他值變成1.
colnames(hyper)[1]<-"target"
colnames(hyper)
##  [1] "target"                    "age"                      
##  [3] "sex"                       "on_thyroxine"             
##  [5] "query_on_thyroxine"        "on_antithyroid_medication"
##  [7] "thyroid_surgery"           "query_hypothyroid"        
##  [9] "query_hyperthyroid"        "pregnant"                 
## [11] "sick"                      "tumor"                    
## [13] "lithium"                   "goitre"                   
## [15] "TSH_measured"              "TSH"                      
## [17] "T3_measured"               "T3"                       
## [19] "TT4_measured"              "TT4"                      
## [21] "T4U_measured"              "T4U"                      
## [23] "FTI_measured"              "FTI"                      
## [25] "TBG_measured"              "TBG"
hyper$target<-ifelse(hyper$target=="negative",0,1)
# 檢查下陽性和陰性的結果
table(hyper$target)
## 
##    0    1 
## 3012  151
prop.table(table(hyper$target))
## 
##       0       1 
## 0.95226 0.04774
# 可見,1僅有5%。這顯然是一個扭曲的數據集，也是罕見事件。
head(hyper,2)
##   target age sex on_thyroxine query_on_thyroxine on_antithyroid_medication
## 1      1  72   M            f                  f                         f
## 2      1  15   F            t                  f                         f
##   thyroid_surgery query_hypothyroid query_hyperthyroid pregnant sick tumor
## 1               f                 f                  f        f    f     f
## 2               f                 f                  f        f    f     f
##   lithium goitre TSH_measured TSH T3_measured   T3 TT4_measured TT4
## 1       f      f            y  30           y 0.60            y  15
## 2       f      f            y 145           y 1.70            y  19
##   T4U_measured  T4U FTI_measured FTI TBG_measured TBG
## 1            y 1.48            y  10            n   ?
## 2            y 1.13            y  17            n   ?
# 這數據都是因子型變量（字符型的值），這些都需要轉換成二值化的數字，以方便建模：
ind<-sapply(hyper,is.factor)
hyper[ind]<-lapply(hyper[ind],as.character)

hyper[hyper=="?"]=NA
hyper[hyper=="f"]=0
hyper[hyper=="t"]=1
hyper[hyper=="n"]=0
hyper[hyper=="y"]=1
hyper[hyper=="M"]=0
hyper[hyper=="F"]=1

hyper[ind]<-lapply(hyper[ind],as.numeric)

replaceNAWithMean<-function(x) {replace(x,is.na(x),mean(x[!is.na(x)]))}

hyper<-replaceNAWithMean(hyper)

模型研究
我們利用caret包中的createDataPartition（數據分割功能）函數將數據隨機分成相同的兩份。

library(caret)
## Loading required package: lattice
## Loading required package: ggplot2
set.seed(1234)
splitIndex<-createDataPartition(hyper$target,time=1,p=0.5,list=FALSE)
trainSplit<-hyper[splitIndex,]
testSplit<-hyper[-splitIndex,]

prop.table(table(trainSplit$target))
## 
##       0       1 
## 0.95006 0.04994
prop.table(table(testSplit$target))
## 
##       0       1 
## 0.95446 0.04554
兩者的分類結果是平衡的，因此仍然有5%左右的代表，我們仍然處于良好的水平。

我們利用caret包中的treebag模型算法，對訓練集數據建立模型，并對測試集數據進行預測。

ctrl<-trainControl(method="cv",number=5)
tbmodel<-train(target~.,data=trainSplit,method="treebag",
               trControl=ctrl)
## Loading required package: ipred
## Loading required package: plyr
predictors<-names(trainSplit)[names(trainSplit)!='target']
pred<-predict(tbmodel$finalModel,testSplit[,predictors])
為了評估模型，我們用pROC包的roc函數算auc得分和畫圖
library(pROC)
## Type 'citation("pROC")' for a citation.
## 
## Attaching package: 'pROC'
## 
## 下列對象被屏蔽了from 'package:stats':
## 
##     cov, smooth, var
auc<-roc(testSplit$target,pred)
print(auc)
## 
## Call:
## roc.default(response = testSplit$target, predictor = pred)
## 
## Data: pred in 1509 controls (testSplit$target 0) < 72 cases (testSplit$target 1).
## Area under the curve: 0.985
plot(auc,ylim=c(0,1),print.thres=TRUE,main=paste('AUC',round(auc$auc`1`,2)))
## 
## Call:
## roc.default(response = testSplit$target, predictor = pred)
## 
## Data: pred in 1509 controls (testSplit$target 0) < 72 cases (testSplit$target 1).
## Area under the curve: 0.985
abline(h=1,col="blue",lwd=2)
abline(h=0,col="red",lwd=2)




auc得分是0.98，已經是非常不錯的結果了（因為它的范圍是在0.5到1之間）。

很難想象SMOTE對此能再有提高了，但接下來我們利用SMOTE對數據處理后再建模，看看auc結果

在R中，SMOTE算法是DMwR軟件包的一部分，主要參數有如下三個：perc.over:過采樣時，生成少數類的樣本個數;k:過采樣中使用K近鄰算法生成少數類樣本時的K值，默認是5；perc.under:欠采樣時，對應每個生成的少數類樣本，選擇原始數據多數類樣本的個數。例如，perc.over=500表示對原始數據集中的每個少數樣本，都將生成5個新的少數樣本；perc.under=80表示從原始數據集中選擇的多數類的樣本是新生的數據集中少數樣本的80%。

library(DMwR)
## Loading required package: grid
## 
## Attaching package: 'DMwR'
## 
## 下列對象被屏蔽了from 'package:plyr':
## 
##     join
trainSplit$target<-as.factor(trainSplit$target)
trainSplit<-SMOTE(target~.,trainSplit,perc.over=100,perc.under=200)
trainSplit$target<-as.numeric(trainSplit$target)
# 我們再次用prop.table()函數檢查結果的平衡性，確定我們已經讓陰性、陽性數據達到相同。
prop.table(table(trainSplit$target))
## 
##   1   2 
## 0.5 0.5
# 再次建立treebag模型
tbmodel<-train(target~.,data=trainSplit,method="treebag",
               trControl=ctrl)
predictors<-names(trainSplit)[names(trainSplit)!='target']
pred<-predict(tbmodel$finalModel,testSplit[,predictors])
auc<-roc(testSplit$target,pred)
print(auc)
## 
## Call:
## roc.default(response = testSplit$target, predictor = pred)
## 
## Data: pred in 1509 controls (testSplit$target 0) < 72 cases (testSplit$target 1).
## Area under the curve: 0.99
哇，達到0.99，比之前的0.985有提高
plot(auc,ylim=c(0,1),print.thres=TRUE,main=paste('AUC',round(auc$auc`1`,2)))
## 
## Call:
## roc.default(response = testSplit$target, predictor = pred)
## 
## Data: pred in 1509 controls (testSplit$target 0) < 72 cases (testSplit$target 1).
## Area under the curve: 0.99
abline(h=1,col="blue",lwd=2)

abline(h=0,col="red",lwd=2)