邏輯迴歸是分類最簡單的演算法。

題目是有一個學校招生，考2門專業課。給出了100個樣本，每一個樣本包含2門課的成績以及是否被錄取。也就是LogiReg_data.csv資料。原始是txt格式，為了方便，人為的加上了表頭。

  1 import numpy as np
  2 import pandas as pd
  3 import matplotlib.pyplot as plt
  4 # 匯入相關庫，並重命名
  5 pdData = pd.read_csv("LogiReg_data.csv")  
  6 #開啟檔案
  7 positive = pdData[pdData['
 
Admitted'] == 1] #將'Admitted'一列為1的行（樣本）提取出來，放入positive
  8 negative = pdData[pdData['Admitted'] == 0] #將'Admitted'一列為0的行（樣本）提取出來，放入negative
  9 fig, ax = plt.subplots(figsize=(10,5))
 10 #subplots固定用法
 11 #subplots是matplotlib中的畫子影象的函式，
 12 #原型：
 13 #plt.subplots(
 14 #    nrows=1,
 15 #    ncols=1,
 16 #
 
    sharex=False,
 17 #    sharey=False,
 18 #    squeeze=True,
 19 #    subplot_kw=None,
 20 #    gridspec_kw=None,
 21 #    **fig_kw,
 22 #)
 23 #nrows，ncols控制幾個子圖，例如：nrows=1，ncols=2就是2個子圖，左一個右一個。sharex、sharey控制
 24 #是否共享x、y軸。返回值是兩個，一個是fig，影象，一個是axes。figsize是圖片大小。
 25 ax.scatter(positive['Exam 1'], positive['
 
Exam 2'], s=30, c='b', marker='o', label='Admitted')
 26 ax.scatter(negative['Exam 1'], negative['Exam 2'], s=30, c='r', marker='x', label='Not Admitted')
 27 #scatter可以用來畫散點圖，裡邊引數positive['Exam 1'], positive['Exam 2']分別是x軸和y軸是什麼。“c”是可選顏色
 28 #“marker”是可選點的形狀，具體可以百度。s是size，為點的大小，label是圖中'o'的點是'Admitted'，'x'是'Not Admitted'
 29 ax.legend()
 30 #ax.legend控制圖例的引數，例如大小、位置、標題等等。https://blog.csdn.net/qq_35240640/article/details/89478439?utm_medium=distribute.pc_aggpage_search_result.none-task-blog-2~all~first_rank_v2~rank_v25-2-89478439.nonecase
 31 ax.set_xlabel('Exam 1 Score')
 32 ax.set_ylabel('Exam 2 Score')
 33 #set_ylabel是x軸和y軸的名稱叫什麼。用plt.xlabel和plt.ylabel也可以。
 34 
 35 def sigmoid(z):
 36     return 1 / (1 + np.exp(-z))
 37 #sigmoid 函式,用來對映成【0,1】的一個概率值
 38 
 39 nums = np.arange(-10, 10, step=1) #np.arange構造一個【-10,10】的一個ndarray，步長為1.（10，-9，-8，-7…………，7,8,9,10）
 40 fig, ax = plt.subplots(figsize=(12,4))
 41 ax.plot(nums, sigmoid(nums), 'r')  #ax.plot畫出sigmoid函式。其中sigmoid(0)=0.5
 42 print(type(nums)) #發現nums就是一個ndarray。
 43 
 44 def model(X, theta):
 45     return sigmoid(np.dot(X, theta.T))  #np.dot 計算點積，需要注意的是，點積計算時，要求第一個引數的列數要和第二個引數的行數一樣才可以。最後的返回值也是一個矩陣，列數為1
 46 
 47 pdData.insert(0, 'Ones', 1,allow_duplicates=True)   #在第一列資料前增加新的一列，變為0、1、2、3列
 48 orig_data = pdData.iloc[:,:].values  #將dataframe轉成numpy.array。
 49 cols = orig_data.shape[1]  #計算出二位陣列有幾列，也可以寫成cols=len(orig_data[0])
 50 X = orig_data[:,0:cols-1]  #將前3列0、1、2列賦值給x
 51 y = orig_data[:,cols-1:cols] #將最後一列賦值給y
 52 theta = np.zeros([1, 3])  #建立一個theta  [[0. 0. 0.]]  這裡寫np.zeros((1,3))也可以
 53 
 54 
 55 
 56 def cost(X, y, theta):   #定義損失函式，損失函式可以評判估計值和實際值的誤差
 57     left = np.multiply(-y, np.log(model(X, theta)))  #np.multiply是引數相乘、乘法。np.log計算數學的log，model返回值是計算sigmoid，結果（0,1）
 58     right = np.multiply(1 - y, np.log(1 - model(X, theta)))
 59     return np.sum(left - right) / (len(X))  #np.sum求和，計算總的損失，再除以樣本數量，就是平均損失。
 60 
 61 cost(X, y, theta)
 62 
 63 def gradient(X, y, theta):  #用來計算梯度下降方向的，這裡需要先對損失函式求導，目的是計算出損失最小的位置，參考數學公式。
 64     grad = np.zeros(theta.shape)  #grad是一個二維陣列，[[0. 0. 0.]]
 65     error = (model(X, theta)- y).ravel() #把error轉成一維陣列，是為後邊的和X第j列相乘做準備
 66     for j in range(len(theta.ravel())):  #theta本身是二維陣列，但是通過theta.ravel() 將theta轉成一維陣列， len(theta.ravel()) 計算陣列長度，在這裡實際長度是3。
 67         term = np.multiply(error, X[:,j])  #注意X是一個二維陣列，所以採用X[:,j]選取所有行的第j列。
 68         grad[0, j] = np.sum(term) / len(X)  #grad是一個二維陣列，[[0. 0. 0.]]，theta有三個引數，所以要計算三個梯度的方向
 69     return grad  #grad是一個二維陣列，[[0. 0. 0.]]
 70 
 71 STOP_ITER = 0
 72 STOP_COST = 1
 73 STOP_GRAD = 2
 74 #設定三種不同的停止策略
 75 def stopCriterion(type, value, threshold):
 76     if type == STOP_ITER:        return value > threshold  #設定迭代次數
 77     elif type == STOP_COST:      return abs(value[-1]-value[-2]) < threshold  #給損失函式cost設定閾值
 78     elif type == STOP_GRAD:      return np.linalg.norm(value) < threshold  #給梯度設定閾值  注意三種情況返回的都是bool值
 79 
 80 import numpy.random
 81 #洗牌，無法保證資料是具有隨機性，通過np.random.shuffle函式打亂資料順序，使資料具有隨機性
 82 def shuffleData(data):  #函式的功能就是將X和y資料重新排列組合
 83     np.random.shuffle(data)
 84     cols = data.shape[1]  #計算出二位陣列有幾列，也可以寫成cols=len(orig_data[0])
 85     X = data[:, 0:cols-1]  #重新選擇X，但是和之前一樣，還是第0、1、2列
 86     y = data[:, cols-1:]  #重新選擇y，但是和之前一樣，還是第3列
 87     return X, y
 88 
 89 
 90 import time
 91 
 92 
 93 def descent(data, theta, batchSize, stopType, thresh, alpha):
 94     # 梯度下降求解
 95 
 96     init_time = time.time()  # 顯示計算時間
 97     i = 0  # 迭代次數  迭代次數越多，損失函式越趨於0，看影象就可以
 98     k = 0  # batch  選取的樣本，當k=1時，一個樣本計算一次，計算速度快，但不一定收斂，這是隨機梯度下降；
 99     # 但是慢當k=10時，10個樣本計算一次，比較實用，這是小批量梯度下降；當k=樣本數量時，計算所有樣本，速度慢，但是容易找到最優解，這是批量梯度下降。
100     X, y = shuffleData(data)  # 呼叫shuffleData函式，打亂資料
101     grad = np.zeros(theta.shape)  # 初始化grad，grad是梯度
102     costs = [cost(X, y, theta)]  # 初始化損失值，這裡打算用appand，所以新生成的costs是二維陣列。
103     # 以上是將變數初始化。
104     while True:  # 死迴圈，碰到break停止
105         grad = gradient(X[k:k + batchSize], y[k:k + batchSize],
106                         theta)  # X和y是二維陣列X[k:k+batchSize]，意思是第k行到第k+batchSize行，y同理
107         k += batchSize  # 取batch數量個數據
108         if k >= n:  # n會在後邊賦值，這個n代表如果計算的樣本數量超出了原來的樣本數量，就重新洗牌，重新計算引數
109             k = 0
110             X, y = shuffleData(data)  # 重新洗牌
111         theta = theta - alpha * grad  # 引數更新，可以參考數學部分
112         costs.append(cost(X, y, theta))  # 計算新的損失
113         i += 1
114 
115         if stopType == STOP_ITER:
116             value = i
117         elif stopType == STOP_COST:
118             value = costs
119         elif stopType == STOP_GRAD:
120             value = grad
121         if stopCriterion(stopType, value, thresh): break
122 
123     return theta, i - 1, costs, grad, time.time() - init_time
124 
125 def runExpe(data, theta, batchSize, stopType, thresh, alpha):
126     #import pdb; pdb.set_trace();
127     theta, iter, costs, grad, dur = descent(data, theta, batchSize, stopType, thresh, alpha)
128 #以下都是為圖的標題做準備的
129     name = "Original" if (data[:,1]>2).sum() > 1 else "Scaled"
130     name += " data - learning rate: {} - ".format(alpha)
131     if batchSize==n: strDescType = "Gradient"
132     elif batchSize==1:  strDescType = "Stochastic"
133     else: strDescType = "Mini-batch ({})".format(batchSize)
134     name += strDescType + " descent - Stop: "
135     if stopType == STOP_ITER: strStop = "{} iterations".format(thresh)
136     elif stopType == STOP_COST: strStop = "costs change < {}".format(thresh)
137     else: strStop = "gradient norm < {}".format(thresh)
138     name += strStop
139     print ("***{}\nTheta: {} - Iter: {} - Last cost: {:03.2f} - Duration: {:03.2f}s".format(
140         name, theta, iter, costs[-1], dur))
141 #以上都是為圖的標題做準備的
142     fig, ax = plt.subplots(figsize=(12,4))
143     ax.plot(np.arange(len(costs)), costs, 'r')
144     ax.set_xlabel('Iterations')
145     ax.set_ylabel('Cost')
146     ax.set_title(name.upper() + ' - Error vs. Iteration')
147     plt.show()
148     return theta
149 n=100
150 runExpe(orig_data, theta, n, STOP_ITER, thresh=5000, alpha=0.000001)  #1.1次
151 #一次選取整個樣本，採用控制迭代次數的方式，5000次，Last cost: 0.63,time=1.69s
152 runExpe(orig_data, theta, n, STOP_ITER, thresh=50000, alpha=0.000001)  #1.2次
153 #一次選取整個樣本，採用控制迭代次數的方式，5000次，Last cost: 0.63,time=17.30s
154 runExpe(orig_data, theta, n, STOP_ITER, thresh=5000, alpha=0.001)   #1.3次
155 #一次選取整個樣本，採用控制迭代次數的方式，5000次，Last cost: 0.61,time=1.69s
156 #1.2和1.1次相比，迭代次數*10，時間也大約*10，但是損失值並沒有降低，所以這種情況下增加迭代次數並沒有用。
157 #1.1和1.3相比，迭代次數不變，學習率降低，但是cost並沒有降低很多，一般的cost取0.001比較合適
158 runExpe(orig_data, theta, n, STOP_COST, thresh=0.000001, alpha=0.001)  #2次
159 #一次選取整個樣本，採用控制cost方式，Last cost: 0.38，time: 38.37s，需要迭代11000次左右
160 #通過增大或者減小cost的值很容易發現，cost越小，迭代次數越多，消耗時間越多；cost增大則相反。
161 runExpe(orig_data, theta, 100, STOP_GRAD, thresh=0.02, alpha=0.001) # 3.1次
162 #一次選取整個樣本，採用控制grad梯度下降速度方式，Last cost: 0.31，time: 78.83s，需要迭代220000次左右
163 runExpe(orig_data, theta, 100, STOP_GRAD, thresh=0.05, alpha=0.001) # 3.2次
164 #一次選取整個樣本，採用控制grad梯度下降速度方式，Last cost: 0.49，time: 14.96s，需要迭代40000次左右
165 #3.2和3.1對比，增大grad，發現迭代次數降低，cost升高。實際上grad為0時為最理想值，但是耗時非常大。
166 runExpe(orig_data, theta, 1, STOP_ITER, thresh=5000, alpha=0.001)  #4次
167 #一次選取1個樣本，採用控制迭代次數方式，影象不穩定，說明1個樣本下模型不收斂，也就是結果不趨於0
168 #即使增加迭代次數也沒用,減小學習率alpha也沒用，影象始終不穩定。
169 runExpe(orig_data, theta, 16, STOP_ITER, thresh=15000, alpha=0.001)  #5次
170 #一次選取16個樣本，採用控制迭代次數方式，影象仍然不穩定，所以在控制迭代次數方式這裡，問題出在了資料
171 #需要對資料進行標準化
172 
173 #以下是資料進行標準化
174 from sklearn import preprocessing as pp
175 scaled_data = orig_data.copy()
176 scaled_data[:, 1:3] = pp.scale(orig_data[:, 1:3])
177 #以上是資料進行標準化
178 runExpe(scaled_data, theta, n, STOP_ITER, thresh=5000, alpha=0.001)  #6次
179 #一次選取整個樣本，採用控制迭代次數的方式，5000次，Last cost: 0.38，time:1.92s
180 #6次和1.1次相比，是資料標準化之後的結果，發現時間上沒有太大改善（因為迭代次數沒有變），但是cost下降了很多
181 #所以資料預處理很重要，要先進行標準化，去量綱，變為純數
182 runExpe(scaled_data, theta, n, STOP_GRAD, thresh=0.02, alpha=0.001)  #7次
183 #一次選取整個樣本，採用控制grad梯度下降速度方式，迭代次數不到60000次，Last cost: 0.22 - time: 22.24s
184 #7次和3.1相比，是資料標準化之後的結果，迭代次數減小，時間減小，cost降到了0.22
185 runExpe(scaled_data, theta, 1, STOP_GRAD, thresh=0.02, alpha=0.001)  #8次
186 #一次選取1個樣本，採用控制grad梯度下降速度方式，迭代次數15000次左右，cost: 0.28，time：1.95s
187 #8次和7次相比，樣本數為1，整體迭代次數縮減了很多，時間縮減了很多
188 
189 theta = runExpe(scaled_data, theta, 1, STOP_GRAD, thresh=0.002/5, alpha=0.001) #9次
190 #一次選取1個樣本，採用控制grad梯度下降速度方式，70000多次，cost降到了0.22 - time: 8.78s
191 runExpe(scaled_data, theta, 16, STOP_GRAD, thresh=0.002*2, alpha=0.001)  #10次
192 #一次選取16個樣本，採用控制grad梯度下降速度方式，6000多次，cost降到了0.21 - time: 1.19s
193 #需要注意的是10次和9次比，雖然看起來只是變了計算一次的樣本數量，但是10的theta初試值並不是0，
194 #而是用的9次的計算結果，也就是說初始值是經過訓練之後的值，所以迭代次數大幅下降。
195 
196 def predict(X, theta):  #函式功能是用原始資料驗證模型的準確性
197     #return [1 if x >= 0.5 else 0 for x in model(X, theta)]  #這是原博主寫的，但是不太好理解，但是習慣矩陣運算很重要
198     a = []
199     for i in model(X, theta):
200         if i>=0.5:
201             a.append(1)
202         else:
203             a.append(0)
204     return a #函式的返回值是一個list，是每一個樣本是否被錄取的結果[1, 0, 1, 1, 0,…………, 1, 1, 1, 0]
205     
206     scaled_X = scaled_data[:, :3]  #選取前三列資料
207 y = scaled_data[:, 3]  #選取最後一列資料
208     predictions = predict(scaled_X, theta)  #predictions是預測值
209 correct = [1 if ((a == 1 and b == 1) or (a == 0 and b == 0)) else 0 for (a, b) in zip(predictions, y)]
210 #correct是一個list，如果真實值和預測值相等，該位置寫1，否則寫0.這裡也可以用for迴圈實現。
211 #zip函式是將(predictions, y)組合成一個可迭代物件，不然for迴圈無法實現。
212 accuracy = (sum(map(int, correct)) % len(correct))  #map是將correct轉成int型list（猜的），可以不用map，分開計算也可以
213 print ('accuracy = {0}%'.format(accuracy))  #計算概率，需要注意的是，這裡是整數除以整數。format很簡單，百度就可以

數學公式：

以上是根據唐宇迪視訊，以及某部落格https://www.cnblogs.com/douzujun/p/8370993.html所寫

程式碼是基於python 3.7，沒有任何問題。