在周二我给精算师上的5小时机器学习速成课结束时,皮埃尔问了我一个有趣问题,是关于不同技术的计算时间的。我一直在介绍各种算法的思想,却忘了提及计算时间。我想在数据集上尝试几种分类算法来阐述这些技术。
> rm(list=ls())
> myocarde=read.table(
"http://freakonometrics.free.fr/myocarde.csv",
head=TRUE,sep=";")
> levels(myocarde$PRONO)=c("Death","Survival")
数据集相当小,包括71个观测值和7个解释变量。因此我决定复制观测结果,并添加一些协变量,
> levels(myocarde$PRONO)=c("Death","Survival")
> idx=rep(1:nrow(myocarde),each=100)
> TPS=matrix(NA,30,10)
> myocarde_large=myocarde[idx,]
> k=23
> M=data.frame(matrix(rnorm(k*
+ nrow(myocarde_large)),nrow(myocarde_large),k))
> names(M)=paste("X",1:k,sep="")
> myocarde_large=cbind(myocarde_large,M)
> dim(myocarde_large)
[1] 7100 31
> object.size(myocarde_large)
2049.064 kbytes
数据集虽然不大...但去跑一个回归至少会大于0.0001秒。实际上,跑逻辑回归用了0.1秒
> system.time(fit< glm(PRONO~.,
+ data=myocarde_large, family="binomial"))
user system elapsed
0.114 0.016 0.134
> object.size(fit)
9,313.600 kbytes
令我吃惊的是回归对象的大小是9M,是数据集大小的四倍多。使用大数据集,大小要大100倍,
> dim(myocarde_large_2)
[1] 710000 31
这花了20秒。
> system.time(fit<-glm(PRONO~.,
+ data=myocarde_large_2, family="binomial"))
utilisateur système écoulé
16.394 2.576 19.819
> object.size(fit)
90,9025.600 kbytes
而这个对象的大小‘只是’之前的十倍。
注意到对于样条函数,计算时间也很相似
> library(splines)
> system.time(fit<-glm(PRONO~bs(INSYS)+.,
+ data=myocarde_large, family="binomial"))
user system elapsed
0.142 0.000 0.143
> object.size(fit)
9663.856 kbytes
如果我们换用另一个函数,更确切地说是用多项式回归,它的耗时会是之前的两倍
> library(VGAM)
> system.time(fit1<-vglm(PRONO~.,
+ data=myocarde_large, family="multinomial"))
user system elapsed
0.200 0.020 0.226
> object.size(fit1)
6569.464 kbytes
然而对象的大小较小。现在,如果我们使用步进式程序,它耗时会长一些:几乎要一分钟,比单次逻辑回归慢500倍
> system.time(fit<-step(glm(PRONO~.,data=myocarde_large,
family="binomial")))
...
Step: AIC=4118.15
PRONO ~ FRCAR + INCAR + INSYS + PRDIA + PVENT + REPUL + X16
Df Deviance AIC
<none> 4102.2 4118.2
- X16 1 4104.6 4118.6
- PRDIA 1 4113.4 4127.4
- INCAR 1 4188.4 4202.4
- REPUL 1 4203.9 4217.9
- PVENT 1 4215.5 4229.5
- FRCAR 1 4254.1 4268.1
- INSYS 1 4286.8 4300.8
user system elapsed
50.327 0.050 50.368
> object.size(fit)
6,652.160 kbytes
我也想尝试caret,这个软件包很适合用来对比模型。在JRSS-A 计算精算科学(R语言)这本书的解读中,Andrey Kosteko注意到这个软件包甚至没有被提及,相关内容也是空白的。所以我尝试了逻辑回归
> library(caret)
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="glm"))
user system elapsed
5.908 0.032 5.954
> object.size(fit)
12,676.944 kbytes
它花费了6秒(是标准调用glm函数的50倍),并且对象大小也大些。而且如果我们尝试步进式调用,那会更糟糕
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="glmStepAIC"))
...
Step: AIC=4118.15
.outcome ~ FRCAR + INCAR + INSYS + PRDIA + PVENT + REPUL + X16
Df Deviance AIC
<none> 4102.2 4118.2
- X16 1 4104.6 4118.6
- PRDIA 1 4113.4 4127.4
- INCAR 1 4188.4 4202.4
- REPUL 1 4203.9 4217.9
- PVENT 1 4215.5 4229.5
- FRCAR 1 4254.1 4268.1
- INSYS 1 4286.8 4300.8
user system elapsed
1063.399 2.926 1068.060
> object.size(fit)
9,978.808 kbytes
花了15分钟,还是在只有30个协变量的情况下......过程如下(我是用microbenchmark来绘制的)
让我们尝试下树。
> library(rpart)
> system.time(fit<-rpart(PRONO~.,
+ data=myocarde_large))
user system elapsed
0.341 0.000 0.345
> object.size(fit4)
544.664 kbytes
跑起来很快,而且大小较小。如果我们改变复杂度参数去得到更深的树,这前后几乎是一样的
> system.time(fit<-rpart(PRONO~.,
+ data=myocarde_large,cp=.001))
user system elapsed
0.346 0.000 0.346
> object.size(fit)
544.824 kbytes
但是同样的,如果我们通过caret调用相同的函数,速度会慢十倍以上,
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="rpart"))
user system elapsed
4.076 0.005 4.077
> object.size(fit)
5,587.288 kbytes
而对象大小是之前的十倍。现在考虑一下随机森林。
> library(randomForest)
> system.time(fit<-randomForest(PRONO~.,
+ data=myocarde_large,ntree=50))
user system elapsed
0.672 0.000 0.671
> object.size(fit)
1,751.528 kbytes
在“只有”50棵树的情况下,它只用了两倍的时间就跑出了结果。但如果是500棵树(默认值)就需要20多倍的时间(从比例上看这也是合理的时间,创建了500棵树而不是50)
> system.time(fit<-randomForest(PRONO~.,
+ data=myocarde_large,ntree=500))
user system elapsed
6.644 0.180 6.821
> object.size(fit)
5,133.928 kbytes
如果在每个节点上改变所使用的协变量的个数,我们可以看到几乎没有影响。用5个协变量(协变量总数的平方根,即默认值),需要6秒,
> system.time(fit<-randomForest(PRONO~.,
+ data=myocarde_large,mtry=5))
user system elapsed
6.266 0.076 6.338
> object.size(fit)
5,161.928 kbytes
但如果我们使用10,它几乎是相同的(甚至更少)
> system.time(fit<-randomForest(PRONO~.,
+ data=myocarde_large,mtry=10))
user system elapsed
5.666 0.076 5.737
> object.size(fit)
2,501.928 bytes
如果我们在caret中使用随机森林算法,则需要10分钟,
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="rf"))
user system elapsed
609.790 2.111 613.515
可视化过程如下
如果我们考虑caret中使用K近邻算法,也需要10分钟(译者注:此处应该是1分钟)
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="knn"))
user system elapsed
66.994 0.088 67.327
> object.size(fit)
5,660.696 kbytes
这与在树上(译者注:决策树)应用bagging算法耗时几乎是相同的
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="treebag"))
Le chargement a nécessité le package : plyr
user system elapsed
60.526 0.567 61.641
> object.size(fit)
72,048.480 kbytes
但是这一次,对象却相当大!
我们也可以考虑应用标准欧几里德距离的SVM算法
> library(kernlab)
> system.time(fit<-ksvm(PRONO~.,
+ data=myocarde_large,
+ prob.model=TRUE, kernel="vanilladot"))
Setting default kernel parameters
user system elapsed
14.471 0.076 14.698
> object.size(fit)
801.120 kbytes
或者设置一些内核参数
> system.time(fit<-ksvm(PRONO~.,
+ data=myocarde_large,
+ prob.model=TRUE, kernel="rbfdot"))
user system elapsed
9.469 0.052 9.701
> object.size(fit)
846.824 kbytes
这两种技术都需要10秒左右,远远超过基本的逻辑回归模型(慢100倍以上)。同样的,如果我们用caret跑,那就需要一段时间了......
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large, method="svmRadial"))
user system elapsed
360.421 2.007 364.669
> object.size(fit)
4,027.880 kbytes
输出如下
我也想尝试一些算法,比如ridge和LASSO。
> library(glmnet)
> idx=which(names(myocarde_large)=="PRONO")
> y=myocarde_large[,idx]
> x=as.matrix(myocarde_large[,-idx])
> system.time(fit<-glmnet(x,y,alpha=0,lambda=.05,
+ family="binomial"))
user system elapsed
0.013 0.000 0.052
> system.time(fit<-glmnet(x,y,alpha=1,lambda=.05,
+ family="binomial"))
user system elapsed
0.014 0.000 0.013
我惊讶地发现它跑得很快。如果我们使用交叉验证来量化耗时
> system.time(fit10<-cv.glmnet(x,y,alpha=1,
+ type="auc",nlambda=100,
+ family="binomial"))
user system elapsed
11.831 0.000 11.831
这需要一些时间......但与其他技术相比,这是合理的。最后,考虑使用boosting包。
> system.time(fit<-gbm.step(data=myocarde_large,
+ gbm.x = (1:(ncol(myocarde_large)-1))[-idx],
+ gbm.y = ncol(myocarde_large),
+ family = "bernoulli", tree.complexity = 5,
+ learning.rate = 0.01, bag.fraction = 0.5))
user system elapsed
364.784 0.428 365.755
> object.size(fit)
8,607.048 kbytes
很长。超过6分钟。通过caret使用glmboost包快得多
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="glmboost"))
user system elapsed
13.573 0.024 13.592
> object.size(fit)
6,717.400 bytes
通过caret调用gbm的时间是之前的十倍多,
> system.time(fit<-train(PRONO~.,
+ data=myocarde_large,method="gbm"))
user system elapsed
121.739 0.360 122.466
> object.size(fit)
7,115.512 kbytes
所有这些都是笔记本电脑完成的。现在我需要在更快的机器上运行相同的代码,来尝试更大的数据集......