Yang Cheng

行万里路，读万卷书

Pegasos算法

| Tags: Machine Learning, Online Learning

摘要： 本文介绍了svm的一种online learning算法pegasos，并基于pegasoso算法实现了一个数字手写识别脚本。

本文参考了博文Online Learning in Clojure和论文Pegasos: Primal Estimated sub-GrAdient SOlver for SVM(PDF)

online learning

Online learning的算法结构是非常简单的，下面的描述是监督的online learning算法框架，其中有经验损失函数$L$，样本流$S$，样本的格式为$(x,y)$:

Initialise a starting model w
While there are more examples in S
    Get the next feature vector x
    Predict the label y' for x using the model w
    Get the true label y for x and incur a penaly L(y,y')
    Update the model w if y ≠ y'

一般来是，训练出来的模型都是一个与样本相同维度的向量。对应二分的分类器，往往涉及到的是计算内积$\langle w,x \rangle$，模型的更新是沿着损失函数的梯度下降方向的。

Pegasos

论文Pegasos: Primal Estimated sub-GrAdient SOlver for SVM是一种svm的online learning算法。

首先来看svm的经验损失函数：

上面式子中，$k$是训练集$S$的大小，$h()$是the hinge loss（合页损失函数），$\langle w, x\rangle$表示$w,x$的内积，$\lambda$是正则化项。

在《统计学习方法》这本书的7.2.4证明了合页损失函数与引入松弛变量后的损失函数是等价的，并证明了$\lambda$与惩罚系数$C$是成反比的。引入松弛变量后的损失函数为:

训练过程中，如果遇到了一个预测错误的样本$(x,y)$, 对模型的更新方法如下：

其中$t$表示已经训练过的样本个数，$ {w_{t + \frac{1}{2}}}$表示训练过$t$个的样本后的模型，${w_{t + \frac{1}{2} }}$ 表示新模型。 根据pegasos算法，新模型的$l_2$范数如果超出了以 $\frac{1}{ {\sqrt \lambda } }$ 为半径的超球面，那么需要将新模型投射到这个超球面上。即：

为什么需要将新的模型投射到以$\frac{1}{ {\sqrt \lambda } }$为半径的超球面上呢？论文证明了svm的最优解是在下面这个集合中的：

而且在pegasos算法的推导，以及模型初始化$w$的时候，都使用了条件

由上面模型的更新公式可以简单分析一下正则化参数$\lambda$的作用，它决定了训练过程中，后面出现的预测错误的样本，对应模型的修正程度。$\lambda$越大，修正程度越小，$\lambda$越小，修正程度越大。同时$\lambda$与惩罚系数$C$是成反比的，所以也可理解为，在训练过称中，出现预测错误样本时，对模型的惩罚程度。$\lambda$越大，惩罚越小，$\lambda$越小，惩罚越大。

Pegasos的算法描述在论文”Pegasos: Primal Estimated sub-GrAdient SOlver for SVM”也是给出了的，可以参考。

但实际上pegasos是一个线性的svm，而且还是一个没有bias的svm，训练出来的线性函数是$y=\langle w,x \rangle$，在上面的论文中的Extensions小节中也讲到了，目前pegasos还没有证明可应用于线性模型$y=\langle w,x \rangle + b$或者是非线性svm模型。

Pegasos的实现例子

实现了一个基于SMO算法的svm，今天就来基于Pegasos实现数字手写识别。svm用于多分类，还是一对多的方式，手写数据还是来自“Machine Learning in Action”的第二章的数据。下面是实现代码

基于pegasos的数字手写识别 (pegasos.py) download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106

#!/usr/bin/env python
# -*- coding: utf-8 -*-

# Pegasos implemented in Python

import os
import sys
import math

G_WEIGHT = []

def parse_image(path):
    img_map = []
    fp = open(path, "r")
    for line in fp:
        line = line[:-2]
        for ch in line:
            img_map.append(int(ch))
    return img_map

def predict(model, data):
    ret = 0
    for i in range(32*32):
        ret += model[i]*data[i]
    return ret

def train_one_model(data, label, sampleNum, modelNum):
    pvalue = predict(G_WEIGHT[modelNum], data)
    # the hinge loss
    if pvalue * label >= 1:
        return

    # update model
    lambd = 0.5
    new_weight = []
    for i in range(32*32):
        # pegasos
        a = G_WEIGHT[modelNum][i] * ( 1 - 1.0/sampleNum) + (1.0 / (lambd * sampleNum))*label*data[i]
        new_weight.append(a)

    # projection
    norm2 = 0
    for i in range(32*32):
        norm2 += math.pow(new_weight[i], 2)
    norm2 = math.sqrt(norm2)
    if norm2 > (1/math.sqrt(lambd)):
        for i in range(32*32):
            G_WEIGHT[modelNum][i] = new_weight[i]/(norm2 * math.sqrt(lambd))
    else:
        G_WEIGHT[modelNum] = new_weight

def train_one_sample(data, num, sampleNum):
    for modelNum in range(10):
        label = -1
        if num == modelNum:
            label = 1
        train_one_model(data, label, sampleNum, modelNum)

if __name__== "__main__":
    for i in range(10):
        G_WEIGHT.append([])
        for j in range(32 * 32):
            G_WEIGHT[i].append(0)

    dirpath = "./trainingDigits/"
    files = os.listdir(dirpath)
    sampleNum = 0
    for file in files:
        print "training:", file
        data = parse_image(dirpath + file)
        num = int(file[0])
        sampleNum += 1
        train_one_sample(data, num, sampleNum)

    # test
    testdir = "./testDigits/"
    files = os.listdir(testdir)
    right = 0
    wrong = 0
    can_not_classify = 0
    total = 0
    for file in files:
        total += 1
        data = parse_image(testdir + file)
        print "testing:", file
        num = int(file[0])
        classify_failed = True
        for i in range(10):
            pvalue = predict(G_WEIGHT[i], data)
            if pvalue > 0:
                classify_failed = False
                print i, "prdict:", 1
                if i == num:
                    right += 1
                else:
                    wrong += 1
            else:
                print i, "prdict:", -1
        if classify_failed:
            can_not_classify += 1

    print "right=", right
    print "wrong=", wrong
    print "can_not_classify=", can_not_classify
    print "total=", total

训练出来的模型测试结果如下：

right= 849
wrong= 46
can_not_classify= 72
total= 946

一共有946个测试样本，其中46个分类错误，72个没有找到分类，849个正确分类，正确分类率89.7%。$\lambda$取值为0.5。我也没有仔细调整$\lambda$的取值，不过看来结果还是慢不错的。但比起SMO算法实现的svm效果要差一些。但是pegasos的优势是快啊，同样的1934个训练样本，基于SMO的svm，花了3、4个小时训练，而pegasos算法只用了30多秒，逆天了。

实现例子的代码和数据可以在github上下载。pegasos有两个版本，pegasos2.py是pegasos.py的升级版，用了numpy库，使得代码更精简好看，同时运行效率更高。这个目录下还包含了论文的pdf文档Pegasos.pdf。

PS：发现numpy和scipy、matplotlib真是好东西啊，python数学运算离不开。另外发现了一个讲numpy/scipy文档翻译为中文的网站用Python做科学计算，好东西啊。

还发现了一个和机器学习相关的网站http://hunch.net/，有很不多不错的学术方面的东西。