连载的上一篇文章，我们实现了数据预处理与绘图的步骤。今天，小鱼将和大家一起实现逻辑回归的 4 个基础模块，为梯度下降策略做准备。

概率映射函数 sigmoid

Sigmoid 函数的公式为：

根据公式实现 sigmoid 方法：

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

Sigmode 函数的自变量取值为任意实数，值域为 [0,1] 。因此 Sigmode 函数可以将任意的值映射到 [0,1] 区间。我们在回归中得到的预测值，将其映射到 Sigmode 函数，即可完成由值到概率的转换，也就是分类任务。

使用上述实现的 sigmoid 函数绘图：

nums = np.linspace(-10, 10, 100)
fig,ax = plt.subplots(figsize=(12,4))
ax.plot(nums, sigmoid(nums))

Sigmoid 函数曲线如下：

计算预测值的模型 model

在线性回归中，我们使用特征参数和偏置项来计算预测值，其中，偏置项 θ₀ 可以视为 θ₀*1 ，即和一个恒为 1 的特征相乘。预测值的计算公式如下：

使用 model 函数实现上述公式：

def model(X, theta):
    """返回预测结果值"""
    return sigmoid(np.dot(X, theta.T))

下面，我们来验证一下 model 函数的实现是否正确。首先，在 DataFrame 中第 0 列的位置插入名为 Ones 的全为 1 的新列，在计算预测值时，该列将和偏置项相乘：

df.insert(0, 'Ones', 1)
df.head()

接下来，定义特征矩阵 X ，真实值 y，以及特征参数 theta 。其中特征 X 为数据集中的前 3 列组成的 ndarray 数组：

>> orig_data = df.values
>> cols = orig_data.shape[1]
>> X = orig_data[:,0:cols-1]
>> X[:5]
array([[ 1.        , -1.60224763,  0.63834112],
       [ 1.        , -1.82625564, -1.2075414 ],
       [ 1.        , -1.53903969,  0.3612943 ],
       [ 1.        , -0.28210129,  1.0863683 ],
       [ 1.        ,  0.69152826,  0.49337794]])

真实值 y 为数据集中的最后一列，但必须是 2 维数组的形式，因此第二个维度也使用切片的形式定位：

>> y = orig_data[:,cols-1:cols]
>> y[:5]
array([[0.],
       [0.],
       [0.],
       [1.],
       [1.]])

特征参数 theta 在初始化时只需要定义为全 0 矩阵即可。注意：theta 必须为 2 维矩阵，并且第一个维度为 1 ，第二个维度的长度与特征的个数相等，此处为 3：

>> theta = np.zeros([1,3])
>> theta
array([[0., 0., 0.]])

X y theta 的形状：

>> X.shape, theta.shape, y.shape
((100, 3), (1, 3), (100, 1))

这样， X 和 theta 的转置（3行1列）计算内积才可以得到 (100,1) 的预测值。

>> model(X, theta).shape
(100, 1)

计算平均损失 cost

在介绍逻辑回归时，我们推导出了如下的对数似然函数（负号表示梯度下降）：

求平均损失的公式就为：

下面，使用代码来实现上述公式：

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)

原始数据集使用默认参数的平均损失如下：

>> cost(X, y , theta)
0.6931471805599453

计算每个参数的梯度方向 gradient

梯度下降的关键步骤是在每个特征方向上求偏导，找到下山最快的方向。损失函数求偏导的公式如下：

根据上述公式，实现计算梯度的函数 gradient：

def gradient(X, y , theta):
    grad = np.zeros(theta.shape)
    error = (model(X, theta) - y).ravel()
    for i in range(len(theta.ravel())):
        term = np.multiply(error, X[:,i])
        grad[0, i] = np.sum(term) / len(X)
    return grad

上述求偏导的函数略显复杂，小鱼来为大家解释一下~

首先，因为我们要在每个特征方向上对 theta 求偏导，因此最终梯度的个数和特征参数的个数是一致的。故而使用 theta.shape 来初始化 gradient 数组：

grad = np.zeros(theta.shape)

在计算误差时，我们使用预测值 - 真实值，相当于将公式中的负号提到了括号里面。ravel 方法则可以将高位数组拉平，返回 1 维数组：

 error = (model(X, theta) - y).ravel()

最后，使用 for 循环计算每个特征参数的偏导数：

for i in range(len(theta.ravel())):
    term = np.multiply(error, X[:,i])
    grad[0, i] = np.sum(term) / len(X)

可以计算一下当前 theta 的梯度：

>> gradient(X, y, theta)
array([[-0.1       , -0.28122914, -0.25098615]])

有了梯度，下节我们就可以开始更新参数了：

即当前特征的参数 - 学习率 alpha * 当前特征参数对应的梯度。其中涉及到了学习率、每次迭代的样本个数以及梯度下降停止策略等，精彩内容下节奉上~