菜鸟入门_Python_机器学习（1）_线性可分的双月实验

系统：win7-CPU;

编程环境：Anaconda2-Python2.7，IDE：pycharm5;

来看我们的任务：

Generate random data in two groups controlled by parameters r, w, d, Na, Nb, and plot out. (ref. p61 on Haykin)
?Implement Perceptron Convergence Algorithm p54 on Haykin in PyCharm, and then apply to your data sets with following different parameters
?r=10, w=2, d=2, Na=100, Nb=100
?r=10, w=2, d=-2, Na=100, Nb=100

《Neural Networks and Learning Machines(Third Edition)》这本书中对所谓的‘双月’做了很多处理，简单来说就是生成这样一个东西：

看起来很简单，但是我第一次编程还是比较蒙蔽的，莫名的恐惧和厌恶……总之多敲敲代码就好了。好歹最后搞了出来：

# -*- coding:gb2312 -*-
import matplotlib.pyplot as plt
import numpy as np
import math
# 生成双月
def point_sample(N,r,w,d):
    x_sample = []
    y_sample = []
    inner_r = r - (w/2)
    outer_r = r + (w/2)

    #point the area A
    while True:
        data_x_A = np.random.uniform(-outer_r,outer_r)
        data_y_A = np.random.uniform(0,outer_r)
        r_A = math.sqrt((data_x_A**2) + (data_y_A**2))
        if r_A >= inner_r and r_A <= outer_r:
            x_sample.append(data_x_A)
            y_sample.append(data_y_A)
            if len(x_sample) == N:
                break
            else:
                continue
        else:
            continue

    #point the area B
    while True:
        data_x_B = np.random.uniform(-outer_r + r,outer_r + r)
        data_y_B = np.random.uniform(-d - outer_r,-d)
        r_B = math.sqrt(((data_x_B - r)**2) + ((data_y_B + d)**2))
        if r_B >= inner_r and r_B <= outer_r:
            x_sample.append(data_x_B)
            y_sample.append(data_y_B)
            if len(x_sample) == 2 * N:
                break
            else:
                continue
        else:
            continue
    plt.figure(1)
    plt.plot(x_sample,y_sample,'b*')

    data_xy = np.array([np.reshape(x_sample,len(x_sample)),np.reshape(y_sample,len(y_sample))]).transpose()

    return data_xy

最后，为了方便之后处理都转置成了二维的数据。
呃……似乎有点挫，给一个网上的版本，简化了很多：

def dbmoon(N=100, d=2, r=10, w=2):
    N1 = 10*N
    w2 = w/2
    done = True
    data = np.empty(0)
    while done:
        #generate Rectangular data
        tmp_x = 2*(r+w2)*(np.random.random([N1, 1])-0.5)
        tmp_y = (r+w2)*np.random.random([N1, 1])
        tmp = np.concatenate((tmp_x, tmp_y), axis=1)
        tmp_ds = np.sqrt(tmp_x*tmp_x + tmp_y*tmp_y)
        #generate double moon data ---upper
        idx = np.logical_and(tmp_ds > (r-w2), tmp_ds < (r+w2))
        idx = (idx.nonzero())[0]

        if data.shape[0] == 0:
            data = tmp.take(idx, axis=0)
        else:
            data = np.concatenate((data, tmp.take(idx, axis=0)), axis=0)
        if data.shape[0] >= N:
            done = False
    print data
    db_moon = data[0:N, :]
    print db_moon
    #generate double moon data ----down
    data_t = np.empty([N, 2])
    data_t[:, 0] = data[0:N, 0] + r
    data_t[:, 1] = -data[0:N, 1] - d
    db_moon = np.concatenate((db_moon, data_t), axis=0)
    return db_moon

可以发现充分利用numpy库可以让代码简化很多。OK，双月得到了，两个区域之间的距离我们定义的较大，让他们线性可分，然后用一个Rosenblatt感知机分开。

很多地方介绍机器学习都会首先介绍单层感知机，作为之后神经网络的基础，它可以看做类似于神经细胞的结构，多个输入，各自有可训练权值，经过求和并激活之后输出

到底为什么要加偏置Bias呢？一种思想是，假设模型的输出始终与真实的数据有差距，我们为了逐步调整这个差距，我们在输入层对每个样本添加了一维特征，相当于在输出端添加了一个始终为一的观察点。从这个角度去看，偏置是可以去掉的，我们可以通过不断调整输入的权值来优化网络。通常为了方便计算将偏置作为输入的一部分：

其中所谓的学习率是自己设定的，我理解的就是给误差一个权值，它不能太大，否则会让权值在更新时跳过最优值而不断震荡，也不能太小，否则会导致寻找最有权值的过程太慢。激活函数的目的是将输入的加权和归一化并更好的体现分类特性，因为我们要分类的双月是线性可分的而且是两类，激活函数用符号函数就好。而统计误差的时候我们用常用的最小均方误差（这个不了解的可以自己查），之后我们用这个误差不断更新权值，达到想要的分类效果。有两种方式来判断网络时候达到我们想要的效果，一种是限定迭代次数，另一种是看代价函数（误差）值时候降低到我们设定的足够低的标准。我这里使用第一种方法，一开始尝试迭代50次，直接上代码：

N = 100
d = -2
r = 10
w = 2
a = 0.1
num_MSE = []
num_step = []
data = dbmoon(N, d, r, w)

x0 = [1 for x in range(1,201)]
x = np.array([np.reshape(x0, len(x0)), np.reshape(data[0:2*N, 0], len(data)), np.reshape(data[0:2*N, 1], len(data))]).transpose()
m = np.array([1, 0, 0])
b_pre = [1 for y in range(1, 101)]
b_pos = [-1 for y in range(1, 101)]
b=b_pre+b_pos

def sgn(v):
    if v >= 0:
        return 1
    else:
        return -1
#compute y(n)
def compute_yn(myw, myx):
    return sgn(np.dot(myw.T, myx))
#Update the weights
def new_w(old_w, myd, myx, a):
    return old_w+a*(myd-compute_yn(old_w,myx))*myx

for ii in range(50):
    i = 0
    sum=0
    for xn in x:
        m = new_w(m, b[i], xn, a)
        sum += (b[i]-compute_yn(m, xn))**2+0.0
        i += 1
    #compute MSE
    mean_square =np.sqrt(sum/N/2)
    num_MSE.append(mean_square)
    print mean_square
    print ii
    num_step.append(ii+1)
plt.subplot(212)
plt.plot(num_step, num_MSE, 'r-')

print m
#draw The decision boundary
testx = np.array(range(-15, 25))
testy = -testx*m[1]/m[2]-m[0]/m[2]
plt.subplot(211)
plt.plot(data[0:N, 0], data[0:N, 1], 'r*', data[N:2*N, 0], data[N:2*N, 1], 'b*', testx, testy, 'g--')
plt.show()

最后的结果：

分得很好，MSE曲线下降的也很快。但是当我们把双月之间的距离调成线性不可分时，分类效果很差，而且MSE曲线不断震荡，如何应对这种情况我们之后讨论。