數(shù)學(xué)建模——邏輯回歸模型詳解Python代碼

程序用到的測試數(shù)據(jù)：

鏈接：https://pan.baidu.com/s/1LGD1MAxk2lxO93smSPNyZg
提取碼：uukr

代碼正文

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=['Exam1','Exam2','Admitted']) pdData.head() print(pdData.head()) print(pdData.shape) positive=pdData[pdData['Admitted']==1]#定義正 nagative=pdData[pdData['Admitted']==0]#定義負(fù) fig,ax=plt.subplots(figsize=(10,5)) ax.scatter(positive['Exam1'],positive['Exam2'],s=30,c='b',marker='o',label='Admitted') ax.scatter(nagative['Exam1'],nagative['Exam2'],s=30,c='r',marker='x',label='not Admitted') ax.legend() ax.set_xlabel('Exam 1 score') ax.set_ylabel('Exam 2 score') plt.show()#畫圖 ##實(shí)現(xiàn)算法 the logistics regression 目標(biāo)建立一個(gè)分類器 設(shè)置閾值來判斷錄取結(jié)果 ##sigmoid 函數(shù) def sigmoid(z):return 1/(1+np.exp(-z)) #畫圖 nums=np.arange(-10,10,step=1) fig,ax=plt.subplots(figsize=(12,4)) ax.plot(nums,sigmoid(nums),'r')#畫圖定義 plt.show() #按照理論實(shí)現(xiàn)預(yù)測函數(shù) def model(X,theta):return sigmoid(np.dot(X,theta.T))pdData.insert(0,'ones',1)#插入一列 orig_data=pdData.as_matrix() cols=orig_data.shape[1] X=orig_data[:,0:cols-1] y=orig_data[:,cols-1:cols] theta=np.zeros([1,3]) print(X[:5]) print(X.shape,y.shape,theta.shape) ##損失函數(shù) def cost(X,y,theta):left=np.multiply(-y,np.log(model(X,theta)))right=np.multiply(1-y,np.log(1-model(X,theta)))return np.sum(left-right)/(len(X)) print(cost(X,y,theta))#計(jì)算梯度 def gradient(X, y, theta):grad = np.zeros(theta.shape)error = (model(X, theta) - y).ravel()for j in range(len(theta.ravel())): # for each parmeterterm = np.multiply(error, X[:, j])grad[0, j] = np.sum(term) / len(X)return grad ##比較3種不同梯度下降方法 STOP_ITER=0 STOP_COST=1 STOP_GRAD=2def stopCriterion(type,value,threshold):if type==STOP_ITER: return value>thresholdelif type==STOP_COST: return abs(value[-1]-value[-2])<thresholdelif type==STOP_GRAD: return np.linalg.norm(value)<thresholdimport numpy.random #打亂數(shù)據(jù)洗牌 def shuffledata(data):np.random.shuffle(data)cols=data.shape[1]X=data[:,0:cols-1]y=data[:,cols-1:]return X,yimport timedef descent(data, theta, batchSize, stopType, thresh, alpha):# 梯度下降求解init_time = time.time()i = 0 # 迭代次數(shù)k = 0 # batchX, y = shuffledata(data)grad = np.zeros(theta.shape) # 計(jì)算的梯度costs = [cost(X, y, theta)] # 損失值while True:grad = gradient(X[k:k + batchSize], y[k:k + batchSize], theta)k += batchSize # 取batch數(shù)量個(gè)數(shù)據(jù)if k >= n:k = 0X, y = shuffledata(data) # 重新洗牌theta = theta - alpha * grad # 參數(shù)更新costs.append(cost(X, y, theta)) # 計(jì)算新的損失i += 1if stopType == STOP_ITER:value = ielif stopType == STOP_COST:value = costselif stopType == STOP_GRAD:value = gradif stopCriterion(stopType, value, thresh): breakreturn theta, i - 1, costs, grad, time.time() - init_time #選擇梯度下降 def runExpe(data, theta, batchSize, stopType, thresh, alpha):#import pdb; pdb.set_trace();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 += strStopprint ("***{}\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')return theta n= 100 runExpe(orig_data,theta,n,STOP_ITER,thresh=5000,alpha=0.000001) plt.show() runExpe(orig_data,theta,n,STOP_GRAD,thresh=0.05,alpha=0.001) plt.show() runExpe(orig_data,theta,n,STOP_COST,thresh=0.000001,alpha=0.001) plt.show() #對(duì)比 runExpe(orig_data, theta, 1, STOP_ITER, thresh=5000, alpha=0.001) plt.show() runExpe(orig_data, theta, 1, STOP_ITER, thresh=15000, alpha=0.000002) plt.show() runExpe(orig_data, theta, 16, STOP_ITER, thresh=15000, alpha=0.001) plt.show() ##對(duì)數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化 將數(shù)據(jù)按其屬性(按列進(jìn)行)減去其均值，然后除以其方差。 #最后得到的結(jié)果是，對(duì)每個(gè)屬性/每列來說所有數(shù)據(jù)都聚集在0附近，方差值為1from sklearn import preprocessing as ppscaled_data = orig_data.copy() scaled_data[:, 1:3] = pp.scale(orig_data[:, 1:3])runExpe(scaled_data, theta, n, STOP_ITER, thresh=5000, alpha=0.001) #設(shè)定閾值 def predict(X, theta):return [1 if x >= 0.5 else 0 for x in model(X, theta)]# if __name__=='__main__':scaled_X = scaled_data[:, :3] y = scaled_data[:, 3] predictions = predict(scaled_X, theta) correct = [1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0 for (a, b) in zip(predictions, y)] accuracy = (sum(map(int, correct)) % len(correct)) print ('accuracy = {0}%'.format(accuracy))

總結(jié)

以上是生活随笔為你收集整理的数学建模——逻辑回归模型Python代码的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。