R語言與機器學習(分類算法)支持向量機
說到支持向量機,必須要提到july大神的《支持向量機通俗導論》,個人感覺再怎么寫也不可能寫得比他更好的了。這也正如青蓮居士見到崔顥的黃鶴樓后也只能嘆“此處有景道不得”。不過我還是打算寫寫SVM的基本想法與libSVM中R的接口。
回到我們最開始討論的KNN算法,它占用的內存十分的大,而且需要的運算量也非常大。那么我們有沒有可能找到幾個最有代表性的點(即保留較少的點)達到一個可比的效果呢?
要回答這個問題,我們首先必須思考如何確定點的代表性?我想關于代表性至少滿足這樣一個條件:無論非代表性點存在多少,存在與否都不會影響我們的決策結果。顯然如果仍舊使用
KNN算法的話,是不會存在訓練集的點不是代表點的情況。那么我們應該選擇一個怎樣的“距離”滿足僅依靠代表點就能得到全體點一致的結果?
我們先看下面一個例子:假設我們的訓練集分為正例與反例兩類,分別用紅色的圓圈與藍色的五角星表示,現(xiàn)在出現(xiàn)了兩個未知的案例,也就是圖中綠色的方塊,我們如何去分類這兩個例子呢?
在
KNN算法中我們考慮的是未知樣例與已知的訓練樣例的平均距離,未知樣例與正例和反例的“距離”誰更近,那么他就是對應的分類。
同樣是利用距離,我們可以換一個方式去考慮:假設圖中的紅線是對正例與反例的分類標準(記為w ? x+b=0),那么我們的未知樣例與紅線的“距離”就成了一個表示分類信度的標準,而w ? y+b(y為未知樣例的數(shù)據(jù))的符號則可以看成是分類的標識。
但是遺憾的是我們不知道這樣的一條分類標準(分類線)是什么,那么我們一個比較自然的想法就是從已知的分類數(shù)據(jù)(訓練集)里找到離分割線最近的點,確保他們離分割面盡可能的遠。這樣我們的分類器會更穩(wěn)健一些。
從上面的例子來看,虛線穿過的樣例便是離分割線最近的點,這樣的點可能是不唯一的,因為分割線并不確定,下圖中黑線穿過的訓練樣例也滿足這個要求:
所以“他們離分割面盡可能的遠”這個要求就十分重要了,他告訴我們一個穩(wěn)健的超平面是紅線而不是看上去也能分離
數(shù)據(jù)的黃線。
這樣就解決了我們一開始提出的如何減少儲存量的問題,我們只要存儲虛線劃過的點即可(因為在w ? x+b=-1左側,w ? x+b=1右側的點無論有多少都不會影響決策)。像圖中虛線劃過的,距離分割直線(比較專業(yè)的術語是超平面)最近的點,我們稱之為支持向量。這也就是為什么我們這種分類方法叫做
支持向量機的原因。
至此,我們
支持向量機的分類問題轉化為了如何尋找最大間隔的優(yōu)化問題。
支持向量機的實現(xiàn)涉及許多有趣的細節(jié):如何最大化間隔,存在“噪聲”的數(shù)據(jù)集怎么辦,對于線性不可分的數(shù)據(jù)集怎么辦等。
我這里不打算討論具體的算法,因為這些東西完全可以參閱july大神的《
支持向量機通俗導論》,我們這里只是介紹遇到問題時的想法,以便分析數(shù)據(jù)時合理調用R中的函數(shù)。
幾乎所有的
機器學習問題基本都可以寫成這樣的數(shù)學表達式:
給定條件:n個獨立同分布觀測樣本(x1 , y1 ), (x2 , y2 ),……,(xn , yn )
目標:求一個最優(yōu)函數(shù)f (x,w* )
最理想的要求:最小化期望風險R(w)
不同的是我們如何選擇f,R。對于
支持向量機來說,f(x,w*)=w? x+b,最小化風險就是最大化距離|w ?x|/||w||,即arg max{min(label ?(w ?x+b))/||w||} (也就是對最不confidence 的數(shù)據(jù)具有了最大的 confidence)
這里的推導涉及了對偶問題,拉格朗日乘子法與一堆的求導,我們略去不談,將結果敘述如下:
我們以鳶尾花數(shù)據(jù)來說說如何利用svm做分類,由于svm是一個2分類的辦法,所以我們將鳶尾花數(shù)據(jù)也分為兩類,“setosa”與“versicolor”(將后兩類均看做一類),那么數(shù)據(jù)按照
特征:花瓣長度與寬度做分類,有分類:
從上圖可以看出我們通過最優(yōu)化原始問題或者對偶問題就可以得到w,b,利用 sign(w.x+b)就可以判斷分類了。
我們這里取3, 10,56, 68,107, 120號數(shù)據(jù)作為測試集,其余的作為訓練集,我們可以看到:
訓練集 setosa virginica
setosa 48 0
virginica 0 96
測試集 setosa virginica
setosa 2 0
virginica 0 4
也就是完全完成了分類任務。
我們來看看鳶尾花后兩類的分類versicolor和virginica的分類,我們將數(shù)據(jù)的
散點圖描繪如下:(我們把第一類“setosa“看做”versicolor“)
不難發(fā)現(xiàn)這時無論怎么畫一條線都無法將數(shù)據(jù)分開了,那么這么辦呢?我們一個自然的辦法就是允許分類有一部分的錯誤,但是錯誤不能無限的大。我們使用一個松弛項來分類數(shù)據(jù)。最優(yōu)化問題轉變?yōu)椋?
當我們確定好松弛項C后,就可以得到分類:
我們還是先來看看分類的效果:(C=10)
訓練集 versicolor virginica
versicolor 93 2
virginica 3 46
測試集 versicolor virginica
versicolor 4 2
virginica 0 0
雖然分類中有一些錯誤,但是大部分還是分開了的,也就是說還是不錯的,至少完成了分類這個任務。
我們再來看看一個更麻煩的例子:
假設數(shù)據(jù)是這樣的:
這時再用直線作為劃分依據(jù)就十分的蹩腳了,我們這時需要引入核的方法來解決這個問題。
在上圖中,我們一眼就能看出用一個S型去做分類就可以把數(shù)據(jù)成功分類了(當然是在允許一點點錯誤的情況下),但是計算機能識別的只有分類器的分類結果是-1還是1,這時,我們需要將
數(shù)據(jù)做出某種形式的轉換,使得原來不可用直線剖分的變得可分,易分。也就是需要找到一個從一個
特征空間到另一個
特征空間的映射。我們常用的映射有:
線性核:u'*v
多項式核:(gamma*u'*v+ coef0)^degree
高斯核:exp(-gamma*|u-v|^2)
Sigmoid核:tanh(gamma*u'*v + coef0)
我們這里使用各種常見的核來看看分類效果:
從圖中我們可以看到正態(tài)核的效果是最好的,用數(shù)據(jù)說話,我們來看看分類錯誤率與折10交叉驗證的結果(報告平均分類正確率):
我們可以看到,無論從存儲數(shù)據(jù)量的多少(支持向量個數(shù))還是分類精確度來看,高斯核都是最優(yōu)的。所以一般情況,特別是在大樣本情況下,優(yōu)先使用高斯核,至少可以得到一個不太壞的結果(在完全線性可分下,線性函數(shù)的支持向量個數(shù)還是少一些的)。
有許多介紹
SVM的書都有類似的表述“由于理解
支持向量機需要掌握一些理論知識,而這對讀者來說有一定的難度,建議直接下載LIB
SVM使用。”
確實,如果不是為了訓練一下編程能力,我們沒有必要自己用前面提到的做法自己實現(xiàn)一個效率不太高的
SVM。R的函數(shù)包e1071提供了lib
SVM的接口,使用e1071的函數(shù)
SVM()可以得到lib
SVM相同的結果,write.svm()更是可以把R訓練得到的結果寫為標準的lib
SVM格式供其他環(huán)境下的lib
SVM使用。
在介紹R中函數(shù)的用法時,我們先簡要介紹一下
SVM的類型,以便我們更好地理解各個參數(shù)的設置。
對于線性不可分時,加入松弛項,折衷考慮最小錯分樣本和最大分類間隔。增加了算法的容錯性,允許訓練集不完全分類,以防出現(xiàn)
過擬合。加入的辦法有以下3類,寫成最優(yōu)化問題形式總結如上圖:
上圖中e為所有元素都為1的列向量,Qij=yiyjK(xi; xj), K(xi; xj) =phi(xi) ?phi (xj), phi(.)為核函數(shù),K(. ;.)表示對應元素核函數(shù)的內積。
現(xiàn)在我們來看看svm()函數(shù)的用法。
## S3 method for class 'formula'
svm(formula, data = NULL, ..., subset,na.action =na.omit, scale = TRUE)
## Default S3 method:
svm(x, y = NULL, scale = TRUE, type = NULL,kernel =
"radial", degree = 3, gamma = if(is.vector(x)) 1 else 1 / ncol(x),
coef0 = 0, cost = 1, nu = 0.5,
class.weights = NULL, cachesize = 40,tolerance = 0.001, epsilon = 0.1,shrinking = TRUE, cross = 0, probability =FALSE, fitted = TRUE, seed = 1L,..., subset, na.action = na.omit)
主要參數(shù)說明:
Formula:分類模型形式,在第二個表達式中使用的的x,y可以理解為y~x。
Data:數(shù)據(jù)集
Subset:可以指定數(shù)據(jù)集的一部分作為訓練集
Na.action:
缺失值處理,默認為刪除數(shù)據(jù)條目
Type:
SVM的形式,使用可參見上面的
SVMformulation,type的選項有:C-classification,nu-classification,one-classification (for novelty detection),eps-regression,nu-regression。后面兩者為利用
SVM做回歸時用到的,這里暫不介紹。默認為C分類器,使用nu分類器會使決策邊界更光滑一些,單一分類適用于所有的訓練數(shù)據(jù)提取自同一個類里,然后
SVM建立了一個分界線以分割該類在
特征空間中所占區(qū)域和其它類在
特征空間中所占區(qū)域。
Kernel:在非線性可分時,我們引入核函數(shù)來做非線性可分,R提供的核介紹如下:
線性核:u'*v
多項式核:(gamma*u'*v + coef0)^degree
高斯核:exp(-gamma*|u-v|^2)
Sigmoid核:tanh(gamma*u'*v + coef0)
默認為高斯核(RBF),lib
SVM的作者對于核的選擇有如下建議:Ingeneral we suggest you to try the RBF kernel first. A recent result by Keerthiand Lin shows that if RBF is used with model selection, then there is no need to consider the linear kernel. The kernel matrix using sigmoid may not be positive definite and in general it's accuracy is not better than RBF. (see thepaper by Lin and Lin. Polynomial kernels are ok but if a high degree is used,numerical difficulties tend to happen (thinking about dth power of (<1) goes to 0 and (>1) goes to infinity).
順帶說一句,在kernlab包中,可以自定義核函數(shù)。
Degree:多項式核的次數(shù),默認為3
Gamma:除去線性核外,其他的核的參數(shù),默認為1/數(shù)據(jù)維數(shù)
Coef0,:多項式核與sigmoid核的參數(shù),默認為0
Cost:C分類的懲罰項C的取值
Nu:nu分類,單一分類中nu的取值
Cross:做K折交叉驗證,計算分類正確性。
由于svm的編程確實過于復雜,還涉及到不少最優(yōu)化的內容,所以在第二部分我的分類都是使用svm函數(shù)完成的(偷一下懶),現(xiàn)將部分R代碼展示如下:
dataSim的函數(shù):
[plain] view plaincopyprint
simData=function(radius,width,distance,sample_size)
{
aa1=runif(sample_size/2)
aa2=runif(sample_size/2)
rad=(radius-width/2)+width*aa1
theta=pi*aa2
x=rad*cos(theta)
y=rad*sin(theta)
label=1*rep(1,length(x))
x1=rad*cos(-theta)+rad
y1=rad*sin(-theta)-distance
label1=-1*rep(1,length(x1))
n_row=length(x)+length(x1)
data=matrix(rep(0,3*n_row),nrow=n_row,ncol=3)
data[,1]=c(x,x1)
data[,2]=c(y,y1)
data[,3]=c(label,label1)
data
}
dataSim=simData(radius=10,width=6,distance=-6,sample_size=3000)
colnames(dataSim)<-c("x","y","label")
dataSim<-as.data.frame(dataSim)
Sigmoid核的分類預測:
[plain] view plaincopyprint
m1 <- svm(label ~x+y, data =dataSim,cross=10,type="C-classification",kernel="sigmoid")
m1
summary(m1)
pred1<-fitted(m1)
table(pred1,dataSim[,3])
核函數(shù)那一小節(jié)作圖的各種東西:
[plain] view plaincopyprint
linear.svm.fit <- svm(label ~ x + y, data = dataSim, kernel ='linear')
with(dataSim, mean(label == ifelse(predict(linear.svm.fit) > 0,1, -1)))
polynomial.svm.fit <- svm(label ~ x + y, data = dataSim, kernel ='polynomial')
with(dataSim, mean(label == ifelse(predict(polynomial.svm.fit) >0, 1, -1)))
radial.svm.fit <- svm(label ~ x + y, data = dataSim, kernel ='radial')
with(dataSim, mean(label == ifelse(predict(radial.svm.fit) > 0,1, -1)))
sigmoid.svm.fit <- svm(label ~ x + y, data = dataSim, kernel ='sigmoid')
with(dataSim, mean(label == ifelse(predict(sigmoid.svm.fit) > 0,1, -1)))
df <- cbind(dataSim,
data.frame(Linear
SVM = ifelse(predict(linear.svm.fit) > 0, 1, -1),
Polynomial
SVM = ifelse(predict(polynomial.svm.fit) > 0, 1, -1),
Radial
SVM = ifelse(predict(radial.svm.fit) > 0, 1, -1),
Sigmoid
SVM = ifelse(predict(sigmoid.svm.fit) > 0, 1, -1)))
library("reshape")
predictions <- melt(df, id.vars = c('x', 'y'))
library('ggplot2')
ggplot(predictions, aes(x = x, y = y, color = factor(value))) +
geom_point() +
facet_grid(variable ~ .)
最后,我們回到最開始的那個手寫數(shù)字的案例,我們試著利用
支持向量機重做這個案例。
運行代碼:
[plain] view plaincopyprint
setwd("D:/R/data/digits/trainingDigits")
names<-list.files("D:/R/data/digits/trainingDigits")
data<-paste("train",1:1934,sep="")
for(i in 1:length(names))
assign(data[i],as.vector(as.matrix(read.fwf(names[i],widths=rep(1,32)))))
label<-rep(0:9,c(189,198,195,199,186,187,195,201,180,204))
data1<-get(data[1])
for(i in 2:length(names))
data1<-rbind(data1,get(data[i]))
m <- svm(data1,label,cross=10,type="C-classification")
m
summary(m)
pred<-fitted(m)
table(pred,label)
setwd("D:/R/data/digits/testDigits")
names<-list.files("D:/R/data/digits/testDigits")
data<-paste("train",1:1934,sep="")
for(i in 1:length(names))
assign(data[i],as.vector(as.matrix(read.fwf(names[i],widths=rep(1,32)))))
data2<-get(data[1])
for(i in 2:length(names))
data2<-rbind(data2,get(data[i]))
pred<-predict(m,data2)
labeltest<-rep(0:9,c(87,97,92,85,114,108,87,96,91,89))
table(pred,labeltest)
模型摘要:
Call:
svm.default(x = data1, y = label, type ="C-classification", cross =10)
Parameters:
SVM-Type: C-classification
cost: 1
gamma: 0.0009765625
Number of Support Vectors: 1139 ( 78 130 101 124 109 122 87 93 135 160 )
Number of Classes: 10
Levels: 0 1 2 3 4 5 6 7 8 9
10-fold cross-validation on training data:
Total Accuracy: 96.7425
Single Accuracies:
97.40933 98.96373 91.7525899.48187 94.84536 94.30052 97.40933 96.90722 98.96373 97.42268
當然,我們還可以通過m$SV查看支持向量的情況,m$index查看支持向量的標簽,m$rho查看分類時的截距b。
訓練集分類結果:
我們拿測試數(shù)據(jù)來看:
分類正確率為:0.9735729,誤差率為2.6%左右,確實達到了開篇提出的可比的目的,而需要儲存的支持向量個數(shù)僅為1139個,比原來的訓練
數(shù)據(jù)1934個要少了近50%,也達到了我們要求的節(jié)約存儲的目的。當然值得一提的是線性分類的效果在實際中也沒有那么糟糕,可以犧牲線性核函數(shù)的正確率來換取分類速度與存儲空間。另外,支持向量的個數(shù)與訓練集的出錯率也沒有特別必然的聯(lián)系,而是與容錯率cost有一定的聯(lián)系。
CDA數(shù)據(jù)分析師考試相關入口一覽(建議收藏):
? 想報名CDA認證考試,點擊>>>
“CDA報名”
了解CDA考試詳情;
? 想學習CDA考試教材,點擊>>> “CDA教材” 了解CDA考試詳情;
? 想加入CDA考試題庫,點擊>>> “CDA題庫” 了解CDA考試詳情;
? 想了解CDA考試含金量,點擊>>> “CDA含金量” 了解CDA考試詳情;