遗传算法--背包问题

背包问题描述：
有n个物品，每个物品有2个属性：重量w和价值v。指定一个背包，该背包能承受的最大重量为W，
对每个物品来说，只有放入背包和不放入背包2种选择。
求解：如何在不超过该背包最大承重W的前提下，使得放入物品的总价值V最大？

曾经用动态规划和分支限界搞定过，换个方法玩一下。
达尔文说过：能够生存下来的往往不是最强大的物种，也不是最聪明的物种，而是最能适应环境的物种。
这个算法就是这么来的，物竞天择，适者生存。
那么这个算法就有四个基本的要点：
种群：由个体组成，每个个体携带有染色体和基因。
选择：优胜劣汰。
交叉：通过交叉来进行基因的重组。
变异：基因概率突变。

1.初始化：
那么背包问题如何与生物这个遗传进化的过程相联系呢。
我们将染色体用二进制，那么0表示物品不装入背包，1表示物品装入背包。
在进化过程当中，0即为基因丢失，1即为基因保留。
个体的适应性用适应度衡量，背包问题中也就是我们背包中的价值总和越高越好，但是这里的总重量
不可以超过背包的预设重量，进化不可逆自然规律而行。
2.选择：
从总体中选择合适的个体，让其进行交配。
这里为了保证生物的多样性，不至于在几代之后差异过小，那么选用轮赌法。
也就是相当于将一个圆分为总体中染色体数量的n个部分，每条染色体所占圆的面积由其适应度按比例来
界定。
那么为这个圆设置一个指针，轮盘旋转，当轮盘停下的时候指针所指即为所选择个体，
那么重复这个过程两次，即可选中两个亲本。
3.交叉：
交叉分为很多中，单点交叉，多点交叉，均匀交叉，这里选择单点交叉。
随机在染色体上选择一个交叉点，然后交换交叉点前后的部分，形成两个新的个体。
表现在这里，就是我们的二进制序列上选择一个分界点，亲本的前半部分和后半部分分别重新组合即可。
4.变异：
生物在进化过程当中，不仅有遗传，还有变异。
体现在我们这里就是二进制序列上的某一位概率的发生改变。

这个过程进行完，那么就基本完成了一个进化的过程，也就是一轮繁殖。
那么我们什么时候找到最好的那一代呢，也就是进化最完全的那一代呢。
也就是如何找到最佳的适应水平。
一般有如下几种适应条件：
1.再进行繁殖不会发生太大的改变，趋于稳定。
2.事先规定迭代次数。
3.事先预定最佳的适应度。

以下即为模拟过程代码：

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
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
PACK_MAX_CAP=1000
MAX_GENERATION=100
import random
def getRand():
    intnum=random.randint(0,65535)
    return intnum
def getRand_float():
    floatnum=random.uniform(0,1)
    return floatnum
class good:
    def __init__(self,weight,value):
        self.weight=weight
        self.value=value
class entity:
    def __init__(self,goods_num):
        self.fitness = 0
        self.sum_weight = 0
        self.generation_id = 0
        self.gene = []
        self.initgene(goods_num)
    def initgene(self,goods_num):
        for i in range(goods_num):
            self.gene.append(0)
class gene_algorithm:
    def __init__(self,grnum,cr,vr,gonum):
        self.group_num=grnum
        self.goods_num=gonum
        self.cross_rate=cr
        self.varia_rate=vr
        self.best_entity = None
        self.group = []  # 群体
        self.goods = []  # 物品
        self.initgoods()
    def get_totalfitness(self,group):
        sum=0.0
        for i in range(len(group)):
            sum+=group[i].fitness
        return sum
    def print_opt(self):
        print("最优解所在代数"+str(self.best_entity.generation_id)+"\n")
        print("总重量"+str(self.best_entity.sum_weight)+"\n")
        print("总价值"+str(self.best_entity.fitness)+"\n")
        print("装入情况为：\n")
        for i in range(len(self.best_entity.gene)):
            print(self.best_entity.gene[i])
    def initgoods(self):
        print("物品集合信息：")
        for i in range(self.goods_num):
            self.goods.append(good(getRand()%150,getRand()%150))
            print("物品号："+str(i)+"重量："+str(self.goods[i].weight)+"价值："+str(self.goods[i].value))
    def initgroup(self):
        self.best_entity=entity(self.goods_num)
        for i in range(self.group_num):
            ent=entity(self.goods_num)
            #基因全零个体；
            vir_capacity=PACK_MAX_CAP*random.uniform(0.5,1)
            weight=0
            count=0
            while(weight<=vir_capacity):
                idx=getRand()%self.goods_num
                #随机取得每个染色体上某个基因的编号；
                if count==3:
                    #连续三次随机到同一个基因位，该个体初始化结束；
                    break;
                if ent.gene[idx]==1:
                    count+=1
                    continue
                ent.gene[idx]=1
                weight+=self.goods[idx].weight
            self.group.append(ent)
            #初始化后加入种群；
    def cal_fitness(self):
        for i in range(len(self.group)):
            ent=self.group[i]
            wei=0;val=0
            for j in range(len(ent.gene)):
                if ent.gene[j]==1:
                    wei+=self.goods[j].weight
                    val+=self.goods[j].value
            if wei>PACK_MAX_CAP:
                ent.fitness=0
                continue
            ent.fitness=val
            ent.sum_weight=wei
    def recordoptimalentity(self,gid):
        sorted(self.group,key=lambda entity:entity.fitness)#通过fitness排序；
        ent=self.group[len(self.group)-1]
        if ent.fitness>self.best_entity.fitness:
            self.best_entity=ent
            self.best_entity.generation_id=gid
        print("遗传代数"+str(gid)+"总重量："
              +str(ent.sum_weight)+" 总价值:"+str(ent.fitness)+"\n")
    def select(self):
        new_group=[]
        selected_rate=[]
        sorted(self.group, key=lambda entity: entity.fitness)
        group_num=len(self.group)
        sum_fitness=self.get_totalfitness(self.group)
        selected_rate.append(self.group[0].fitness/sum_fitness)
        for i in range(1,group_num):
            selected_rate.append(selected_rate[i-1]+self.group[i].fitness/sum_fitness)
            #构造轮盘；
        left_group_num=group_num*0.5
        for i in range(int(left_group_num)):
            rand_rate=getRand_float()
            for idx in range(len(selected_rate)):
                if rand_rate<=selected_rate[idx]:
                    new_group.append(self.group[idx])
                    break;
        #更新群体
        self.group.clear()
        self.group=new_group
    def Willcross(self):
        return (getRand_float()<=self.cross_rate)
    def cross(self):
        group_num=len(self.group)
        for i in range(0,(group_num-1),2):
            father=self.group[i]
            mother=self.group[i+1]
            for j in range(len(father.gene)):
                if self.Willcross():
                    temp=father.gene[j]
                    father.gene[j]=mother.gene[j]
                    mother.gene[j]=temp
            self.group.append(father)
            self.group.append(mother)
    def Willvaria(self):
        return (getRand_float()<=self.varia_rate)
    def varia(self):
        group_num=len(self.group)
        for i in range(group_num):
            if self.Willvaria():
                ent=self.group[i]
                for j in range(len(ent.gene)):
                    if self.Willvaria():
                        if ent.gene[j]==1:
                            ent.gene[j]==0
                        else:
                            ent.gene[j]==1
    def run(self):
        self.initgroup()
        for i in range(MAX_GENERATION):
            self.cal_fitness()
            self.recordoptimalentity(i)
            self.select()
            self.cross()
            if i%5==0&i!=0:
                self.varia()
if __name__=='__main__':
    ga=gene_algorithm(20,0.8,0.5,40)
    ga.run()
    ga.print_opt()

展开全文 >>

tp中CRUD add()方法添加数据主键自增，插入为空解决

thinkphp中操作数据库使用add()方法添加数据，然后查询，发现只有主键，
其他字段值并没有返回，然后去看了一下，发现字段值确实为空，说明add执行失败了。
经过一番请教，原来是缓存的原因，

如果在使用表过后，又手动更改了表结构，Runtime中缓存了表的结构信息，
然后手动修改表，并未与Runtime一致，就会导致这个原因，
只要将缓存清空即可，也就是将Runtime文件清空，问题就会解决。
我们也可以编写自己的脚本，来自动化这个Runtime清空过程。

展开全文 >>

轮廓检测入门

人脸检测完成后，我们需要将检测到的人脸进行框选，那么就需要轮廓检测了。
不讨论过多的算法方面的东西(因为实在不太懂)，我这次直接使用了opencv中的
findContours()函数。

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
def getcontours(img):
    best_rate=2
    best_w=0
    best_h=0
    best_x=0
    best_y=0
    rate=1
    height=img.shape[0]
    width=img.shape[1]
    min_area=height*width/1000
    want_op=img.copy()
    want_op=cv2.cvtColor(want_op,cv2.COLOR_BGR2GRAY)
    ret,bina=cv2.threshold(want_op,0,255,cv2.THRESH_BINARY)
    _,contours,her=cv2.findContours(bina,cv2.RETR_TREE,cv2.CHAIN_APPROX_SIMPLE)
    #连通区域，比例最接近，坐标限制；
    for contour in contours:
        area=cv2.contourArea(contour)
        x, y, w, h = cv2.boundingRect(contour)
        if x>width/5 and x<width*0.8 and y<height/4:
            if area > min_area:
                my_rate=w/h
                if abs(my_rate-rate)<abs(best_rate-rate):
                    best_rate=my_rate
                    best_h=h
                    best_w=w
                    best_x=x
                    best_y=y
    return best_x,best_y,best_w,best_h

直接说这段代码吧,我们会在其中看到这么一句:

1
_,contours,her=cv2.findContours(bina,cv2.RETR_TREE,cv2.CHAIN_APPROX_SIMPLE)

这个就是opencv内置的轮廓检测的函数，
这里一定要注意，在网上看到的大多数教程这个函数的返回值都是两个，
但是由于版本问题，我使用的其实返回了三个值,即使第一个值很少用，
但是如果不接就会报错: ValueError: too many values to unpack
返回值:
(1)基本用不到，处理的图像；
(2)轮廓的点集，在opencv中，轮廓就是一系列点的集合；
(3)各层轮廓的索引；
主要参数:
(1)二值图，一定是二值图，所以我们需要对图像进行预处理。
我这里其实之前传进去需要处理的就是一副黑白图，但是不知道为什么还会报错，所以
就把黑白图依然进行下述处理。

1
2
want_op=cv2.cvtColor(want_op,cv2.COLOR_BGR2GRAY)
ret,bina=cv2.threshold(want_op,0,255,cv2.THRESH_BINARY)

这个处理过程是原地处理，所以如果后续需要原图，那么最好copy()一下。

(2)第二个参数是轮廓的检索方式
CV_RETR_EXTERNAL:
表示只检测外轮廓；
CV_RETR_LIST:
检测的轮廓不建立等级关系；
CV_RETR_CCOMP:
建立两个等级的轮廓，上面的一层为外边界，里面的一层为内孔的边界信息。
如果内孔内还有一个连通物体，这个物体的边界也在顶层；
CV_RETR_TREE:
建立一个等级树结构的轮廓；

(3)第三个参数是轮廓的近似方法
CV_CHAIN_APPROX_NONE:
存储所有的轮廓点，相邻的两个点的像素位置差不超过1;
CV_CHAIN_APPROX_SIMPLE:
压缩水平方向，垂直方向，对角线方向的元素，只保留该方向的终点坐标，
例如一个矩形轮廓只需4个点来保存轮廓信息;
CV_CHAIN_APPROX_TC89_L1，CV_CHAIN_APPROX_TC89_KCOS使用teh-Chinl chain 近似算法。

获取轮廓后，我们只要对返回的点集，也就是第二个参数进行遍历，每次遍历得到一个contour，
如果想要对其进行处理，比如标记等等，或者是在其中选出需要的轮廓(比如最大的轮廓)。
在循环中处理逻辑即可。

下面的代码是对获取到的轮廓区域内的像素点进行统计投票，
如果满足条件的像素点数大于总像素点数的八分之三，返回0，否则返回1。
通过分析，还发现rectangle函数的框选区域把标记框也标记进去了，
而且轮廓线上左宽为2像素，下右宽为1像素,需要精确内部区域的一定要把框所占数量去掉。

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
def has_cap(bin_img,img):
    red_cap=[[110,220],[0,100],[0,100]]
    yellow_cap=[[205,255],[110,220],[0,100]]
    white_cap=[[205,255],[205,255],[205,255]]
    x,y,w,h=getcontours(bin_img)
    cv2.rectangle(img, (x, y), (x + w, y + h), (255, 0, 0), 2)
    area=(w-3)*(h-3)
    if w-3<=0 or h-3<=0:
        return 2
    #检测不到人头
    #轮廓线上左宽为2，下右宽为1;
    is_cap=np.zeros([h-3,w-3],np.int)
    print(str(h)+','+str(w))
    for i in range(y+2,y+h-1):
        print(' ')
        for j in range(x+2,x+w-1):
            bgr=img[i,j]
            print(bgr)
            if bgr[0]>=0 and bgr[0]<=100:
                if bgr[1]>=0 and bgr[1]<=100 and bgr[2]>=110 and bgr[2]<=220:
                    is_cap[i-y-2,j-x-2]=1
                elif bgr[1]>=110 and bgr[1]<=220 and bgr[2]>=205 and bgr[2]<=255:
                    is_cap[i-y-2,j-x-2]=1
            elif bgr[0]>=205 and bgr[0]<=255:
                if bgr[1]>=205 and bgr[1]<=255 and bgr[2]>=205 and bgr[2]<=255:
                    is_cap[i-y-2,j-x-2]=1
    area2=0
    for i in range(is_cap.shape[0]):
        for j in range(is_cap.shape[1]):
            area2+=is_cap[i,j]
    print(str(area)+','+str(area2))
    if area2/area>3/8:
        return 0
    else:
        return 1

展开全文 >>

肤色检测

这几天做安全帽检测器，需要定位人脸的部位，就需要用到皮肤检测。
在经过行人检测过后，定位图像中人的位置以后，标记出人的位置，然后将标记过后的行人框当作
皮肤检测的目标范围。
在网上学习了一下，目前有很多检测的技术，这篇文章里做个总结。
1.基于RGB通道检测
在检测范围内对每个像素点通过RGB值来进行判断:
经过大量研究，肤色在RGB模型中的范围满足一下判别式：

1
2
3
4
#均匀光照中：
R>95 AND G>40 B>20 AND MAX(R,G,B)-MIN(R,G,B)>15 AND ABS(R-G)>15 AND R>G AND R>B
#侧光情况中：
R>220 AND G>210 AND B>170 AND ABS(R-G)<=15 AND R>B AND G>B

只要在像素通道值判断中通过上述表达式进行检测即可。

2.基于YCrCb空间椭圆模型检测
同样是经过大量的研究，将皮肤模型映射到YCrCb空间中，那么在CrCb的二维空间当中，满足肤色的
值将会落在一个椭圆的区域内:

1
2
#椭圆参数如下：
ellipse=cv2.ellipse(skincbcr,[113,155.6],[23.4,15.6],43.0,0.0,360.0,[255,255,255],-1)

那么通过代入(Cr,Cb)坐标的然后判断是否落在椭圆区域内即可。

3.YCrCb空间中Cr分量+Otus阈值法
(1)YCrCb空间:
Y即为明亮度，CrCb即为色度，分别为色调与饱和度。
更多的关于YCrCb的信息网上都有，就不具体介绍了。
(2)最大类间方差法(Otus):
使用一个阈值将整个数据空间分为两类，通过两个类之间的方差值来获得最优的阈值。
我们需要做的就是将RGB空间映射到YCrCb空间，然后提取Cr分量;
然后对Cr进行二值化阈值分割,利用上面Otus法判断即可。

4.基于YCrCb空间CrCb分量范围筛选
类似法一，也是通过前人的大量统计，
黄种人的Cr范围在133和173之间，Cb分量在77和127之间。
然后根据项目的实际情况灵活调控这个范围在进行范围检测即可。

5.HSV空间H分量范围筛选
HSV空间：分别为H色调，S饱和度，V明度。根据颜色的直观特性建立的色彩空间模型，
也被称为六角锥模型。


(1)色调H：用角度度量，取值范围为0°～360°，
从红色开始按逆时针方向计算，红色为0°，绿色为120°,蓝色为240°。
它们的补色是：黄色为60°，青色为180°,品红为300°；
(2)饱和度S：取值范围为0.0～1.0；
(3)亮度V：取值范围为0.0(黑色)～1.0(白色)。
满足皮肤模型的范围:

1
2
(H>=0 AND H<=25) OR (H>=335 AND H<=360)
(S>=0.2 AND S<=0.6) AND (V>=0.4)

通过上面的范围判断即可。

6.基于YCrCb空间高斯模型检测
高斯模型主要利用了统计学原理，肤色这样符合正态分布的统计样本，在特征空间中的分布应当
满足高斯分布。
关于高斯分布不做过多的介绍。
(1)我们要做的就是在YCrCb空间当中，
得到Cr和Cb分量的值，通过样本得到均值Mean(MCb,MCr)和协方差矩阵C(这是一个
统计过程，所以说这种检测方法和机器学习结合可能实际效果好一些吧)。
(2)将均值和协方差矩阵代入二维高斯函数:

1
2
P(Cb,Cr)=exp(-0.5*(x-Mean)*inv(C)*(x-Mean).T)
#其中x=(Cb,Cr)，是一个1*2的矩阵；

下面是一个python实现的版本。

1
2
3
4
5
6
7
8
9
cbcr = np.array([[Cb], [Cr]])
#得到一个2*1的矩阵;
FaceProbImg[y, x] =
    np.exp(-0.5 * ((cbcr - Mean).T * np.linalg.inv(C) * (cbcr - Mean))[0, 0])
    #得到的是一个1*1的矩阵，取[0,0]即取出其值;
    if FaceProbImg[y, x] > threshhold:
        img_bin[y, x] = [255,255,255]
    else:
        img_bin[y,x]=[0,0,0]

上面得到的FaceProbImg是像素属于皮肤区域的概率，我们设置一个阈值threshold,
比较得到的概率和阈值即可。

经过上面的对比，最后选择了椭圆模型，感觉效果最好,相比来说实现也比较简单。
关键代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
img_ycbcr=img.copy()
    img_bin=img.copy()
    img_ycbcr=cv2.cvtColor(img_ycbcr,cv2.COLOR_BGR2YCrCb)
    height=img_ycbcr.shape[0]
    width=img_ycbcr.shape[1]
    for y in range(height):
        for x in range(width):
            Y=img_ycbcr[y,x,0]
            CrY=img_ycbcr[y,x,1]
            CbY=img_ycbcr[y,x,2]
            Cb,Cr=Before_Model(Y,CbY,CrY)
            #椭圆模型；
            ellipse=ellipse2()
            cbcr=np.array([[Cb,Cr]])
            if ellipse[cbcr[0,0]]>0:
                img_bin[y, x] = [255,255,255]
            else:
                img_bin[y,x]=[0,0,0]

上面这段代码可以更换结构，来使用不同的模型进行检测。

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
Kl, Kh = 125, 188
Ymin, Ymax = 16, 235
Wlcb, Wlcr = 23, 20
Whcb, Whcr = 14, 10
Wcb, Wcr = 46.97, 38.76
#这些数据是通过85万多个肤色点统计得出的;
def Before_Model(Y,CbY,CrY):
    if Y<Kl or Y>Kh:
        if Y<Kl:
            CbY_evr=108+(Kl-Y)*(118-108)/(Kl-Ymin)
            CrY_evr=154-(Kl-Y)*(154-144)/(Kl-Ymin)
            WcbY=Wlcb+(Y-Ymin)*(Wcb-Wlcb)/(Kl-Ymin)
            WcrY=Wlcr+(Y-Ymin)*(Wcr-Wlcr)/(Kl-Ymin)
        elif Y>Kh:
            CbY_evr=108+(Y-Kh)*(118-108)/(Ymax-Kh)
            CrY_evr=154+(Y-Kh)*(154-132)/(Ymax-Kh)
            WcbY=Whcb+(Ymax-Y)*(Wcb-Whcb)/(Ymax-Kh)
            WcrY=Whcr+(Ymax-Y)*(Wcr-Whcr)/(Ymax-Kh)
        CbKh_aver = 108 + (Kh - Kh) * (118 - 108) / (Ymax - Kh)
        CrKh_aver = 154 + (Kh - Kh) * (154 - 132) / (Ymax - Kh)
        Cb=(CbY-CbY_evr)*Wcb/WcbY+CbKh_aver
        Cr=(CrY-CrY_evr)*Wcr/WcrY+CrKh_aver
    else:
        Cb=CbY
        Cr=CrY
    return Cb,Cr

上面这段代码是
为了使肤色聚类不受亮度Y 的影响将YCbCr 颜色空间中的色度Cb、Cr进行非线性变换。
参考文献:
《基于分裂式K均值聚类的肤色检测方法》

1
2
3
4
def ellipse2():
    skincbcr=np.zeros([256,256],cv2.CV_8UC1)
    ellipse=cv2.ellipse(skincbcr,[113,155.6],[23.4,15.6],43.0,0.0,360.0,[255,255,255],-1)
    return ellipse

目前学习到的就这些了。

展开全文 >>

pyinstaller遇到的问题

今天花了一早晨的时间，终于将最近做的一个安全帽检测器收了尾，做的虽然不是很满意，但是
期间踩了很多坑，涨了一波姿势，同时也是自己第一次真正的用机器学习做项目，感受了这个漫长
的过程，学了不少东西，不亏不亏。
因为当初觉得python写机器学习的东西，方便一点，所以就用python开始了。快收尾的时候，发现，
其实用c++写问题也不大啊，这个项目主要用的wx和opencv依然可以随便用啊。。但是都快写完了，
这么多代码移植一下真的太费时间了，但是python写完之后，就是要面对这么一个问题，我们写成
的程序是多个py文件，如果换一台没有解析器，第三方模块没装全的机器，那就废了，所以就要想
办法把他们打包成一个exe文件，开始吧。
首先两种选择，py2exe和pyinstaller，对比了一下，py2exe也太老了，而且在官网好像找不到
支持py3的了，而且网上的一致观点是py2exe效果很不好，打包成的exe存在包括兼容性的一系列
问题，于是便选择了pyinstaller。
直接

1
pip install pyinstaller

装好了，那么开始做吧，于是第一步就遇到了问题，cmd下进入要打包的py文件所在的目录，
执行pyinstaller -F 指令，然后结果只有一句：
failed to create process
这是什么鬼，于是到网上一搜，这样解决：
去第三方库中找到文件:

1
pyinstaller-script.py

会看到第一行是一个python.exe路径,把这行的”路径”用引号括起来，保存，再次执行命令，
这个信息不会在出现了，开始打包过程，然而在打包过程就失败了，出现错误。
这里就要说一个问题了，目前来说，以pip方式安装的pyinstaller对于py3.5以上的版本
支持并不是很完善，因此就会出现各种各样的错误。
所以到这里，依然决然的从github上下载源zip的方式进行安装。
下载以后，然后解压后，里面会有一个pyinstaller文件夹，进入
执行setup.py文件即可：

1
python setup.py install

然后再次执行打包命令就不会出现错误了…有点坑…
这次成功打包了，生成了exe文件，于是迫不及待的运行它。
然后跳出了控制台窗口，然后就报错了，说找不到我导入的模块，这个找不到的文件夹就是我的
的打包目录，因为我打算打包的是几个py文件，并不是只有一个，我以为只要将程序入口设在一个
py文件当中，然后将目录里其他模块import就好了，但是报错找不到。
上网了寻求了一下帮助，并没有找到打包多个py的正确方式，所有的打包都是针对单个py的，对于
其他导入的第三方模块，比如cv2，wx之类的并未报错，说明打包过程至少是可以识别python的
sit-packages文件夹的，所以可能将我们的文件都放进去打包应该没啥问题了吧。想了一下，
那我的第三方文件夹就乱套了啊。果断放弃，把所有的文件的代码全部整合到一个py文件当中，
重新打包，ok了，这个时候我的UI窗口可以执行了。
开始运行，导入视频，然后又出问题了，我在控制台中捕捉的视频帧数是0，说明程序中这行代码
并未执行：

1
2
cv2.VideoCapture(filename)
#通过文件读取视频

然而这句代码却执行了

1
2
cv2.VideoCapture(0)
#通过摄像头读取视频

解决这个问题的方法：
在cv2库中找到如下dll文件：

1
2
3
opencv_ffmpegXYZ.dll
opencv_ffmpegXYZ_64.dll
#XYZ为opencv版本号;

32位选择第一个，64位选择第二个，然后将这个文件夹放到exe文件所在的同级目录当中，
在我们编写的py文件开头加入一下代码：

1
2
3
4
5
6
import ctypes
try:
    temp=ctypes.windll.LoadLibrary('opencv_ffmpegXYZ.dll')
except:
    print("load binary failed")
    pass

然后重新打包，一定要将动态链接库和打包的exe一直放在一起。
同理，我们训练完成的模型bin文件，音频文件以及其他exe中要使用到的文件，
都要放在和exe同级文件夹当中，并且所有类似load()的函数都直接写成文件名即可，比如：

1
2
3
4
5
6
ctypes.windll.LoadLibrary('opencv_ffmpeg341.dll')
hog = cv2.HOGDescriptor()
hog.load("myHogDector1.bin")
winsound.PlaySound("901028.wav",winsound.SND_FILENAME)

打包完成，能够正常运行了,大功告成。

展开全文 >>

<学习笔记> gxx_personality_sj0错误

今天在codeblocks中编译写好的坦克大战的代码，结果报了这个错误：

1
2
3
undefined reference to 'gxx_personality_sj0'
undefined reference to '_Unwind_SjLj_Resume'
undefined reference to '_Unwind_SjLj_Register'

看不太懂，于是上网搜了一下，这是编译器的问题，可能是编译和链接不一致导致的。
于是自己好好检查了一遍文件结构，在src文件当中有一个CommonClass.a压缩文件，
里面是一个CommonClass.o对象文件，相当于一个编译的中间文件，将这样的多个对象文件
连接在一个就可以构成我们的可执行文件。
也就是说，接下来我们rebuild要做的工作是先将我们自己写好的cpp文件编译成一系列的obj文件，
然后再整合成为exe文件，但是当初编译它里面自带的obj文件的编译器和我们现在链接的编译器
可能不同，于是我们换了一个版本的gcc，再次进行编译，这次成功了。
所以导致这个错误的原因就在于编译器的问题，可能是链接和编译版本不一致导致的。

展开全文 >>

<学习笔记> HOG+SVM人体识别基本框架

先来说一下HOG吧，全称方向梯度直方图，是一种用来图像物体检测的描述子。
结合支持向量机，当前广泛应用于行人检测。
生成过程：

(1)图像归一化，为了减小光照等各种因素对特征目标提取的影响，对图片进行预处理。
(2)计算像素梯度，个人理解这里算是为什么这个算法对行人检测有着良好的效果，像素梯度我认为可以看作是
描述某个像素点与其周围像素点差异性的大小度量。
就像我们地理学里所知的登高设色地形图一样，我们的图像中的点的像素值就可以看作地形图的海拔高度值。

如图所示，

肉眼大致估计，向量2所指的方向就是梯度方向，这也是为什么梯度这个量对于图像边缘识别那么有用了。
像素梯度的计算公式就不多说了(天知地知，你我都知)。
(3)这一步算是开始落实我们的HOG算法的核心了，我们采取“投票制度”。
这里有一个三层级别
Pixel级(像素)，Cell级(单元)，Block级(块)。
Block是划分图像的第一次，通常使用的HOG结构有三种，矩形，圆形和环绕。

在进行下一步的细分，每一个块由若干个单元组成，每个单元又由若干个像素组成。

如下图所示，我们的Cell就相当于一个投票空间，每个像素都有投票的权利，每个像素都拿这个权利选择一个
方向，这个方向就是像素点自己的梯度方向。

在每个Cell中，统计票数，也就是梯度统计，以梯度方向为横坐标，梯度大小和为纵坐标绘制梯度方向直方图。
方向可以是180度(无符号),也可以是360度(有符号)，然后将这个角度划分为若干个区间，

如图中所示，采用360度，每20度为一个区间，区间边界度数对应一个直方柱。
就如上图中，红色向量属于20-40度区间，那么假设它是35度，它的梯度大小为40，那么如果采取线性插值法
进行投票，若不对大小进行处理(有时候会进行权重处理，比如梯度幅值的平方等等)，那么就应该向40度投票
30，向20度投票10。
在每个cell中进行上述统计绘制完整的直方图。
这里选取合适的空间大小和区间大小是最大化检测能力的关键，网上关于这方面的说法比较多，不多说了。
(4)梯度归一化
cell投票完成后，然后将cell整合成block，进行特征向量的归一化，这一步的目的是使向量空间对
光照，阴影，边缘变化具有鲁棒性。
归一化的方式有一下几种：

在人体检测当中，方法a的表现是最好的。
还有就是cell不是某个block私有的，他可能同时属于多个block，也就是block之间可能是重叠的。
(5)整合特征向量
最终得到的向量是一个Block_numCell_numOrientation_num维的向量。
以上就是对HOG特征提取的一般过程，更多理论知识网上都有，自行学习吧。

接下来，就来说一下利用Opencv自带的HOG和SVM进行行人检测的基本框架，由于本次项目有特定需求，
所以不使用它自带的已训练好的模型，决定自己进行训练。
训练基本模块组织：
(1)加载并读取正负样本
正样本：含有人体的图片若干；
负样本：不含人体的图片若干；
手动裁剪样本：hog检测窗口的大小为128*64，所以我们需要将样本进行手动裁剪，比如要进行行人检测，
就需要将样本中的行人裁剪出来，负样本不必做处理，但是要比上面这个尺寸大才行。
处理完成后分别存放在两个对应的文件夹当中。
(2)分别提取正样本和负样本的HOG特征，这里相当于原料精化的过程。
同时对正样本和负样本做标签，比如将正负样本分别标记为1和0。
提取的特征向量和标签就是分类器的输入了。
计算特征：hog.compute()
(3)分类器设置与训练
设置SVM的初始参数，包括各种核函数和损失函数等等。
一定要选择线性SVM，因为cv.HOGDescriptor()检测函数只支持线性检测。
保存参数：svm.save();
加载参数：cv.ml.SVM_load();
训练：svm.train();
(4)训练之后，保存训练结构为一个文本文件。
这个文本文件中存储着几个量：
support vector：支持向量；
alpha:约束条件参数(拉格朗日乘子，支持向量机的训练可以看作找最优alpha的过程)；
rho:判决函数偏置量；
(5)自举法(bootstrap)处理难样本(hardexample)。
利用上面习得的分类器对我们的负样本进行检测，检测方法:多尺度检测
hog.detectMultiScale()
我们图片里的人有大有小，所以要给图像一个多层次的大小度量，也就是我们的图像金字塔，
在每一层中进行检测，然后做标注，将结果综合。

将负样本中错误识别出来的行人区域继续进行上述(2)过程，将得到的特征和我们第一次
得到的特征综合，再次重复(3)的训练过程。
(6)存储模型，大功告成。

展开全文 >>

<学习笔记> 线性回归

之前我们构造了各类分类器，目的是分类，目标变量是标称型数据，从这一篇开始，开始
对连续型数据进行预测，那么我们通过什么方法进行呢，那就是回归。
我们利用回归寻找两种或两种以上变量间的相互依赖关系。
回归在化学分析，遗传学等我们常见的学科当中有着非常广泛的应用。
这篇文章主要对线性回归进行学习总结。

一.利用线性回归找到最佳拟合直线
回归的目标是预测数值型的目标值。
最直接的方式是通过输入写出一个目标的表达式：

1
y=ax+b;

这就是所谓的回归方程，其中的a,b可以被称作回归系数。
我们求得回归系数的过程就是回归。
一旦我们知道了a,b两个值的大小，如果给定输入x，那么结果y就很容易得到了。
那么怎样从一大堆的数据中求出回归方程呢，假设输入矩阵放入矩阵x中，而回归系数存放在向量
w中，那么对于给定的数据X1，预测结果将会通过

1
Y1=(X1^T)*w

给出，现在的问题是，我们给出的是一系列的x和y，要得到的是w，一个常用的方法是得到误差最小的
w。这里的误差是指预测y和实际y之间的差距，直接使用该误差的简单累计会使得结果相互抵消，所以
我们将其平方化处理，采用平方误差。
那么平方误差即可写作：

如果是用矩阵表示的话，那么可以写作：

都w求导，得到X^T(Y-Xw),令其等于0，接触w如下：

w头上有个小标记，代表这是当前可以估计出的w的最优解。从现有数据估计出的w可能与真实的w
值有一定的差异，这个解只是代表一种当前的最优估计，因此要加以区分。
还要注意上面的计算涉及到矩阵求逆，那么这个方程只有在逆矩阵存在的条件下才适用，
所以我们编写代码的时候要进行判断。
下面是用python实现的标准回归函数以及数据读取函数：

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
from numpy import *
def loadDataSet(filename):
    numfeat=len(open(filename)
                .readline().split('\t'))-1
    dataMat=[]
    labelMat=[]
    fr=open(filename)
    for line in fr.readlines():
        lineArr=[]
        curLine=line.strip().split('\t')
        for i in range(numfeat):
            lineArr.append(float(curLine[i]))
        dataMat.append(lineArr)
        labelMat.append(curLine[-1])
    return dataMat,labelMat
def standRegres(xArr,yArr):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    xTx=xMat.T*xMat
    if linalg.det(xTx)==0.0:
        print("this matrix is singular,cannot do inverse")
        return
    ws=xTx.I*(xMat.T*yMat)
    return ws

上面的数据载入不用过多介绍，第二个函数只用几行代码就实现了我们上面要求的结果，
先将数据向量矩阵化，然后计算平方，通过调用linalg模块求行列式，如果行列式为0，
那么矩阵不可逆，反之，计算结果。
几乎任意的数据点都可以用上述的方法建立模型，但是如何评价我们模型的好坏呢。
对于下面两个图，我们在这两个数据集上分别做线性回归，将得到完全一样的拟合直线。
显然两个数据是不一样的，那么模型分别在二者的效果比较起来如何。
有一种方法就是比较预测值和真实值之间的匹配程度，计算两个序列的相关系数。

上面的两个图具有相同回归系数(0,2.0)，但是上图的相关系数是0.58，
下图的相关系数达到了0.99。
numpy库中提供了计算相关系数的方法：
通过corrcoef(yEstimate,yActual)来计算预测值和真实值之间的相关性。

上面我们将数据当作直线来建模，寻找了最佳拟合曲线，我们还可以发现数据当中存在的潜在模式。

二.局部加权线性回归
线性回归会遇到的问题是有可能出现欠拟合现象，因为它要求的是具有最小均方误差的无偏估计。
显然，欠拟合会使得我们得不到好的结果。所以我们要想办法引入一些偏差，来降低我们的均方误差。
其中的一种方法就是局部加权线性回归(LWLR)。
在这个算法当中，我们基于预测点附近的的每个点一定的权重，然后在这个子集上，进行普通的回归。
该算法解出的回归系数w的形式如下：

其中w是一个权重，用来给每个数据点赋予一个权重，
LWLR使用类似SVM核的“核”来对附近的点基于更高的权重，核的类型可以自由选择，
高斯核对应的权重如下：

那么很容易看出来，点x与x(i)越近，w(i,i)将会越大。
我们可以通过调节上式中的k值，来决定我们对附近的点给予多大的权重。
下图展示了参数k与权重的关系。

从图中可以观察得出，当k值越大的时候，越多的数据用于训练回归模型。
下面给出局部加权线性回归函数的代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
def lwlr(testPoint,xArr,yArr,k=1.0):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    m=shape(xMat)[0]
    weights=mat(eye((m)))
    #创建对角权重矩阵，随着样本点与待估测的点的距离的递增，权重将指数级递减。
    for j in range(m):
        diffMat=testPoint-xMat[j,:]
        weights[j,j]=exp(diffMat*diffMat.T/(-2.0*k**2))
    xTx=xMat.T*(weights*xMat)
    if linalg.det(xTx)==0.0:
        print("cannot inverse")
        return
    ws=xTx.I*(xMat.T*(weights*yMat))
    return testPoint*ws
def lwlrTest(testArr,xArr,yArr,k=1.0):
    m=shape(testArr)[0]
    yHat=zeros(m)
    for i in range(m):
        yHat[i]=lwlr(testArr[i],xArr,yArr,k)
    return yHat

得到的拟合直线如下图：

第一幅图取k=1.0,所有数据等权重，最佳拟合直线与标准回归一致。
第二幅图取k=0.01，抓住了数据的潜在模式，得到了很好的结果。
第三幅图取k=0.003,纳入了很多的噪声点，拟合的直线与数据点过于靠近了，
这就是个过拟合的例子。而第一幅图是一个欠拟合的例子。
所有我们在选择k的时候，并不是越大或者越小越好，而是要通过对比不同的模型，来综合
考虑我们选取的参数值，找到最佳模型。

三.缩减系数来理解数据
有时候我们会遇到数据特征比样本点还多的情况，如果使用之前的方法来进行计算，会出错，
特征比样本点还多，说明输入数据的矩阵不是满秩矩阵，非满秩矩阵求逆会出现很多问题。
为了解决这个问题，引入了“岭回归”，这是其中一种缩减方法，
除此之外，还有lasso法，这种方法效果最好但是最复杂。
最后还有一种“前向逐步回归”,既能得到满意的效果，也容易实现。

1.岭回归
岭回归就是在矩阵X.TX加上一个λE使得矩阵非奇异。当X非满秩，或者是矩阵的列与列之间的
相关性太高的话，也就是X.TX接近奇异的时候，计算(X.T*X)的逆的时候误差会很大。
在这种情况下，回归系数的计算公式为:

岭回归最初用来处理特征数多于样本数的情况(直观上就是矩阵不满秩)，现在也用于在估计
中加入偏差，来得到更好的估计，通过引入λ来限制所有w的和，通过引入该惩罚项，
减少不重要的参数，这个技术就是缩减。

训练方法类似，通过预测误差最小化得到λ：
得到数据集，一部分测试，剩余的用来训练参数w。
训练后用用测试集预测回归性能。
通过选取不同的λ重复上述过程，最终选取使得误差最小的值。
代码实现:

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
def ridgeRegres(xMat,yMat,lam=0.2):
    xTx=xMat.T*xMat
    denom=xTx+eye(shape(xMat)[1])*lam
    if linalg.det(denom)==0.0:
        print("this matrix is singular,cannot do inverse")
        return
    ws=denom.I*(xMat.T*yMat)
    return ws
def redgeTest(xArr,yArr):
    xMat=mat(xArr)
    yMat=mat(yArr).T
    yMean=mean(yMat,0)
    yMat=yMat-yMean
    xMeans=mean(xMat,0)
    xVar=var(xMat,0)
    xMat=(xMat-xMeans)/xVar
    #标准化,减均值，初方差；
    numTestPts=30
    三十个不同lambda值进行测试；
    wMat=zeros((numTestPts,shape(xMat)[1]))
    for i in range(numTestPts):
        ws=ridgeRegres(xMat,yMat,exp(i-10))
        wMat[i,:]=ws.T
    return wMat

2.前向逐步回归
属于一种贪心算法，因为每一步都尽可能的减小误差。
一开始，所有的权重都设置为1，然后每一步所做的决策是对某个权重增加或减少一个很小的值。
代码实现:

1
2
3
def rssError(yArr,yHatArr):
    return ((yArr-yHatArr)**2).sum()
#预测误差大小；

1
2
3
4
5
6
def regularize(xMat):
    xVar=var(xMat,0)
    xMeans=mean(xMat,0)
    xMat=(xMat-xMeans)/xVar
    return xMat
#标准化；
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
def stageWise(xArr,yArr,eps=0.01,numIt=100):#迭代步长，迭代次数；
    xMat=mat(xArr)
    yMat=mat(yArr).T
    yMean=mean(yMat,0)
    yMat=yMat-yMean
    xMat=regularize(xMat)
    m,n=shape(xMat)
    returnMat=zeros((numIt,n))
    ws=zeros((n,1))
    wsTest=ws.copy()
    wxMax=ws.copy()
    for i in range(numIt):
        print(ws.T)
        lowestError=inf
        for j in range(n):
            for sign in [-1,1]:
                wsTest=ws.copy()
                wsTest[j]+=eps*sign
                yTest=xMat*wsTest
                rssE=rssError(yMat.A,yTest.A)
                if rssE<lowestError:
                    lowestError=rssE
                    wsMax=wsTest
        ws=wsMax.copy()
        returnMat[i,:]=ws.T
    return returnMat

当我们应用缩减方法时，模型也就增加了偏差，与此同时减少了模型的方差。
在这之后，我们就要权衡偏差与方差，是我们的模型得到更进一步的优化。

展开全文 >>

<leetcode> Complex Number Multiplication

因为这道题让我想起了编译里的词法分析，所以记下来吧。
题目描述：
Given two strings representing two complex numbers.
You need to return a string representing their multiplication. Note i2 = -1 according to the definition.
Example 1:

1
2
3
Input: "1+1i", "1+1i"
Output: "0+2i"
Explanation: (1 + i) * (1 + i) = 1 + i2 + 2 * i = 2i, and you need convert it to the form of 0+2i.

Example 2:

1
2
3
Input: "1+-1i", "1+-1i"
Output: "0+-2i"
Explanation: (1 - i) * (1 - i) = 1 + i2 - 2 * i = -2i, and you need convert it to the form of 0+-2i.

Note:
The input strings will not have extra blank.
The input strings will be given in the form of a+bi, where the integer a and b will both belong to the range of [-100, 100]. And the output should be also in this form.
代码：

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
#include<sstream>
#include<string>
class Solution {
public:
    string complexNumberMultiply(string a, string b) {
        int a1,ai,b1,bi;
        int *pa=getnum(a);
        int *pb=getnum(b);
        a1=pa[0];ai=pa[1];
        b1=pb[0];bi=pb[1];
        int ans1=a1*b1-ai*bi;
        int ansi=a1*bi+ai*b1;
        stringstream ss;
        ss<<ans1;
        string str1;
        ss>>str1;
        ss.clear();
        ss<<ansi;
        string stri;
        ss>>stri;
        return str1+'+'+stri+'i';
    }
    int* getnum(string s){
        int s1,si;
        const char *p=s.data();
        int sf=0;
        int temp=0;
        bool flag1=false,flagi=false;
        while(*p!='\0'){
            if(*p=='-'&&sf==0){
                flag1=true;
            }
            if(*p=='-'&&sf!=0){
                flagi=true;
            }
            if(*p>='0'&&*p<='9'){
                temp=temp*10+(*p-48);
            }
            if(*p=='+'){
                s1=temp;
                temp=0;
            }
            if(*p=='i'){
                si=temp;
            }
            sf++;
            p++;
        }
        if(flag1){
            s1=-s1;
        }
        if(flagi){
            si=-si;
        }
        int *ans=new int[2];
        ans[0]=s1;
        ans[1]=si;
        return ans;
    }
};

展开全文 >>

<学习笔记> AdaBoost元算法

元算法就是对其他算法进行组合的一种方式。
bagging方法：自举汇聚法，从原始数据中选择s次之后得到s个新数据的技术。新数据集和原数据集的
大小相等。每个新数据集中的数据都是通过在原数据集随机选择来进行替换得到的。
boosting方法：类似bagging，boosting不同分类器权重不同，而bagging相同。
adaboosting运行过程：
训练数据中的每个样本，并赋予一个权重，这些权重构成了向量D。一开始，这些权重都初始化为相等的值。
首先在训练中训练出一个弱分类器并计算该分类器的错误率，然后再同一个数据集上再次训练弱分类器。
在分类器的第二次训练当中，将会调整每个样本的权重，第一次分对的样本的权重降低，分错的权重升高。
为了从所有的弱分类器中获得最终的分类结果，adaboosting为每个分类器都分配了一个权重alpha，
这些alpha值是基于每个弱分类器的错误率进行计算的。
错误率e就是分类错样本数目的和总样本数目的比值。
alpha的计算公式：

1
α=1/2ln(1/e-1)

adaboosting算法的流程可以用下面的图表示：

计算出alpha值之后，就要对我们的向量D进行更新，D的计算方式如下：

计算出D之后，算法开始下一轮迭代，就这样不断重复训练调整权重，边界就是训练错误率为0或者
弱分类器的个数达到指定的值为止。

基于单层决策树构建弱分类器
之前已经接触过了决策树，那么单层决策树就是仅基于单个特征，只有一个分叉过程，也叫做决策树桩。

1
2
3
4
5
6
7
8
9
def loadimpdata():
    datamat=matrix([
        [1.,2.1],
        [2.,1.1],
        [1.3,1.],
        [1.,1.],
        [2.,1.]])
    classlabels=[1.0,1.0,-1.0,-1.0,1.0]
    return datamat,classlabels

上面是我们构造的数据集，将datamat中的数对放入我们的坐标轴，如果想要找到一条平行与坐标轴的直线来将所有的
1型和-1型分开，显然不可能。AdaBoost需要将多个单层决策树结合才能对该数据集进行分类。
上面是数据载入的函数，接下来，要构建多个函数lai建立单层决策树。
第一个函数用于测试是否存在某个值小于或者大于我们正在测试的阈值。
第二个函数用来找到错误率最低的单层决策树。

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
def stumpclassify(datamatrix,dimen,threshval,threshIneq):
    retarray=ones((shape(datamatrix)[0],1))
    if threshIneq =='lt':
        retarray[datamatrix[:,dimen]<=threshval]=-1.0
    else:
        retarray[datamatrix[:,dimen]>threshval]=1.0
    return retarray
def buildstump(dataarr,classlabels,D):
    datamatrix=mat(dataarr)
    labelmat=mat(classlabels).T
    m,n=shape(datamatrix)
    numsteps=10.0
    beststump={}
    bestclassest=mat(zeros(m,1))
    minerr=inf
    for i in range(n):
        rangemin =datamatrix[:,i].min()
        rangemax =datamatrix[:,i].max()
        stepsize=(rangemax-rangemin)/numsteps
        for j in range(-1,int(numsteps)+1):
            for inequal in ['lt','gt']:
                threshval=(rangemin+float(j)*stepsize)
                predictedvals=stumpclassify(datamatrix,i,threshval,inequal)
                errarr=mat(ones(m,1))
                errarr[predictedvals==labelmat]=0
                weightederror=D.T*errarr#计算加权错误率
                if weightederror<minerr:
                    minerr=weightederror
                    bestclassest=predictedvals.copy()
                    beststump['dim']=i
                    beststump['thresh']=threshval
                    beststump['ineq']=inequal
    return beststump,minerr,bestclassest

第一个函数的结构很简单，不做过多解释。
第二个函数最重要的部分就在三层循环的部分，
第一层循环在数据集的所有特征上遍历，通过计算最大值和最小值我们来得知需要多少步长。
第二层循环是在上面选取的值上进行遍历。
可以将阈值设置为整个取值范围之外可以，所以需要额外的步骤。
第三层循环是在大于和小于之间切换不等式。
在最内层循环之中，在数据集以及三个循环变量上面使用第一个函数，基于循环变量，返回分类结果，
接下来构造一个列向量errarr，如果预测值和真实的label值不同，那么填1，然后将错误向量和权重D
相乘，那么就得到了加权错误率，这里就是adaboost和分类器交互的地方，D就是用来评价分类器的。

完整adaboost实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
def adaboosttrainDS(dataarr,classlabels,numIt=40):
    weakclassarr=[]
    m=shape(dataarr)[0]
    D=mat(ones((m,1))/m)
    aggclassest=mat(zeros((m,1)))
    for i in range(numIt):
        beststump,error,classest=buildstump(dataarr,classlabels,D)
        print("D:",D.T)
        alpha=float(0.5*log((1.0*error)/max(error,1e-16)))
        beststump['alpha']=alpha
        weakclassarr.append(beststump)
        print("classest:",classest.T)
        expon=multiply(-1*alpha*mat(classlabels).T,classest)
        D=multiply(D,exp(expon))
        D=D/D.sum
        #为进行下一次迭代计算D
        aggclassest+=alpha*classest
        print("aggclassest:",aggclassest)
        aggerrors=multiply(sign(aggclassest)!=mat(classlabels).T,ones((m,1)))#sign()函数：正负分离；
        errorrate=aggerrors.sum()/m
        print("total error:",errorrate,"\n")
        if errorrate==0.0:
            break
    return weakclassarr

上面的函数的参数包括数据集，类别标签以及迭代次数，迭代次数由用户指定。
如果我们在预设的迭代次数未到达之前就达到errorrate为0，那么就直接退出迭代，不需要指定的那么多了。
上面我们融合的是单层决策树，也就是将其作为基分类器。
我们可以通过修改代码将其修改为其他的弱分类器，前面的任何一种都可以。
向量D非常重要，包含了每个数据点的权重，一开始，赋予了每个权重值相同的值，在后续的迭代中，adaboost
会在增加错误分类的数据的权重的时候，降低分类正确的数据的权重。
D是一个概率分布的向量，因此所有的元素之和为1.0。为了满足这个需求，那么一开始就将他切为1/m。
同时，程序还会建立一个列向量aggclassest，来记录每个数据点的类别估计累积值。
adaboost算法的核心在于for循环，该循环运行指定次数或者直到错误率为0退出。
循环中的第一件事是建立单层决策树。随后是计算alpha的值，这个值告诉总分类器这次的单层决策树的
结果的权重，其中的max函数用来避免除零异常。
然后，将这个值加入到beststump字典中，这个字典又加入到列表中。
多个分类器得出结果后，就要进行“汇总”了：
下面的函数就可以利用基于adaboosttrainDS()中的弱分类器对数据进行分类了。

1
2
3
4
5
6
7
8
9
10
11
12
def adaclassify(dattoclass,classifierarr):
    datamatrix=mat(dattoclass)
    m=shape(datamatrix)[0]
    aggclassest=mat(zeros((m,1)))
    for i in range(len(classifierarr)):
        classest=stumpclassify(datamatrix,
                               classifierarr[i]['dim'],
                               classifierarr[i]['thresh'],
                               classifierarr[i]['ineq'])
        aggclassest+=classifierarr[i]['alpha']*classest
        print(aggclassest)
    return sign(aggclassest)

文件载入与自检特征：

1
2
3
4
5
6
7
8
9
10
11
12
def loadDataSet(fileName):
    numFeat = len(open(fileName).readline().split('\t'))
    dataMat = []; labelMat = []
    fr = open(fileName)
    for line in fr.readlines():
        lineArr = []
        curLine = line.strip().split('\t')
        for i in range(numFeat - 1):
            lineArr.append(float(curLine[i]))
        dataMat.append(lineArr)
        labelMat.append(float(curLine[-1]))
    return dataMat, labelMat

其他分类度量指标
之前的分类器我们基本都是基于错误率来衡量成功程度的。
错误率也就是分类错误的占总测试样例中的比例。错误率越低，我们的分类器越成功。
实际上，这样的度量遮盖了样例如何被分错的事实。
接下来，我们使用一个名叫混淆矩阵的工具来更好的理解分类中的错误。
下面这个图是关于在房子周围可能发现动物类型的预测。

很显然，一个完美的分类器应当是使得这个矩阵的非对角元素全为0。
接下来，考虑一个二分类的混淆矩阵：

通过这个矩阵，我们可以定义一些概念：
将一个正例判为正例，即会产生一个真正例，也称真阳。
将一个反例正确地判定为一个反例，则认为产生了一个真反例，也称真阴。
在分类中，某个类别的重要性高于其他类别的时候，就可以定义出多个比错误率更优秀的
的新指标。
第一个指标称为正确率，表示真阳占所有预测为正例的比例。
第二个指标称为召回率，表示真阳站所有真实结果为正例的比例。
同时想要达到高正确率和高召回率是困难的。

另一个用于度量分类中非均衡性的工具是ROC曲线，ROC代表接收者操作特征，
ROC曲线是根据一系列不同的二分类方式，以真阳率为纵坐标，假阳率为横坐标绘制的。
一般的评价方法是将我们的结果分为两类，然后进行评估，ROC曲线的评价方式与一般的
最大的区别在于可以根据实际情况，我们可以把实验成果划分为多个有序的分类序列，再
进行分析。
通过ROC曲线图像，我们可以比较不同分类器的优劣，将不同算法得到的曲线绘制到同一个坐标
系当中，很清楚就能判断出孰优孰劣，越是靠近左上角部分的曲线那么肯定更优秀咯。
同样，这个曲线是根据阈值变化来绘制曲线的，因此我们同时可以根据不同的阈值来使用表现最
好的分类器，将他们组合起来可能会得到更好的结果。
我们可以根据曲线下面积(AUC)来比较不同的ROC曲线，一个完美的分类器AUC值为1，而随机猜测
的则为0.5，似乎符合我们对概率的一般印象。

上图显示了蛋白质水解的酶体消化模式的提取，上图中有四个重要的点：
(0.0,1.0):超级完美的分类器。
(1.0,0.0):成功避开了所有的正确答案的分类器。。自不必多言了。。
(1.0,1.0):分类器预测所有样本都为正样本。
(0.0,0.0):分类器预测所有样本都为负样本。
还有一条虚线，就是随机猜测咯，大概有一半猜对，一半猜错。
其实相较混淆矩阵，我们对ROC可以有一种感受，就是其实相当于将离散的布尔值连续化了，
就像平时我们做出决断并不是之后好坏一样，可能在不同决策水平上，会得到不大一样的结果，
这种结果会依赖于决策水平。
那么如何绘制一条ROC曲线呢。
最简单直白的思路，取点描线。
因此必须知道每个样本被判断为阳性或者阴性的可信赖程度。
分类器有一个重要的功能就是概率输出，也就是描述一个样本有多大的概率属于正样本，
有多大的概率属于负样本，例如朴素贝叶斯告诉了我们一个可能性，Logistic回归给出一个
数值输入到Sigmoid()函数当中，SVM中也会计算出一个这样的值输入到sign()函数当中。
这些值就可以用于衡量这个分类器的得出的结果值不值得信赖。
为了创建ROC曲线，我们会先做一个排序，从排名最低的样例开始，排名更低的判定为反例，
排名更高的被判定为正例。然后将这个最低的移到次低的样例中去，如果属于正例，修改
真阳率，输入反例，对假阴率进行修改。
其实相当于我们设计了一个cursor，从排名最低的开始，在cursor指在某个被排列项的时候，
都将其作为一个阈值，然后比它排名高的自然可看作正样本，反之。
下面举了这样一个例子，

其中class为p即为正样本，为n即为负样本，score是每个测试样本属于正样本的概率，也就是
我们相信其为正样本的可信赖度。
那么以第5个样本为例，其概率为0.55，那么很显然如果我们将其作为阈值，那么样例1，2，3，
4都可以看作正样本，其他的都是负样本，cursor到达一个样例的时候，都会得到一组FPR和
TPR，那么上面这个例子就可以得到20组值对。
绘制到图上：

接下来，我们将绘制ROC曲线和计算AUC的代码加入：

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
def plotROC(predstrengths,classlabels):
    import matplotlib.pyplot as plt
    cur={1.0,1.0}
    ySum=0.0
    numPosClas=sum(array(classlabels)==1.0)
    yStep=1/float(numPosClas)
    xStep=1/float(len(classlabels)-numPosClas)
    sortedIndicies=predstrengths.argsort()
    fig=plt.figure()
    fig.clf()
    ax=plt.subplot(111)
    for index in sortedIndicies.tolist()[0]:
        if classlabels[index]==1.0:
            delX=0;delY=yStep
        else:
            delX=xStep;delY=0
            ySum+=cur[1]
        ax.plot([cur[0], cur[0] - delX], [cur[1], cur[1] - delY], c='b')
        cur = (cur[0] - delX, cur[1] - delY)
    ax.plot([0, 1], [0, 1], 'b--')
    plt.xlabel('False Positive Rate');
    plt.ylabel('True positive Rate')
    plt.title('ROC curve for AdaBoost Horse colic Detection System')
    ax.axis([0, 1, 0, 1])
    plt.show()
    print("the Area Under the Curve is:", ySum * xStep)

基于代价函数的分类器决策控制
我们上面介绍了通过调节分类器的阈值来处理非均匀代价问题的方法，
除此之外，我们还有一些其他的方法，其中的一种被称作代价敏感学习。
有如下两个代价矩阵：

第一张表是目前为止分类器的代价矩阵，那么总代价可以计算：

1
TP*0+FN*1+FP*1+TN*0

同理第二张表可以计算：

1
TP*(-5)+FN*1+FP*50+TN*0

采用第二张表，分类错误的代价也是不同的。那么在我们构建分类器的时候，如果知道这些代价值，
就可以选择代价最小的分类器了。
在分类算法中，有很多方法可以引入代价信息。
在adaboost中，可以依据代价函数对权重向量D进行调整。
在朴素bayes中，可以选择具有最小期望代价而不是最大概率的类别作为结果。
在SVM中，可以在代价函数中对不同类别选择不同的参数C。
这些做法会给予较小的类更大的权重，即在训练中，小类只允许有更少的错误。

处理非均匀问题的数据抽样
另一种对非均匀问题调节分类器的方法，就是从输入原料着手，对训练数据进行改造加工。
可以通过两种方式实现：
过抽样：复制样例。
欠抽样：删除样例。
我们在分类过程中会有一些罕见的样例需要识别，
而这些比较特殊的罕见的样例反而会裹挟着更重要的东西，因此我们需要保留更多的信息。
那么我们对罕见样例进行过抽样处理，对那些相反的采取欠抽样。
这其中也存在缺点，我们需要了解哪些是可以剔除的，因为在某些剔除样例中也可能携带了
剩余样例中不包含的有价值信息。
我们可以采取一种方法，就是选择那些离决策边界较远的进行剔除。
例如在信用卡欺诈案例当中，可能有50例属于欺诈交易(正例，罕见样例)，
有5000例属于正常交易(反例)，如过为了达到平衡对反例进行欠抽样，那么就要删掉4950例，
而这些被删掉的样例中可能含有很多有价值信息。那么为了避免这种极端的产生，我们最好还是
反例欠抽样和正例过抽样混合的方法更好。
对正例进行过抽样，我们可以复制已有样例，也可以加入已有样例的相似样例，一种方式是加入
已有数据点的插值点，但是这可能会导致过拟合的问题。
感谢博主孔明的文章：http://alexkong.net/2013/06/introduction-to-auc-and-roc/

展开全文 >>