用Python处理多分类的线性可分问题

我们接着Logistics Regression模型的思路，考虑如何处理多分类的线性可分问题。

下面我们介绍Softmax Regression模型，这一模型是对Logistics 模型的扩展，用于处理多分类的线性可分问题

1. 理论

对于有k个分类的样本，我们规定每个样本在k个分类上的概率为

上面这个式子中，j为单个分类，y为样本类别，X为特征向量，θ为系数，我们的任务就是利用样本数据估计θ的取值。

2. 损失函数

由此我们构建损失函数：

I{}为boolean函数

3. 更新方程

再对损失函数运用梯度下降算法，得到系数矩阵θ的更新方程

有了更新方程我们就可以写程序进行迭代啦。下面是个例子

4.举个栗子

我们的训练样本为下图

                                                                                                                                                                （数据见文章末尾）

# -*- coding: utf-8 -*-

#%% 导入包
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np

#%% 导入数据
data=pd.read_table('Softmax_Regression/data.txt',sep='\t')
x=data.iloc[:,0] #特征向量一
y=data.iloc[:,1] #特征向量二
label=data.iloc[:,2] #标识

#%%画图
plt.figure('训练样本')
plt.scatter(x[label==0],y[label==0],label='第一类')
plt.scatter(x[label==1],y[label==1],label='第二类')
plt.scatter(x[label==2],y[label==2],label='第三类')
plt.scatter(x[label==3],y[label==3],label='第四类')
plt.xlabel('特征向量一')
plt.ylabel('特征向量二')
plt.legend()
plt.title('训练样本')

#%% 调整参数
feature=np.mat((x,y)).T
label=np.mat(label).T
print(feature)
print(type(feature))
print(np.shape(feature))

#%% 训练模型
w=gradientAscent(feature,label,4,50000,0.2)

#%% 输出权重
print('权重为：',w)

#%% 生成预测样本
xx=(np.random.random(10000)-0.5)*6
yy=(np.random.random(10000))*15
b=np.ones(10000)
print(np.shape(np.mat((xx,yy,b)).T))
test=np.mat((xx,yy,b)).T
re=test*w
predict=re.argmax(axis=1)

#%%
plt.figure('测试样本')
print(np.column_stack((test,predict))[:,1])
merge=pd.DataFrame(np.column_stack((test,predict)))
plt.scatter(merge[merge.iloc[:,3]==0].iloc[:,0],merge[merge.iloc[:,3]==0].iloc[:,1],5)
plt.scatter(merge[merge.iloc[:,3]==1].iloc[:,0],merge[merge.iloc[:,3]==1].iloc[:,1],5)
plt.scatter(merge[merge.iloc[:,3]==2].iloc[:,0],merge[merge.iloc[:,3]==2].iloc[:,1],5)
plt.scatter(merge[merge.iloc[:,3]==3].iloc[:,0],merge[merge.iloc[:,3]==3].iloc[:,1],5)
plt.title('测试样本')

#%% 梯度更新函数
def gradientAscent(feature,label_data,k,maxCycle,alpha):
         #input: feature_data(mat):特征
          #     label_data(mat):标签
          #     k(int):类别的个数
          #     maxCycle(int):最大迭代次数
          #     alpha(float):学习率
          #     m(int)样本个数
          #     n(int)变量特征
        #output: weights(mat):权重
     feature_data=np.column_stack((feature,np.ones(np.shape(feature)[0])))
     m=np.shape(feature_data)[0]
     n=np.shape(feature_data)[1]
     weights=np.mat(np.ones((n,k)))
     i=0
     while i<=maxCycle:
         err=np.exp(feature_data*weights)
         if i%100==0:
             print ("\t-------iter:",i,\
             ",cost:",cost(err,label_data))
         rowsum=-err.sum(axis=1)
         rowsum=rowsum.repeat(k,axis=1)
         err=err/rowsum
         for x in range(m):
             err[x,label_data[x,0]]+=1
         weights=weights+(alpha/m)*feature_data.T*err
         i=i+1
     return weights

#%% 计算损失值函数
def cost(err,label_data):
    # input: err(mat):exp的值
     #          label_data:标签的值
      #  output: sum_cost/ m(float):损失函数的值
    
    m=np.shape(err)[0]
    sum_cost=0.0
    for i in range(m):
        if err[i,label_data[i,0]]/np.sum(err[i,:])>0:
            sum_cost -= np.log(err[i,label_data[i,0]]/np.sum(err[i,:]))
        else:
            sum_cost -=0
    return sum_cost/m

训练样本，得到参数θ的矩阵

5. 模型预测

将得到的参数θ代入随机生成的数据进行测试，得到的结果如下图

是不是很精确！学到了吗，那就给小哥顶一个吧，thx

 

数据：

http://pan.baidu.com/s/1hrLUCss