几个分类器和对机器学习的一些思考

最近一周，下班后都在kaggle上玩。kaggle现在做的越来越好，对于我们这种搞机器学习的新手来说，这个网站提供了很多实际数据。关注这个网站也有快一年时间了，这里不仅传递着数据和算法，更多的是传递着全球数据爱好者们对数据的信仰。

最近在做的是手写数字识别的那道题。手写数字识别，是机器学习用于解决实际问题的一个相当成功的案例。几乎所有机器学习的书籍中都会提到这个例子，国内甚至还有人专门利用各种算法解决该问题而写成了一本不薄的书。这个问题的解决，涉及到很多内容，除去机器学习中的一些分类算法之外，更多的是工程上的一些问题。比如如何收集大量的手写数字样本，如何提取图像中的信息，在诸多的学习算法中如何找出一个最合适的（未必是预测准确度最高的，也许还要考虑时间，成本等其他的因素）。工程师的思维和一个理科学生的思维之间的差别，在于，工程师往往需要综合考虑很多因素，得出一个最适宜的解决问题的方案，而一个理科学生，可能会对自己所掌握的算法数量和难度引以为傲，但是却会忽视算法之外的一些事。

关于手写数字识别这个例子，数据情况是，训练样本量73.2M，测试集48.7M。每一个观测，是一个28*28像素的手写数字的图片，该图片用长度为28*28的向量来表示，在训练集中，人工标识出了每一个观测对应的数字标签。这是一个典型的有指导学习问题。目的是从训练集中学习出数字识别的规则，然后在测试集上进行检验。目前kaggle上的最好成绩（测试集上的识别准确率）已经达到了99.47%。

我也尝试了几种不同算法来做这个问题，这个过程中也遇到了很多不可避免的工程方面的问题。

我依然是使用最为顺手的R作为工具，但是由于32位版本最多只能分配2g内存，而很多算法需要储存大矩阵，因此总会出现“无法分配xx内存的矢量”这样的错误提示，令人吐血不止。关于R语言对大数据的处理，目前已经有些解决方案，R本身也有bigmemory包之类的。但是内存依然无法满足需要。这也是R语言最为人所诟病的一点。

我所采用的解决方法比较简单，每次随机从训练样本中抽出三分之一的数据作为训练集，训练完成得到模型之后便把所用的训练数据删掉（腾出内存来），然后用模型在检测集上进行分类，将分类结果输出写到硬盘上，然后删除内存上的相应结果。重复多次这样的过程，得到多组结果，然后再按照少数服从多数的原则在多组结果中进行投票，得到最终的预测结果。这样的方法，显得很不靠谱。也很令人心痛，因为明明有那么多数据可以用，但是由于计算机和自身水平的限制无法利用，只能提取一小部分来利用。这是令人气郁的。然而，即便只用了三分之一的数据，可是当你还是发现，训练依然需要消耗大量的时间时，你也只能遗憾的接受这样残酷的现实了。

算法方面，我选取了naive bayes，k-nearest neighborhood（k=10），generalized linear model（with lasso regularization），randomForest以及support vector machine。这些方法都是最为基本（常用）的机器学习算法。除了naive bayes算法预测准确度较低，60%左右之外，剩下几种算法的预测准确度都在92%到96%。其中SVM算法精度最高，这个结果是预料之中的。令我感到比较意外的是加了L1范数惩罚项的广义线性模型的效果也比较出色。这个模型使用了cross-validation，因此整个过程非常耗时，这个模型仅使用了三分之一的样本，并不似其他算法，我都做了model average，但是还是得到了比较满意的预测效果。让我又一次对Efron，hastie，freedman他们充满了敬仰。

照例，代码将共享。

#####prepare the data#####
setwd("F:/kaggle/numRecognition")
training=read.csv("train.csv",header=TRUE)
test=read.csv("test.csv",header=TRUE)
#Check the data set
#head(training)
#head(test)
label=training[,1]
label.f <- as.factor(label)
train=as.matrix(training[,-1])
test=as.matrix(test)
###in order to run the code on my computer,i have to split the train set and make
###model average,what a fuck!!




#####KNN#####
library(FNN)
idx.1=sample(c(1:42000),size=42000/3,replace=TRUE)
train.1=train[idx.1,]
label.1=label[idx.1]
res.KNN.1 =knn(train.1, test, label.1, k = 10, algorithm="cover_tree")
write(res.KNN.1, file="RESKNN1.csv", ncolumns=1)
rm(train.1)
rm(label.1)
rm(res.KNN.1)
rm(idx.1)

idx.2=sample(c(1:42000),size=42000/3,replace=TRUE)
train.2=train[idx.2,]
label.2=label[idx.2]
res.KNN.2 =knn(train.2, test, label.2, k = 10, algorithm="cover_tree")
write(res.KNN.2, file="RESKNN2.csv", ncolumns=1)
rm(train.2)
rm(label.2)
rm(res.KNN.2)
rm(idx.2)

idx.3=sample(c(1:42000),size=42000/3,replace=TRUE)
train.3=train[idx.3,]
label.3=label[idx.3]
res.KNN.3 =knn(train.3, test, label.3, k = 10, algorithm="cover_tree")
write(res.KNN.3, file="RESKNN3.csv", ncolumns=1)
rm(train.3)
rm(label.3)
rm(res.KNN.3)
rm(idx.3)

idx.4=sample(c(1:42000),size=42000/3,replace=TRUE)
train.4=train[idx.4,]
label.4=label[idx.4]
res.KNN.4 =knn(train.4, test, label.4, k = 10, algorithm="cover_tree")
write(res.KNN.4, file="RESKNN4.csv", ncolumns=1)
rm(train.4)
rm(label.4)
rm(res.KNN.4)
rm(idx.4)

idx.5=sample(c(1:42000),size=42000/3,replace=TRUE)
train.5=train[idx.5,]
label.5=label[idx.5]
res.KNN.5 =knn(train.5, test, label.5, k = 10, algorithm="cover_tree")
write(res.KNN.5, file="RESKNN5.csv", ncolumns=1)
rm(train.5)
rm(label.5)
rm(res.KNN.5)
rm(idx.5)

res.knn.1=read.csv("RESKNN1.csv",header=FALSE)
res.knn.2=read.csv("RESKNN2.csv",header=FALSE)
res.knn.3=read.csv("RESKNN3.csv",header=FALSE)
res.knn.4=read.csv("RESKNN4.csv",header=FALSE)
res.knn.5=read.csv("RESKNN5.csv",header=FALSE)
res=data.frame(res.knn.1,res.knn.2,res.knn.3,res.knn.4,res.knn.5)
res=as.matrix(res)
res.knn=vector()
for(i in 1:length(res.knn.1[,1]))
{
res.knn[i]=names(sort(table(as.vector(res[i,])),decreasing=TRUE)[1])
}
res.knn=as.numeric(res.knn)
write(res.knn, file="resknn_submit.csv", ncolumns=1)

#####randomForest#####
library(randomForest)
idx.1=sample(c(1:42000),size=42000/3,replace=TRUE)
train.1=train[idx.1,]
label.1=label.f[idx.1]
model1.randomforest = randomForest(train.1, label.1, xtest=test, ntree=1000)
pred1.randomforest = levels(label.1)[model.randomforest$test$predicted]
write(pred1.randomforest, file="predRandF1.csv", ncolumns=1)
rm(train.1)
rm(label.1)
rm(model1.randomforest)
rm(pred1.randomforest)


idx.2=sample(c(1:42000),size=42000/3,replace=TRUE)
train.2=train[idx.2,]
label.2=label.f[idx.2]
model2.randomforest = randomForest(train.2, label.2, xtest=test, ntree=1000)
pred2.randomforest = levels(label.2)[model.randomforest$test$predicted]
write(pred2.randomforest, file="predRandF2.csv", ncolumns=1)
rm(train.2)
rm(label.2)
rm(model2.randomforest)
rm(pred2.randomforest)

idx.1=sample(c(1:42000),size=42000/3,replace=TRUE)
train.1=train[idx.1,]
label.1=label.f[idx.1]
model1.randomforest = randomForest(train.1, label.1, xtest=test, ntree=1000)
pred1.randomforest = levels(label.1)[model.randomforest$test$predicted]
write(pred1.randomforest, file="predRandF1.csv", ncolumns=1)
rm(train.1)
rm(label.1)
rm(model1.randomforest)
rm(pred1.randomforest)

idx.1=sample(c(1:42000),size=42000/3,replace=TRUE)
train.1=train[idx.1,]
label.1=label.f[idx.1]
model1.randomforest = randomForest(train.1, label.1, xtest=test, ntree=1000)
pred1.randomforest = levels(label.1)[model.randomforest$test$predicted]
write(pred1.randomforest, file="predRandF1.csv", ncolumns=1)
rm(train.1)
rm(label.1)
rm(model1.randomforest)
rm(pred1.randomforest)

idx.1=sample(c(1:42000),size=42000/3,replace=TRUE)
train.1=train[idx.1,]
label.1=label.f[idx.1]
model1.randomforest = randomForest(train.1, label.1, xtest=test, ntree=1000)
pred1.randomforest = levels(label.1)[model.randomforest$test$predicted]
write(pred1.randomforest, file="predRandF1.csv", ncolumns=1)
rm(train.1)
rm(label.1)
rm(model1.randomforest)
rm(pred1.randomforest)

#####naive bayes######
library(e1071)
rm(training)
rm(label)
model.nb=naiveBayes(train,label.f)
pred.nb=predict(model.nb,test)
head(pred.nb) #####lasso#####library(glmnet)
idx.1=sample(c(1:42000),size=42000/3,replace=TRUE)
train.1=train[idx.1,]
label.f.1=label.f[idx.1]
fit.glmnet.1=cv.glmnet(train.1,label.f.1,family="multinomial")
pred.glmnet.1=predict(fit.glmnet.1,newx=test,type="class",s="lambda.min")
#print(fit.glmnet.1)
#plot...