逻辑回归小案例

先又一组学生数据,要求通过这个学生的两门成绩来预测这个学生是否能被录取。

数据格式:




先引入三件套:


import numpy as np import pandas as pd import matplotlib.pyplot as plt import
os 按照表格样式打印出数据,方便可视化:



path = 'data' + os.sep + 'LogiReg_data.txt' pdData = pd.read_csv(path,header =
None ,names = ['Exam 1','Exam 2','admin']) pdData.head()


admin就是是否被录取。


接着画出散点图:


positive = pdData[pdData['admin'] == 1] negative = pdData[pdData['admin'] ==
0] fix,ax = plt.subplots(figsize = (10,5)) ax.scatter(positive['Exam
1'],positive['Exam 2'],s = 30,c = 'b',marker = 'o',label = 'Admitted')
ax.scatter(negative['Exam 1'],negative['Exam 2'],s = 30,c = 'r',marker =
'x',label = 'Not Admintted') ax.legend() ax.set_xlabel('Exam 1 score')
ax.set_ylabel('Exam 2 socre') plt.show()






接下来就是分析的步骤了:

logistic regression

目标:建立一个分类器。求解出三个参数,θ0  θ1  θ2,θ0是一个偏移量,其他的就是两个特征值所以对应两个参数。

设定一个阈值,比如大于百分之50就录取,小于就不录取。

1,Sigmoid,一个可以把值映射成概率的函数。

2,model,返回函数预测结果值。

3,cost,根据参数计算损失。

4,gradient,根据每个参数计算梯度方向。

5,descent,进行参数更新。

6,accuracy,计算精度

Sigmoid函数模块:


def sigmoid(z):     return 1 / (1+np.exp(-z)) 可以画出Sigmoid函数的图像:



nums = np.arange(-10,10,step = 1) fig,ax = plt.subplots(figsize = (12,4))
ax.plot(nums,sigmoid(nums),'r') plt.show()
model模块:


model模块是用于计算参数和特征值的乘积的,也就是构建一个模型。


def model(X,theta): return sigmoid(np.dot(X,theta.T))
np.dot函数就是矩阵乘法,theta.T就是矩阵的转置。



建立特征值X,参数theta,Y的矩阵:


pdData.insert(0,'ones',1) matrix = pdData.as_matrix() cols = matrix.shape[1] X
= matrix[:,0:cols-1]#特征值 Y = matrix[:,cols-1:cols]#真实值 theta = np.zeros([1,3])
因为theta0是一个偏置值,所以没有X对应,所以在X里面添加一行作为X0和theta0对应,就设为1。然后把这个表格转化为一个矩阵,因为这个矩阵的最后一行就是真实值,也就是分割的值(逻辑回归就是适合二分类的算法),所以取前一行。theta是参数值,参数值是未知的,所以可以随机,但为了简单所以取0。


cost损失函数:



把负号提进去,先算里面的。



def cost(X,y,theta): left = np.multiply(-y,np.log(model(X,theta)))
print(left[1:5]) right = np.multiply(1-y,np.log(1-model(X,theta)))
print(right[1:5]) return np.sum(left-right)/(len(X))


计算梯度:

求偏导之后的最终结果。


def gradient(X,y,theta): grad = np.zeros.shape(theta.shape) error =
np.multiply(model(X,theta) - y).ravel() for j in range (len(theta.ravel)): term
= np.multiply(error,X[:,j]) grad[0,j] = np.sum(term)/(len(X)) return grad
先建立一个结果grad,大小是和参数的theta是一样的。把里面的负号提进去,就可以求得里面的。ravel函数是把一个矩阵编程链表,这个矩阵本来就是一维的,所以直接就变成了一个链表。Xij是第i个样本j个特征值。而样本是一起乘上去的,就是使用矩阵乘法自己乘,所以只需要取j就可以了。





迭代停止条件:

因为theta值的计算是没有终点的,所以要设置一个计算终点。有三种方法:1,按照迭代次数。2,按照损失函数的精度。3,按照梯度下降的精度。


STOP_ITER = 0 STOP_COST = 1 STOP_GRAD = 2 def
stopCriterion(type,value,threshold): if type == STOP_ITER: return value >
threshold elif type == STOP_COST: return abs(value[-1]-value[-2])<threshold
elif type == STOP_GRAD: np.linalg.norm(value) < threshold
np.linalg.norm()是求范式的意思,如果不加参数就是默认。



洗牌:


import numpy.random #洗牌 def shuffleData(data): np.random.shuffle(data) cols =
data.shape[1] X = data[:, 0:cols-1] y = data[:, cols-1:] print('X',X[1:5])
print('Y',Y[1:5]) return X, y使用np.random.shuffle(矩阵)来洗牌,再重新获取函数。



梯度下降求解:
每一次都要洗牌后再迭代,然后计算梯度下降的幅度,学习率乘于幅度就是每个增长的参数值了。
import time def
descent(data,theta,batchsize,stopType,thresh,alpha):#数据,参数,梯度下降的类型,停止的类型,停止的边缘,学习率
init_time = time.time() i = 0#迭代次数 k = 0#每一次梯度下降的数据 X,y = shuffleData(data)
grad = np.zeros(theta.shape) costs = [cost(X,y,theta)] while True: grad =
gradient(X[k:k+batchsize],y[k:k+batchsize],theta)#梯度下降 k += batchsize if k >=
100:#大于总数据 k = 0 X,y = shuffleData(data)#重新再来 theta = theta - alpha * grad
costs.append(cost(X,y,theta)) i += 1 if stopType == STOP_ITER: value = i elif
stopType == STOP_COST: value = costs elif stopType == STOP_GRAD: value = grad
if stopCriterion(stopType,value,thresh): break return
theta,i-1,costs,grad,time.time() - init_time执行函数:

def runExpe(data,theta,batchSize,stopType,thresh,alpha): n = 100
theta,iter,costs,grad,dur =
descent(data,theta,batchSize,stopType,thresh,alpha)#核心函数 name = "Original" if
(data[:,1]>2).sum() > 1 else "Scaled" name += " data - learning rate: {} -
".format(alpha) if batchSize==n: strDescType = "Gradient" elif batchSize==1:
strDescType = "Stochastic" else: strDescType = "Mini-batch
({})".format(batchSize) name += strDescType + " descent - Stop: " if stopType
== STOP_ITER: strStop = "{} iterations".format(thresh) elif stopType ==
STOP_COST: strStop = "costs change < {}".format(thresh) else: strStop =
"gradient norm < {}".format(thresh) name += strStop print ("***{}\nTheta: {} -
Iter: {} - Last cost: {:03.2f} - Duration: {:03.2f}s".format( name, theta,
iter, costs[-1], dur)) fig, ax = plt.subplots(figsize=(12,4))
ax.plot(np.arange(len(costs)), costs, 'r') ax.set_xlabel('Iterations')
ax.set_ylabel('Cost') ax.set_title(name.upper() + ' - Error vs. Iteration')
plt.show() return theta只有descent求梯度下降才是核心,其他都是画图用的。


可以看得出损失函数到达0.63就维持不变了,三个参数也给了出来。

接下来调换停止策略:

我们使用按照损失函数的精度来进行测试:


result = runExpe(matrix, theta, 100,STOP_COST, thresh=0.000001, alpha=0.001)



哎  电脑垃圾,花了249秒,放到同学电脑里20秒。。。。

可以看到不同的策略就有不同的效果了,比之前的好很多,到0.4了。

当然其实也可以调换其他的策略,比如梯度下降的策略,改变为小批量等等。