ESL-chapter10-boosting

boosting是一个很强大的方法。现介绍如下：

     假定有一个比较弱的分类器，如一层的决策树。分类的结果只比随机猜测要好一些。boosting 就是不断的用这个分类器去拟合加权过的原始数据，得到多个分类器，最后取多个分类器结果的加权和。注意，这里的关键是 “加权过的原始数据”，在boosting的每一轮，我们都要对原始数据进行加权处理，给每个数据赋值一个权值，权值修改的原则是：误分类的数据权值逐渐加大。具体的算法可见书上的算法10.1。

     现在按照书上的示例来生成数据，测试boosting的威力。采用rpart来生成弱分类器decision stump。这里稍微插一句，R里面生成决策树根的方法有四种，tree，ctree, rpart和decisionstump。前三种是通过调整参数来生成decision stump，decisionstump函数是直接生成，但是decisionstump函数不支持加权，所以不好用。前三种我拟合了很多遍，只有rpart才能够生成我想要的decision stump，tree和ctree都不好用（主要是参数不好调整）。最后生成了书上的图10.2.

这个图没有书上的图光滑，但是大概能够说明趋势。

library(rpart)#引入rpart包
x.train <- matrix(rnorm(2000*10,0,1),2000,10)
y.train <- ifelse(rowSums(x.train^2)>9.34,1,-1)
y.train <- as.factor(y.train)
x.test <- matrix(rnorm(10000*10,0,1),10000,10)
y.test <- ifelse(rowSums(x.test^2)>9.34,1,-1)
y.test <- as.factor(y.test)
X.train <- data.frame(x.train,y.train)
X.test <- data.frame(x.test,y.test)#按照书上的要求生成数据，测试数据实际上没有用。


AdaBoost.M1 <- function(x,y,M,...)#
    {
        N <- length(y)
        w <- rep(1/N,N)#每轮迭代的权值
        alpha <- rep(0,M)
        G <- c()#初始化一些必要的变量，G用来加权最后生成的模型
      
        c.error <- c()
        for(m in 1:M)
            {
                 rpart.train <- rpart(y~.,data=data.frame(x,y),weights=w,control=rpart.control(minsplit = 2,minbucket=0,maxdepth=1))
                 rpart.pred <- predict(rpart.train,newdata=data.frame(x),type="class")#用rpart拟合数据，并预测
                indicator <-  ifelse(rpart.pred!=y,1,0)#接下来都是对书上算法的每一步实现。注释掉的地方是别的书上给出的算法，我实现并比较了一下。
                error <- sum(w*indicator)/sum(w)
                 alpha[m] <-log((1-error)/error)
               # alpha[m] <-0.5*log((1-error)/error)
                #in.y <- as.numeric(as.character(y))*as.numeric(as.character(rpart.pred))
                #Z <- sum(w*exp(-1*alpha[m]*in.y))
                #w <- w*exp(-1*alpha[m]*in.y)/Z
                w <- w*exp(alpha[m]*indicator)
                G <- cbind(G,alpha[m]*as.numeric(as.character(rpart.pred)))
                
                labels <- sign(rowSums(G))
                c.error <- c(c.error,(table(labels,y)[2]+table(labels,y)[3])/N)
             }
       
        return(c.error)
    }


error.collection <- AdaBoost.M1(x.train,y.train,400)#最后，画图
plot(1:length(error.collection),error.collection,col="orange",ylim=c(0,0.5),type="l")

我们同样绘制了了图10.3，但是和书上的结果不太一样。指数损失函数的值大很多，原因待查~