Logistic回归从头开始

您可以在我的GitHub上的实验室存储库中找到工作代码示例(包括这个示例)。

有时，为了预测新的、看不见的数据，有必要将现有数据分成几类。这个问题叫做分类，可以用来从数据中学习这些类的算法之一叫做Logistic回归。

在本文中，我们将深入研究Logistic回归模型，以了解它与其他回归模型(如线性回归或多元线性回归)的不同之处，如何从直观的角度考虑它，以及如何在从头开始实现它的同时将所学知识转换为代码。

如果你读过关于线性回归和多元线性回归的文章，你可能还记得我们的算法的主要目标是分别找到最合适的直线或超平面。

简单回顾一下，可以通过斜截形式表示一条线，如下所示：

在线性回归中，我们使用现有的数据找到一条斜率截距形式(a\(m\)和\(b\)组合)的直线，这条直线通过这样的数据拟合得最好。在线性回归中，我们使用了现有的数据来找到一条斜率截距形式的直线(a\(m\)和\(b\)组合)。

略微扩展坡度截取表单以支持多个\(x\)值和多个坡度(我们将使用\(\beta_n\)而不是\(m_n\))会产生以下结果：

这个放大的斜率截距公式被用在多元线性回归模型中，以求出与数据拟合最好的超平面的\(\β\)和\(b\)值。一旦找到，我们就可以通过插入\(x\)值来获得相应的\(y\)值来使用它进行预测。

线性回归模型总是将一组x值映射到连续尺度上得到的y值。这意味着\(y\)值可以是\(0\)、\(42\)或\(5.023.212\)。如果我们的\(y\)值是分类的，比如二进制值是\(0\)或\(1\)，我们将如何使用这样的回归模型？是否有办法定义阈值，以便将诸如\(42\)之类的值分配给类别\(1\)，而将诸如\(0.002\)之类的小值分配给类别\(0\)？

这就是Logistic回归发挥作用的地方。使用Logistic回归，我们可以将任何结果(Y)值映射到\(0\)和\(1\)之间的值，无论其大小如何。

让我们仔细看看我们需要做哪些修改，才能将线性回归模型转变为Logistic回归模型。

Logistic回归的核心是所谓的Sigmoid函数。Sigmoid函数是一类绘制时遵循S形的函数。

最突出的Sigmoid函数是所谓的Logistic函数，它是由Pierre Francois Verhulst发展起来的，用于模拟人口增长。它是通过以下公式进行数学描述的：

不要被数学吓倒！现在您需要知道的是，该函数接受任何\(x\)值，并将其映射到\(y\)值，范围从\(0\)到\(1\)。

绘制一系列\(x\)值的函数证明了这一说法，并产生了前述的S型曲线：

请注意，随着\(x\)值分别变小或变大，函数会越来越接近\(y\)值\(0\)或\(1\)。还要注意，\(x\)值\(0\)产生\(y\)值\(0.5\)。

这正是我们需要的。使用此函数，我们可以将任何数字(无论其大小)压缩为从\(0\)到\(1\)范围内的值。这使得函数结果是可预测的，这在我们稍后定义阈值以将函数输出与类相关联时非常有用。

注意：虽然有很多不同的Sigmoid函数可供选择，但很多人在谈论Logistic函数时都会使用Sigmoid函数这个名称。我们将遵循这一惯例，使用术语Sigmoid函数作为逻辑函数的同义词。

既然我们已经了解了Sigmoid函数的映射功能，我们应该能够在其中包装一个线性回归模型(如多重线性回归)，将回归原始输出转换为一个从\(0\)到\(1\)的值。

让我们把这个想法转化为数学吧。回想一下，我们的多元线性回归模型如下所示：

"；包装"；这在Sigmoid函数中(我们使用\(\sigma\)表示Sigmoid函数)导致以下结果：

我们需要做的第一件事是实现底层的多元线性回归模型。看一下数学，似乎可以使用点积来计算\(\beta\)和\(x\)部分，然后将单个\(b\)值相加。

为了使一切都更容易计算和实现，我们将使用一个小技巧。将值与标识符\(1\)相乘得到该值，因此我们将\(1\)添加到\(x\)值中，并将\(b\)添加到\(\beta\)值中。这样，我们可以只使用点积计算，而不需要稍后单独添加\(b\)。下面是这个把戏的数学公式：

\[\vec{x}=\Begin{pMatrix}1\\x_1\\.\\x_n\end{pMatrix}\vec{\beta}=\Begin{pMatrix}b\beta_1\\.\beta_n\end{pMatrix}\]

\[y=\vec{x}\cdot\vec{m}=\sum_{i=1}^n x_i\beta_i=x_1\Times\beta_1+.+x_n\Times\beta_n\]。

一旦我们计算了点积，我们需要将其传递给Sigmoid函数，以便将其结果转换为介于\(0\)和\(1\)之间的值。

def dot(a：List[Float]，b：List[Float])->；Float：assert len(A)==len(B)return sum([a_i*b_i for a_i，b_i in zip(a，b)])assert dot([1，2，3，4]，[5，6，7，8])==70。

这里的squish函数将\(x\)和\(\beta\)值作为参数(请记住，我们已将a\(1\)附加到\(x\)值，将\(b\)附加到\(\beta\)值)，使用dot函数计算\(x\)和\(\beta\)的点积，然后将此结果传递给Sigmoid函数，以将其映射到\(0\)和\(\beta\)之间的值。

def squish(beta：LIST[FLOAT]，x：LIST[FLOAT])->；FLOAT：Assert len(Beta)==len(X)#计算点积dot_result：Float=dot(beta，x)#使用sigmoid得到0到1之间的结果返回sigmoid(Dot_Result)assert squish([1，2，3，4]，[5，6，7，8])==1.0。

我们已经谈了很多关于Sigmoid函数是如何使函数结果可预测的解决方案，因为所有的值都映射到\(0\)-\(1\)范围。但是，该范围内的值代表什么呢？让我们来看一个例子。

下面是一个数据集，描述了学生为考试学习了多长时间，以及他们在给定学习时间的情况下是否通过了考试。

看一眼数据，似乎学生学习的时间越长，通过考试的可能性就越大。直觉上，这是有道理的。

看着绘制的数据，我们可以立即发现，这些值似乎贴在了图表的底部或顶部。鉴于使用线性回归模型来找到最能描述数据的直线似乎是不可行的。如果我们期望此行生成的值是\(o\)或\(1\)，那么如何通过数据拟合此行呢？

让我们来做一个思维实验。如果我们以某种方式找到了线性回归模型的一些系数，它最好地描述了数据，并通过Sigmoid函数传递它计算的结果，会发生什么呢？下面是添加了Sigmoid函数的上图：

查看上面的图表，我们可以看到Sigmoid函数确保来自底层线性回归模型的结果被映射到\(0\)和\(1\)之间的刻度上，这反过来又使得例如在\(0.5\)处定义阈值成为可能，即大于\(0.5\)的值可能是学生通过考试的预测因素，而小于\(0.5\)的值可能意味着她将考试不及格。

请注意，最后一句话的措辞不是巧合。西格诺德函数产生的值可以解释为概率，其中\(0\)表示\(0%\)概率，\(1\)表示\(100%\)概率。

事实证明，我们可以将上一节中的发现转换为一个称为概率密度函数或(简称PDF)的函数。

特别地，我们可以定义一个条件概率，它规定给定某个\(\β\)和\(x_i\)，每个对应的\(y_i\)应等于\(1\)概率\(\σ(\βx_i)\)和\(0\)概率\(1-\σ(\βx_i)\)：

看看上面的公式，我们是如何从上面的口头描述中推导出来的，这可能是一个谜。下面是我希望您尝试的内容：请通过将\(y_i\)设置为\(0\)来应用公式，然后将其设置为\(1\)，看看会发生什么情况。您将注意到的是，根据您设置的值，公式中只有一部分保持不变，而另一部分被取消。

通过Logistic回归，我们的主要目标是找到模型的β参数，使得对于一对x值，我们的模型计算的y值尽可能接近实际y值的可能性最大。

为了找到最佳参数，我们需要以某种方式计算我们的模型预测与当前设置的错误程度。

在上一节中，我们讨论了概率密度函数(PDF)，它似乎正好捕捉到了这一点。我们只需要稍微调整一下这个函数，这样我们以后就可以更容易地使用它进行计算。

我们将应用的主要调整是在\(\log\)函数中包装\(y_i=0\)和\(y_i=1\)的单独PDF计算。对数具有严格递增的良好特性，这使得以后对其数据进行计算变得更容易。此外，我们的PDF Main属性仍然保留，因为使预测正确(y\)的可能性最大化的任何\(\beta\)值集也会最大化\(\log\)可能性。

我们只需要解决一个小问题。总体而言，我们正试图最大限度地减少我们的模型产生的错误预测的数量，但再次查看对数图，我们看到函数正在严格增加。我们通过将对数与\(-1\)相乘来镜像\(x\)轴上的对数(将其倒置)。因此，我们的两个函数现在看起来如下所示：

现在，我们要做的最后一件事是将两个函数放回到一个等式中，就像我们上面对复合PDF函数所做的那样。执行此操作将导致以下结果：

\[\log L(\beta\midx_i y_i)=-(y_i\log(\sigma(\beta x_i))+(1-y_i)\log(1-\sigma(\beta x_i)\]。

此函数称为对数损失(或对数损失/对数似然率)，稍后我们将使用它来确定我们的模型与其预测结果之间的关系。同样，您可能希望将\(y_i\)设置为\(0\)或\(1\)，以查看等式的一部分已被抵消。

def neg_log_lisilience(y：Float，y_pred：Float)->；FLOAT：return-((y*log(Y_Pred))+((1-y)*log(1-y_pred)assert 2.30<；neg_log_lisilience(1，0.1)<；2.31assert 2.30<；neg_log_lisilience(0，0.9)<；2.31assert 0.10<；0.11assert 0.10<；neg_log_lisilience(0，0.1)<；0.11

接下来，让我们使用我们的编码版Log Lost创建\(y=0\)和\(y=1\)的曲线图：

xs_nll：LIST[FLOAT]=[x/10000 for x in range(1,10000)]图，(ax1，ax2)=plt.subplots(1，2)ax1.lot(xs_nll，[neg_log_lisilience(1，x)for x in xs_nll])ax1.set_title('；y=1'；)ax2.lot(xs_nll，[neg_log_lisilience(0，x)for x in xs_nll])ax1.set_title('；y=1'；)ax2.lot(xs_nll，[neg_log_lisilience(0，x)for x in x in。)；

正如我们所看到的，预测的错误越大，计算的误差就越大。因此，日志丢失功能惩罚的是错误行为，而不是奖励正确的行为。

为了计算整个数据集的总体误差，我们对每个单独的原木损失计算求和并求平均值：

定义错误(ys：List[Float]，ys_pred：List[Float])->；Float：assert len(Ys)==len(Ys_Pred)Num_Items：int=len(Ys)sum_nll：Float=sum([neg_log_lisilience(y，y_pred)for y，y_pred in zip(ys，ys_pred)])return(1/num_item)*sum_nllassert 2.30<；error。错误([0]，[0.9])<；2.31assert 0.10<；Error([1]，[0.9])<；0.11assert 0.10<；Error([0]，[0.1])<；0.11。

现在，我们需要实现的最后一个缺失部分是优化步骤。我们最终想要的是一个带有β参数的Logistic回归模型，这些参数与x值相结合，对任何y值都能产生最准确的预测。

既然我们现在可以计算当前模型及其\(\beta\)参数产生的误差，我们就可以迭代地更改\(\beta\)参数，直到我们的模型无法再改进(无法降低误差值)。

注意：我已经写了一篇专门的帖子，非常详细地解释了这个算法，所以我不会在这里讲太多细节，如果你不熟悉梯度下降，我会鼓励你阅读这篇文章。

它的要点是，给定我们的Log Lost函数，我们可以找到Log Lost函数计算的误差最小的一组\(\beta\)参数。因此，如果错误不能再减少，我们已经找到了局部(或全局)最小值。为了找出最小值位于何处，我们将使用误差函数梯度，它是一个向量，并引导我们到达该位置。

因为我们一次又一次地重复这样的计算，所以我们反复向下迭代地向下递减误差函数表面，因此我们将其命名为梯度下降。

梯度的计算是通过计算对数损失函数相对于我们的(Xi)矢量中的各个(x_{ij})值的偏导数来完成的。

grad：List[Float]=[0 for_in range(len(Beta))]，用于zip(xs，ys)中的x，y：Err：Float=squish(beta，x)-y for i，x_i in Enumerate(X)：grad[i]+=(err*x_i)grad=[1/len(X)*g_i for g_i in grad](X)：grad[i]+=(err*x_i)grad=[1/len(X)*g_i for g_i in grad]。

再说一次，如果你不熟悉梯度和偏导数，你可能会想阅读我关于梯度下降算法的文章。

我们终于把所有的东西都准备好了！让我们抓取一些数据，并使用Logistic回归模型对其进行分类！

我们将要使用的数据集与我们在上面的例子中已经看到的类似，在上面的例子中，我们试图根据学生学习的时间来预测他是否会通过考试。

您可以在此处下载数据集。检查数据，我们可以看到有2个浮点值和一个整数值，表示标签和(可能)表示问题学生是否通过考试的整数值。

我们的任务是使用此数据集训练Logistic回归模型，该模型将帮助我们分配标签\(0\)或\(1\)(即，分类)新的、不可见的数据点。

接下来，我们需要解析文件并提取\(x\)和\(y\)值：

MARKS_DATA_PATH：path=data_dir/'；marks.txt'；xs：LIST[LIST[FLOAT]]=[]ys：LIST[FLOAT]=[]WITH OPEN(MARKS_DATA_PATH)AS FILE：对于文件中的行：DATA_POINT：LIST[STR]=line.strie().Split('；，'；)x1：FLOAT=FLOAT(DATA_POINT[0])x2：FLOAT=FLOAT(DATA_POINT[1])y：INT=INT(DATA_POINT[2])xs.append([x1，x2])ys.append(Y)。

绘制数据图表，看看我们是否需要处理任何异常值或其他意外情况，这总是一个好主意：

x1s：List[Float]=[x[0]for x in xs]x2s：List[Float]=[x[1]For x in xs]plt.dister(x1s，x2s，c=ys)plt.axt([min(X1s)，max(X1s)，min(X2s)，max(X2s)])；

看起来我们这里(几乎)不错。我们只有一个方面需要进一步检查。如果您正在查看轴，您会看到值的范围从\(\约35\)到\(\约95\)。使用如此大的值可能会导致稍后使用\(\log\)和\(\exp\)时出现问题。我们将立即解决这个问题。

您可能还记得在这篇文章的开头，我们应用了一个技巧，我们在每个\(x\)向量前面加上\(1\)，并在\(\beta\)向量前面加上\(b\)，以便可以使用点积(这是一种更简单的计算)。

下面的代码在每个\(x\)向量前面加上一个\(1\)，这样我们以后就可以利用该计算技巧：

对于xs中的x：x.insert(0，1)xs[：5]#[1，34.62365962451697，78.0246928153624]，#[1，30.28671076822607，43.89499752400101]，#[1，35.84740876993872，72.90219802708364]，#[1，60.18259938620976，86.30855209546826]，#[1，79.0327360507101，75.3443764369103]]。

下一步，我们将解决我们在上面的段落中提到的伸缩性问题。为了解决这个问题，我们要做的是标准化(通常通过z-Score)整个数据集：

xs=z_SCORE(Xs)xs[：5]#[1，-1.5942162646576388，0.6351413941754435]，#[1，-1.8171014180340745，-1.2014885239142388]，#[1，-1.531325157335502，0.3594832875590465]，#[1，-0.28068723821760927，1.0809228071415948]，#[1，0.6880619310375534，0.4909048515228952]]。

如果您对z_core函数的作用感到好奇，请查看我在GitHub上的实验室存储库中的整个实现。

现在我们终于在一个位置上，我们可以训练我们的Logistic回归模型通过梯度下降。我们将从对模型参数的一些随机猜测开始，通过使用误差函数计算模型当前的总体错误行为，然后使用误差函数梯度来更新模型的值，以便在下一次迭代中产生较小的误差，从而迭代优化我们的模型。最终，我们将收敛到局部最小值，这将导致\(\beta\)值计算出最小的误差。

beta：list[Float]=[Random()/10 for_in range(3)]print(f'；以"；beta"；：{beta}'；开始)POECHS：INT=5000 Learning_rate：FLOAT=0.01对于范围内的PORCH(POECHS)：#计算"；预测"；(`beta`和`x`的扁平化点积)基于我们当前的`beta`向量ys_pred：list[Float]=[squish(beta，x)for x in xs]#计算并打印错误if pech%1000==True：Lost：Float=error(ys，ys_pred)print(f'；Epoch{Epoch}-->；Loss：{Loss}'；)#为zip(xs，ys)中的x，y计算梯度grad：list[浮点]=[0 for_in range(len(Beta))]：err：浮点=i的squish(beta，x)-y，枚举中的xi(X)：grad[i]+=(err*xi_i)grad=[1/len(X)*g_i for g_i in grad]#朝最大的减少方向迈出一小步β=[b+。Grad)]打印(f'；#34；Beta"；：{Beta}'；)#从"；Beta"；开始：[0.06879018957747185，0.060750489548129484，0.08122488791609535]#纪元1--&>损失：0.6091560801945126#纪元1001--&>损失：0.2037432848849053#2001年纪元--&>损失：0.20350230881468107#纪元3001-&>损失：0.20349779972872906#纪元4001-&>损失：0.20349770371660023#对&#的最佳估计。：[1.7184091311489376、4.01281584290694、3.7438191715393083]。

既然我们的模型已经过训练，我们可以计算一些统计数据来看看它有多精确。对于这些计算，我们将把阈值设置为\(0.5\)，这意味着我们的模型产生的每个高于\(0.5\)的值都被认为是\(1\)，每个小于\(0.5\)的值都被认为是\(0\)。

TOTAL：INT=len(Ys)Thresh：Float=0.5true_Positions：Int=0true_Positions：Int=0false_Positions：Int=0false_Negative：int=0 for i，x in Enumererate(Xs)：y：int=ys[i]pred：Float=Squish(beta，x)y_pred：int=1 if pred<；阈值：y==1且y_pred==1时y_pred=0：TRUE_PRODIES+=1 elif y==0且y_pred==0：TRUE_NECTIONS+=1 ELIF y==1且y_pred==0：FALSE_NECTIONS+=1 ELIF y==0且y_pred==1：FALSE_PORTIONS+=1print(f'；TRUE PRECTIONS：{TRUE_PRECTIONS}'；)

..