NKUMachineLearning-lab2.2-非线性logistics回归

题设

题设存于文章所附带的pdf文件中，数据集不予提供。

---

代码实现

# Lab.2 Logistic回归算法(非线性)
## 1.读取数据
import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as op
#使用loadtxt函数读取数据存于变量ex2data2
#使用变量X储存ex2data2的前两列数据（芯片的两项测试结果）
#使用变量y储存ex2data2的第三列数据（标签，1表示通过测试，0表示不能通过）
#使用变量m储存样本数量
ex2data2=np.loadtxt('D:\Softwares\Git\ML-2024Spring\Lab2-2-Logistic-with-reg\ex2data2.txt',delimiter=',')
X=ex2data2[:,0:2]
y=ex2data2[:,[2]]
m=y.shape[0]
print(X.shape, y.shape)
## 2.可视化数据
# 使用plt.plot（）函数绘制散点图，使用不同颜色绘制正例和负例
plt.xlabel('Microchip Test 1')
plt.ylabel('Microchip Test 2')
for i in range(m):
    if(y[i]==1):
        plt.plot(X[i,0],X[i,1],'xr')
    elif(y[i]==0):
        plt.plot(X[i,0],X[i,1],'ob')

plt.show()
## 3.训练Logistic模型

### 3.1数据预处理：多项式特征构造，初始化θ，初始化正则化参数lambda
#完成函数mapFeature，输入特征值x1,x2后，返回一个数组,该数组存储维度扩增后的特征值，特征最高次为6次
#将X的前两列分别作为参数X1,X2传入mapFeatrure，并用变量x接收完成构造的特征值数组
#使用np.zeros初始化theta,规模为特征数量*1（特这个数量为np.shape(x)[1]）
#初始化正则化参数lbd为1
def mapFeature(X1,X2):
    degree = 6 # 每个Featuer的最高次
    x_mapped=[]
    for i in range(0,degree+1):
        for j in range (i+1):
            x_mapped.append((X1**(i-j)*X2**j))
    return x_mapped

X1=X[:,0]
X2=X[:,1]
x=mapFeature(X1,X2)
x=np.array(x)
x=x.T
n=x.shape[1]
theta=np.zeros((n,1))
lbd=1
print(x.shape, theta.shape)
### 3.2定义损失函数、梯度
# 完成正则化的损失函数costFunction(),返回值为一维数组
# 完成正则化的梯度计算gradient()，返回值为一维数组
def sigmoid(z):
    g = np.zeros(np.shape(z))
    g = 1/(1+np.exp(-1*z))
    return g

def costFunction(theta, X, y, lbd):
    theta = theta.reshape(len(theta), 1)# 为适应后面训练过程的维度要求，可先忽略
    J=0
    first=np.multiply(-y,np.log(sigmoid(X@theta)))
    second=np.multiply((1-y),np.log(1-sigmoid(X@theta)))
    reg=0

    for i in range(theta.shape[0]):
        reg+=theta[i]**2
    
    reg=reg*lbd/(2*len(X))
    J=np.sum(first-second)/(len(X))+reg
    J = J.flatten()
    return J

def gradient(theta, X, y, lbd):
    theta = theta.reshape(len(theta), 1)# 为适应后面训练过程的维度要求，可先忽略
    # grad必须赋予一个
    m=len(X)
    grad=np.zeros(theta.shape)

    for i in range (m):
        xi = X[i].reshape(1, -1)
        yi = y[i]
        grad += xi.T * (sigmoid(xi @ theta) - yi)
        
    grad=grad/len(X)
    reg=theta*lbd/len(X)

    grad=grad+reg

    grad = grad.flatten() # 返回一维数组
    return grad

cost = costFunction(theta, x, y, lbd)
# print(cost)
grad = gradient(theta, x, y, lbd)
print('对初始零向量theta求得的cost为',cost,'\n梯度为:\n',grad)
### 3.3使用scipy.optimize.minimize求最小损失对应的参数theta

在Python中使用Scipy包下的scipy.optimize.minimize(fun, x0, args=(), method, jac)方法参数的维度要求很严格:  
1.fun为进行优化的目标函数，传入需调用的函数名（不需要加括号），此处为fun=costFunction。需注意调用的函数第一个参数（theta）和返回值（J）必须为一维数组  
2.x0即theta需传入一维数组  
3.args传入fun需要的其他参数，需用tuple传入  
4.method指定优化算法，此处我们使用method='TNC'  
5.jac调用梯度计算函数传入参数需与fun调用函数完全相同，且返回值为一维数组  
高维数组a调整为一维可以使用a.flatten() ，该函数会产生一个副本，不会直接改变a的维度
# 验证参数维度是否符合要求：代价函数和梯度函数的返回值为一维数组、传入x0的参数theta为一维数组
print(cost.ndim, grad.ndim,np.ndim(theta.flatten()))
theta=theta.flatten()
result=op.minimize(fun=costFunction,x0=theta.flatten(),args=(x,y,lbd),method='TNC',jac=gradient)
print(result)
theta_star = result.x
## 4.评估Logistic回归模型

### 4.1绘制决策边界

plt.xlabel('Microchip Test 1')
plt.ylabel('Microchip Test 2')
for i in range(m):
    if y[i] == 1: #正例
        plt.plot(X[i,0], X[i,1], 'or')
    elif y[i] == 0: #反例
        plt.plot(X[i,0], X[i,1], 'ob')

grid_x1 = np.linspace(-1, 1.5, 50)
grid_x2 = np.linspace(-1, 1.5, 50)
x_theta = np.zeros((grid_x1.size,grid_x2.size))
for i in range(grid_x1.size):
    for j in range(grid_x2.size):
        x_theta[i][j] = np.dot(mapFeature(grid_x1[i], grid_x2[j]), theta_star)
x_theta = x_theta.T
plt.contour(grid_x1, grid_x2, x_theta,[0], colors='green' )

plt.savefig('test2.svg', dpi=300, format='svg')

### 4.2计算模型准确率
count=0
for i in range(m):
    res=sigmoid(theta_star@x[i])
    if res>0.5:
        if(y[i]==1):
            count+=1
    if res<0.5:
        if(y[i]==0):
            count+=1

print(count/m)