1 import tt  2 import numpy as np  3 import random  4 import time  5 st = time.time()  6   7 def wTx(alpha, trainx, trainy): return np.dot(trainx, np.dot(trainx.T, alpha * trainy))   8   9 def EI(svmst, i): 10     wTtrainx = wTx(svmst.alpha, svmst.trainx, svmst.trainy) 11     return wTtrainx[i] - svmst.trainy[i] 12  13 def UI(svmst, i): 14     wTtrainx = wTx(svmst.alpha, svmst.trainx, svmst.trainy) 15     return float(svmst.trainy[i] * (wTtrainx[i]+ svmst.b) - 1)  16  17 def randomchoosej(i, m): 18     j = i 19     while(j==i): 20         j = random.randint(0, m-1) 21     return j 22  23 def deltacal(v1, v2): return 2.0*np.dot(v1, v2)-np.dot(v1, v1)-np.dot(v2,v2) 24  25 def LH(y1, y2, a1, a2, C): 26     if y1 == y2: L, H = max(0, a1+a2-C), min(a1+a2, C) 27     else: L, H = max(0, a2-a1), min(C+a2-a1, C) 28     return L,H 29  30 def newalpha(L,H,a): 31     if a>H: a = H 32     elif a<L: a = L 33     return a 34  35 def svm(trainx, trainy, C=0.6, toler=0.001, Maxiters=40): 36     st1 = time.time() 37     print "step 2: training..." # step 2: training... 38     iter = 0; p=0; m, n = trainx.shape; alpha = np.zeros((m, 1)); b = 0.0; wTtrainx = wTx(alpha, trainx, trainy) # m*1 wTxi 39     while (iter < Maxiters): 40         alphachange = 0 41         for i in xrange(m): 42             ui = float(trainy[i] * (wTtrainx[i]+ b) - 1) # yi(w.Txi+b) 43             if (alpha[i] < C and ui < -toler) or (alpha[i] > 0 and ui > toler): # choose i 44                 j = randomchoosej(i, m) # randon choose j 45                 alphajo = alpha[j].copy(); alphaio = alpha[i].copy() 46                 L, H = LH(trainy[i], trainy[j], alphaio, alphajo, C)# update alpha or quit 47                 delta = deltacal(trainx[i], trainx[j]) 48                 alpha[j] += trainy[j]*(wTtrainx[j]-trainy[j]-wTtrainx[i]+trainy[i])/delta; alpha[j] = newalpha(L,H,alpha[j]) 49                 if L == H: print "L=H,quit."; continue 50                 if delta >= 0.0: print "delta>=0, quit."; continue 51                 if (abs(alpha[j] - alphajo)[0] < 0.00001):print "j not moving enough, quit." ; continue 52                 alpha[i] += trainy[i] * trainy[j]*(alphajo - alpha[j]) 53                 wTtrainx = wTx(alpha, trainx, trainy) # update b 54                 b1 = trainy[i] - wTtrainx[i]; b2 = trainy[j] - wTtrainx[j] 55                 if alpha[i]>0 and alpha[i]<C: b = b1 56                 elif alpha[j]>0 and alpha[j]<C: b = b2 57                 else: b = (b1+b2)/2.0 58                 alphachange += 1;p+=1 59                 #print "iter: %d i: %d, j: %d pairs changed %d" % (iter, i, j, alphachange) 60         if alphachange == 0: iter += 1 61         #print "iteration number: ", iter 62     print "p",p 63     for i in alpha[alpha>0]: 64         print "%.8f" % i 65     print ‘b‘, b 66     return alpha,b,(time.time()-st1) 67  68 def suportvector(alpha, trainx, trainy): 69     return trainx[np.nonzero(alpha>0)], trainy[np.nonzero(alpha>0)] 70  71 class svmsture(object): 72     def __init__(self, trainx, trainy, C, toler): 73         self.trainx = trainx 74         self.m = len(trainx) 75         self.trainy = trainy 76         self.C = C 77         self.toler = toler 78         self.alpha = np.zeros((self.m,1)) 79         self.b = 0.0 80         self.ek = np.zeros((self.m,2)) 81  82 def choosej(svmst, i, ei): 83     svmst.ek[i] = [1,ei] 84     nzeiid = np.nonzero(svmst.ek)[0] 85     if len(nzeiid) > 1: 86         maxid = 0; maxediff = 0.0 87         for id in nzeiid: 88             if id == i: continue 89             eiid = EI(svmst, id) 90             absdiff = abs(eiid - ei) 91             if absdiff > maxediff: 92                 maxid = id; maxediff = absdiff 93     else: maxid = randomchoosej(i, svmst.m) 94     return maxid 95  96 def update(svmst, i): 97     ei = EI(svmst, i) 98     svmst.ek[i] = [1, ei] 99 100 def svminner(svmst, i):101     ei = EI(svmst, i); ui = UI(svmst, i); wTtrainx = wTx(svmst.alpha, svmst.trainx, svmst.trainy)102     if (svmst.alpha[i] < svmst.C and ui < -svmst.toler) or (svmst.alpha[i] > 0 and ui > svmst.toler): # choose i103         j = choosej(svmst, i, ei)104         alphajo = svmst.alpha[j].copy(); alphaio = svmst.alpha[i].copy()105         L, H = LH(svmst.trainy[i], svmst.trainy[j], alphaio, alphajo, svmst.C)106         delta = deltacal(svmst.trainx[i], svmst.trainx[j])107         if L == H: print "L=H,quit."; return 0108         if delta >= 0.0: print "delta>=0, quit."; return 0109         svmst.alpha[j] += svmst.trainy[j]*(wTtrainx[j]-svmst.trainy[j]-wTtrainx[i]+svmst.trainy[i])/delta; 110         svmst.alpha[j] = newalpha(L,H,svmst.alpha[j]); update(svmst, j)111         if (abs(svmst.alpha[j] - alphajo)[0] < 0.00001):print "j not moving enough, quit." ; return 0112         svmst.alpha[i] += svmst.trainy[i] * svmst.trainy[j]*(alphajo - svmst.alpha[j])113         update(svmst, i); wTtrainx = wTx(svmst.alpha, svmst.trainx, svmst.trainy)114         b1 = svmst.trainy[i] - wTtrainx[i]; b2 = svmst.trainy[j] - wTtrainx[j]115         if svmst.alpha[i]>0 and svmst.alpha[i]<svmst.C: svmst.b = b1116         elif svmst.alpha[j]>0 and svmst.alpha[j]<svmst.C: svmst.b = b2117         else: svmst.b = (b1+b2)/2.0 118         return 1119     else: return 0120 121 def svmouter(trainx, trainy, C = 0.6, toler = 0.001, Maxiters = 40):122     svmst = svmsture(trainx, trainy, C, toler)123     alphachange = 0; iter = 0124     while (alphachange == 0 and iter < Maxiters):125         for i in xrange(svmst.m):126             alphachange += svminner(svmst, i)127         print "FullSet, iter: %d, pairs changed %d" % (iter, alphachange)128         iter += 1; alphachange = 1129         while (alphachange > 0):130             nonboundid = np.nonzero((svmst.alpha>0) * (svmst.alpha<svmst.C))[0]; alphachange = 0 ###131             for i in nonboundid:132                 alphachange += svminner(svmst, i)133             print "non-bound, iter: %d, pairs changed %d" % (iter, alphachange)134             iter += 1135         print "iteration number: ", iter136     return svmst.b, svmst.alpha137 138 def svmclassify(testx):139     trainx, trainy = tt.loaddata(‘testSet‘)140     b, alpha = svmouter(trainx, trainy)141     wTtrainx = wTx(alpha, trainx, trainy)142     testy = np.dot(testx, np.dot(trainx.T, alpha * trainy)) + b143     print ‘b‘,b144     print ‘alpha‘, alpha[alpha>0]145     print ‘w‘,np.dot(trainx.T, alpha * trainy)146     print ‘testy=wTtest + b‘,testy147             148 trainx, trainy = tt.loaddata(‘testSet‘)    149 i=0150 svmclassify(trainx[i])151 print ‘trainy[%d] %d‘ % (i, trainy[i])152 """153 svmst = svmsture(trainx, trainy, 0.6, 0.001)154 for i in xrange(svmst.m):155     print i156     print "svminner",svminner(svmst, i)157     print "type(svm)",type(svminner(svmst, i))158 """159 160 161 162 print "cost time: ", time.time()-st