#coding=utf-8
import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import classification_report#逻辑回归结果评估
from sklearn.preprocessing import PolynomialFeatures
from sympy.logic.algorithms.dpll2 import Level
#用*替换空格和^
def replaceSpace(string,char,char1):
if string == None:return None
if len(string) == 0:return 0
return string.replace(char,char1)
#获取多项式方程。
def get_func(arr,str1):
z=""
for i in range(len(str1)):
if i==len(str1)-1:
z+=str(arr[i])+"*"+str(str1[i])
else:
z+=str(arr[i])+"*"+str(str1[i])+"+"
new_z = replaceSpace(z," ","*")
return replaceSpace(new_z,"^","**")
#关闭科学计数法
np.set_printoptions(suppress=True)
#两个特征的二元分类
data = np.loadtxt('logistic_2d_data.txt',dtype=np.float)
model = LogisticRegression()
model.fit(data[:,:2],data[:,2])
pre_data = model.predict(data[:,:2])
#画出图形
raw_data1 = data[data[:,2]==1.0]#1类别的特征
raw_data0 = data[data[:,2]==0.0]#0类别的特征
plt.scatter(raw_data1[:,0],raw_data1[:,1],c='r',label='feature1',marker="o")
plt.scatter(raw_data0[:,0],raw_data0[:,1],c='yellow',label='feature2',marker="*")
plt.xlabel("x0")
plt.ylabel("x1")
#分类后的图形
#numpy.ogrid()函数生成一定范围一定个数,可用来广播计算的数组
#还需注意x1,x0的顺序。生成的x1是100*1 生成的x0是1*100 两个结构是不一样的。x0代表x轴x1代表y轴
x1,x0 = np.ogrid[0:15:100j,-5:3:100j]
#也可以如下生成
#x0 = np.linspace(-5,3,100).reshape(1,-1)
#x1 = np.linspace(0,15,100).reshape(-1,1)
#print x1
#这个模型只有两个特征,所以关系式是 f(x0,x1) = w0*x0+w1*x1+c
#f(x0,x1)以0为界限,把样本分为两类。所以在接下来的图中直线方程为 f(x0,x1)=0
#我们该如何得到这种画隐函数的图像呢?可以利用绘制等高线的函数来绘制。
#因为f(x0,x1)=0 。是高维函数f(x0,x1)在函数值为0处的等高线图
#下面的get_func()函数就是获取函数解析式f(x0,x1) = w0*x0+w1*x1
#第一个参数是自变量的系数,第二个参数自变量的名称列表。
#因为get_func()函数得到的只是一个字符串,所以需要eval转换为可执行的语句。
#这里得到的z是100*100的数组,而不是100*1的数组,因为前面说过了 x0,x1是可以广播的。
z = eval(get_func(model.coef_[0],["x0","x1"]))+model.intercept_
#注意extent的范围要和x0,x1保持一致。
extent=[-5,3,0,15]
#levels=[0]表示画出等高线为零的线。
plt.contour(z,extent=extent,levels=[0],colors=['black'])
print classification_report(data[:,2],pre_data)+"\n"
print "参数: "+str(model.coef_)
print "截距: "+str(model.intercept_)
#二阶回归目标函数的形式为 f(x0,x1) =
ploy2 = PolynomialFeatures(2,include_bias=False)
#构造特征
new_x = ploy2.fit_transform(data[:,:2])
model.fit(new_x,data[:,2])
#ploy2.get_feature_names()函数返回 ['x0', 'x1', 'x0^2', 'x0 x1', 'x1^2']
#这个就是二阶新特征。原来的一阶特征为['x0','x1']
#这些特征需要传入函数get_func()得到函数解析式.
#print ploy2.get_feature_names()
z2=eval(get_func(model.coef_[0],ploy2.get_feature_names()))+model.intercept_
plt.contour(z2,extent=extent,levels=[0])
plt.xlim(-3,3)
plt.legend(loc="upper right")
plt.show()
运行结果