excel二元回归分析怎么做

题目

在训练的第二部分，我们将要通过加入正则项提升逻辑回归算法。简而言之，正则化是成本函数中的一个术语，它使算法更倾向于“更简单”的模型（在这种情况下，模型将更小的系数）。这个理论助于减少过拟合，提高模型的泛化能力。这样，我们开始吧。

设想你是工厂的生产主管，你有一些芯片在两次测试中的测试结果。对于这两次测试，你想决定是否芯片要被接受或抛弃。为了帮助你做出艰难的决定，你拥有过去芯片的测试数据集，从其中你可以构建一个逻辑回归模型。

预备知识

（1）特征映射

（2）正则化损失函数

（3）正则化梯度下降

答案

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import scipy.optimize as opt

#画出散点图，观察是否呈线性
def plot_data():
    positive = data2[data2['Accepted'].isin([1])]   # 将Accepted列中是1的数据提出来，单独呈一个表格
    negative = data2[data2['Accepted'].isin([0])]   # 同上，提出Accepted是0的表格
    # 接下来画图
    fig, ax = plt.subplots(figsize=(8, 5))
    ax.scatter(positive['Test 1'], positive['Test 2'], s=50, c='b', marker='o', label='Accepted')
    ax.scatter(negative['Test 1'], negative['Test 2'], s=50, c='r', marker='x', label='Rejected')
    ax.legend(loc=2)    # 画出图例
    ax.set_xlabel('Test 1 Score')
    ax.set_ylabel('Test 2 Score')


#定义一个特征映射函数：直到六次幂(power)
def feature_mapping(x1, x2, power):
    data = {}    #定义一个字典，从而在之后可以直接添加元素的名称和内容
    for i in np.arange(power + 1):    #i=n时，就是n次幂的所有组合。一共最多power次幂
        for p in np.arange(i + 1):
            data["f{}{}".format(i - p, p)] = np.power(x1, i - p) * np.power(x2, p)  #format是格式化函数。
    return pd.DataFrame(data)  #将字典转化为一种二维表，类似excel

# 定义代价函数：正则化的代价函数，解决过拟合现象
def cost(theta, X, y):  # 这是未正则化的代价函数（之前的那个）
    first = (-y) * np.log(sigmoid(X @ theta))  # 注意array中有三种乘法
    second = (1 - y) * np.log(1 - sigmoid(X @ theta))
    return np.mean(first - second)

def sigmoid(z):  # 之前那个代价函数要用到这个logist函数
    return 1 / (1 + np.exp(- z))

def costReg(theta, X, y, l=1):  # 这是正则化之后的代价函数 cost + (l/2m)(θ·θ)
    _theta = theta[1:]  # 不惩罚第一项θ0(常数)
    reg = (l / (2 * len(X))) * (_theta @ _theta)
    return cost(theta, X, y) + reg  # 即之前那个代价函数再加一部分

# 计算梯度,定义批量的梯度函数：同样是正则化后的，是那个导数
def gradient(theta, X, y):  # 这是之前未正则化的那个
    return (X.T @ (sigmoid(X @ theta) - y)) / len(X)

def gradientReg(theta, X, y, l=1):  # 这是正则化后的
    reg = (l / len(X)) * theta
    reg[0] = 0  # 同样不惩罚第一个θ
    return gradient(theta, X, y) + reg
    # 这个结果和网页上的那个结果顺序不一样，但是本质相同

def predict(theta, X):
    probability = sigmoid(X @ theta)
    return [1 if x >= 0.5 else 0 for x in probability]

data2 = pd.read_csv('ex2data2.txt', names=['Test 1', 'Test 2', 'Accepted'])
#初始化数据
plot_data()
plt.savefig("4.png")
# 注意到其中的正负两类数据并没有线性的决策界限。
x1 = data2['Test 1'].values   #定义x1，x2的数据
x2 = data2['Test 2'].values
_data2 = feature_mapping(x1, x2, power=6)  #得到特征映射后的二维表
# _data2.head()
x = _data2.values   #将二维表转化为array格式
y = data2['Accepted'].values  #初始化y
theta = np.zeros(x.shape[1])
# print(X.shape, y.shape, theta.shape)  #检查维度，期待输出((118, 28), (118,), (28,))

# 通过高级优化来进行梯度下降（计算速度很快，且不需要人为设定α）
# 只需传入cost函数，已经所求的变量theta，和梯度。cost函数定义变量时变量tehta要放在第一个，若cost函数只返回cost，则设置fprime=gradient。
result2 = opt.fmin_tnc(func=costReg, x0=theta, fprime=gradientReg, args=(x, y, 2))  # 预测出theta值
# 判断拟合的模型结果的精确性
final_theta = result2[0]
predictions = predict(final_theta, x)
correct = [1 if a==b else 0 for (a, b) in zip(predictions, y)]  #通过模型预测的结果与真实结果相比
accuracy = sum(correct) / len(correct)
# print(accuracy)  #输出0.8多
"""还可以使用sklearn来解决
from sklearn import linear_model#调用sklearn的线性回归包
model = linear_model.LogisticRegression(penalty='l2', C=1.0)
model.fit(X, y.ravel())"""

# 画出决策边界：不能像第一题那样ax.plot(x,y)那样画了，因为那样要列出y等于多少多少x,把y看成x的函数。这里幂太高式子复杂，很难变形。
# xx，yy是指二维图上的所有点（网格线）。
# 思路是：1.初始化网格点 2.所有网格点变成高次幂的形式、再×theta。这个式子命名为z，z等于0时，y=0.5。contour函数是默认z=0
x = np.linspace(-1, 1.5, 250)  # x维度是250 numpy.linspace(start, stop, num=50, endpoint=True, retstep=False, dtype=None)在指定的间隔范围内返回均匀间隔的数字
# 把xx,yy数据生成mesh网格状的数据，因为等高线的显示是在网格的基础上添加上高度值
xx, yy = np.meshgrid(x, x)  # 生成网格点坐标矩阵，xx维度是250,250  xx.ravel()的维度是62500，ravel即多维数组转化为一维数组,后面featuremapping的参数是要一维的。
z = feature_mapping(xx.ravel(), yy.ravel(), 6).values  # z的维度是(62500,28)
z = z @ final_theta  # z的维度是(62500,)  xx.shape = (250,250)
z = z.reshape(xx.shape)  # 将z的维度从一维的（62500，1）重新生成为(250,250)的矩阵

plot_data()  # 这个是一开始定义的那个函数，即执行后画出散点图.
# plt.contour([X, Y,] Z, [levels], ** kwargs) - 绘制轮廓（等高线）。
# X，Y ：值Z的坐标。X和Y必须都是2-D，且形状与Z相同，或者它们必须都是1-d，这样len（X）== M是Z中的列数，len（Y）== N是Z中的行数。
# Z ： 绘制轮廓的高度值。
# levels： int或类数组，确定轮廓线/区域的数量和位置。
plt.contour(xx, yy, z, 0)  # 这里z是一个关于xx、yy的高次幂的函数，默认z=0画图
plt.ylim(-.8, 1.2)  # 设置绘制的上下限
plt.savefig("5.png")
plt.show()  #最后显图，注意要所有元素都定义过后才能显图