逻辑回归算法原理及其python实现

1.逻辑回归原理

1.决策边界

梯度下降要做的就是在使损失函数尽量小的情况下求的一组Θ

2.损失函数

这里只给出了损失函数的梯度下降求解，之所以使用对数损失，其实是一个极大似然估计的参数估计问题



为何此处损失函数要用对数损失？
https://www.jianshu.com/p/b6bb6c035d8c

2.线性逻辑回归python实现

数据集展示

from matplotlib import pyplot as plt
import numpy as np

# 载入数据
data = np.genfromtxt("LR-testSet.csv", delimiter=",")
x_data = data[:, :-1]
y_data = data[:, -1]

def plot():
    x0 = []
    x1 = []
    y0 = []
    y1 = []
    # 切分不同类别的数据
    for i in range(len(x_data)):
        if y_data[i] == 0:
            x0.append(x_data[i, 0])
            y0.append(x_data[i, 1])
        else:
            x1.append(x_data[i, 0])
            y1.append(x_data[i, 1])

    # 画图
    scatter0 = plt.scatter(x0, y0, c='b', marker='o')
    scatter1 = plt.scatter(x1, y1, c='r', marker='x')
    # 画图例
    plt.legend(handles=[scatter0, scatter1], labels=['label0', 'label1'], loc='best')

plot()
plt.show()

逻辑回归

# 数据处理，添加偏置项
x_data = data[:,:-1]
y_data = data[:,-1,np.newaxis]

print(np.mat(x_data).shape) # (100, 2)
print(np.mat(y_data).shape) # (100, 1)

# 给样本添加偏置项
# x0=1 x1 x2
x_data = np.concatenate((np.ones((100,1)),x_data),axis=1)
print(x_data.shape) # (100, 3)

def sigmoid(x):
    return 1.0/(1+np.exp(-x))

# 损失函数
def cost(xMat, yMat, ws):
    # y*log(h(x))
    left = np.multiply(yMat, np.log(sigmoid(xMat*ws)))
    # (1-y)*log(h(x))
    right = np.multiply(1 - yMat, np.log(1 - sigmoid(xMat*ws)))
    return np.sum(left + right) / -(len(xMat))

# 梯度下降
def gradAscent(xArr, yArr):
    xArr = preprocessing.StandardScaler().fit_transform(xArr)
    xMat = np.mat(xArr)
    yMat = np.mat(yArr)
    
    lr = 0.001
    epochs = 10000
    costList = []
    # 计算数据行列数
    # 行代表数据个数，列代表权值个数
    m,n = np.shape(xMat)
    # 初始化权值 (3,1) Θ1 Θ2 Θ3
    ws = np.mat(np.ones((n,1)))
    
    for i in range(epochs+1):             
        # xMat和weights矩阵相乘 (100,3)*(3,1)=(100,1)
        h = sigmoid(xMat*ws)   
        # 计算误差   (3,100)*(100,1)=(3,1) 得到Θ1 Θ2 Θ3的梯度
        ws_grad = xMat.T*(h - yMat)/m
        # 重新计算梯度
        ws = ws - lr*ws_grad 
        
        if i % 50 == 0:
            costList.append(cost(xMat,yMat,ws))
    return ws,costList

# 训练模型，得到权值和cost值的变化
ws,costList = gradAscent(x_data, y_data)
print(ws)
# [[ 1.        ]
#  [ 0.30816757]
#  [-1.76171512]]

# 预测
def predict(x_data, ws):
    x_data = preprocessing.StandardScaler().fit_transform(x_data)
    xMat = np.mat(x_data)
    ws = np.mat(ws)
    return [1 if x >= 0.5 else 0 for x in sigmoid(xMat*ws)]
predictions = predict(X_data, ws)
print(classification_report(y_data, predictions))

#               precision    recall  f1-score   support

#          0.0       0.92      1.00      0.96        47
#          1.0       1.00      0.92      0.96        53

#     accuracy                           0.96       100
#    macro avg       0.96      0.96      0.96       100
# weighted avg       0.96      0.96      0.96       100

3.非线性逻辑回归python实现

数据集展示

import matplotlib.pyplot as plt
import numpy as np
from sklearn.metrics import classification_report
from sklearn import preprocessing
from sklearn.preprocessing import PolynomialFeatures

# 载入数据
data = np.genfromtxt("LR-testSet2.txt", delimiter=",")
x_data = data[:,:-1]
y_data = data[:,-1,np.newaxis]
    
def plot():
    x0 = []
    x1 = []
    y0 = []
    y1 = []
    # 切分不同类别的数据
    for i in range(len(x_data)):
        if y_data[i]==0:
            x0.append(x_data[i,0])
            y0.append(x_data[i,1])
        else:
            x1.append(x_data[i,0])
            y1.append(x_data[i,1])

    # 画图
    scatter0 = plt.scatter(x0, y0, c='b', marker='o')
    scatter1 = plt.scatter(x1, y1, c='r', marker='x')
    #画图例
    plt.legend(handles=[scatter0,scatter1],labels=['label0','label1'],loc='best')
    
plot()
plt.show()

# 定义多项式回归,degree的值可以调节多项式的特征
poly_reg  = PolynomialFeatures(degree=3) 
# 特征处理
x_poly = poly_reg.fit_transform(x_data)
print(x_poly.shape) # (118,10)
# x0=1 x1 x2 x1^2 x1x2 x2^2 x1^3 x1^2*x2 x2^2*x1 x2^3

def sigmoid(x):
    return 1.0/(1+np.exp(-x))

# 损失函数
def cost(xMat, yMat, ws):
    left = np.multiply(yMat, np.log(sigmoid(xMat*ws)))
    right = np.multiply(1 - yMat, np.log(1 - sigmoid(xMat*ws)))
    return np.sum(left + right) / -(len(xMat))

# 梯度下降
def gradAscent(xArr, yArr):
    xArr = preprocessing.StandardScaler().fit_transform(xArr)
    xMat = np.mat(xArr)
    yMat = np.mat(yArr)
    
    lr = 0.03
    epochs = 50000
    costList = []
    # 计算数据列数，有几列就有几个权值
    m,n = np.shape(xMat)
    # 初始化权值 (10,1) Θ1 Θ2...Θ10
    ws = np.mat(np.ones((n,1)))
   
    for i in range(epochs+1):             
        # xMat和weights矩阵相乘 (118,10)*(10,1)=(118,1)
        h = sigmoid(xMat*ws)
        # 计算误差 (10,118)*(118,1)=(10,1) 得到Θ1 Θ2...Θ10的梯度
        ws_grad = xMat.T*(h - yMat)/m
        ws = ws - lr*ws_grad 
        
        if i % 50 == 0:
            costList.append(cost(xMat,yMat,ws))
    return ws,costList

# 训练模型，得到权值和cost值的变化
ws,costList = gradAscent(x_poly, y_data)
print(ws)
# [[ 4.16787292]
#  [ 2.72213524]
#  [ 4.55120018]
#  [-9.76109006]
#  [-5.34880198]
#  [-8.51458023]
#  [-0.55950401]
#  [-1.55418165]
#  [-0.75929829]
#  [-2.88573877]]

# 预测
def predict(x_data, ws):
    x_data = preprocessing.StandardScaler().fit_transform(x_data)
    xMat = np.mat(x_data)
    ws = np.mat(ws)
    return [1 if x >= 0.5 else 0 for x in sigmoid(xMat*ws)]
predictions = predict(x_poly, ws)
print(classification_report(y_data, predictions))

#               precision    recall  f1-score   support

#          0.0       0.89      0.83      0.86        60
#          1.0       0.84      0.90      0.87        58

#     accuracy                           0.86       118
#    macro avg       0.87      0.86      0.86       118
# weighted avg       0.87      0.86      0.86       118