知也无涯

吾生也有涯,而知也无涯,能学一点算一点…

  • 在2015年,由 OpenAI 的 DP Kingma 等发布了 《ADAM: A METHOD FOR STOCHASTIC OPTIMIZATION》算法后,由于其迭代效率提升非常明显,所以 ADAM(或其变种)就被广泛的采用。本文将继续对上一篇介绍的梯度下降算法进行优化,并介绍 ADAM 算法(一种对随机梯度下降算法的优化算法)的实现以及效果。

    Stochastic Gradient Descent 或者说 mini-batch解决了样本量巨大时,梯度下降迭代的问题。但是,也带了一些新的问题。最为主要的是,因为样本数据的波动,而导致每次梯度下降计算时,梯度方向的波动,从而降低了梯度下降迭代的效率。

    在前面的《Mini-batch Gradient Descent和随机梯度下降(SGD)》文章中,我们对比了 mini-batch 和 batch gradient descent 的在迭代时,目标函数下降的速度。

    可以看到,batch gradient descent 的目标函数下降非常稳定,而 Mini-batch 的实现则会有明显的波动。为了尝试修正这个问题,从而提高迭代效率,在神经网络算法上,逐渐探索出了一些较为高效的优化算法:Adam SGD。该算法将 RMSprop 和 “Exponential smoothing”的想法结合在一起,形成了一个较为高效的算法,在实践中被广为使用。


    Stochastic Gradient Descent 与 Momentum

    SGD 会在每次迭代时根据样本的偏差,展现出不同的偏差,所以,在使用SGD进行迭代时,观察其 cost函数下降,应该会有更加明显的波动(后续吧自己实现的程序改造后,尝试观察一下)。

    为了加快迭代的速度,一个折中的思路是,引入一个均值替换当前的梯度方向。该如何引入这个均值呢?梯度是一个随时计算推进,不断推进的变量,常用的均值计算可以参考:Moving average。最为常见的实现是使用“Exponential moving average”,这种平均值的计算,在迭代计算时实现非常简单。

    Momentum 就是 “Exponential moving average”实现时的参数“smoothing factor”,在神经网络中,经常使用 \( \beta \)表示(原因是 \( \alpha \) 已经表示学习率了 )。

    而这里的 Momentum ,也是 TensorFlow 在构造 SGD 算法时需要的另一个参数。

    关于Exponential moving average

    或者叫“Exponential smoothing”。我们看看这个算法的具体实现是怎样的?

    原始的迭代:\( w = w – \alpha \frac{\partial J}{\partial w} \)

    使用 “Exponential smoothing” 后的迭代:

    $$
    \begin{align}
    v_0 & = 0 \quad \partial{w}_t = \frac{\partial J}{\partial w}|_{(for \, sample \, t)} \\
    v_{t} & = \beta*v_{t-1} + (1-\beta)\partial{w}_{t} \\
    w & := w – \alpha v_t
    \end{align}
    $$

    考虑 \( \beta = 0.9 \),如果数学直觉比较好的话,可以看出,原本使用梯度\( \partial{w} \)进行迭代的,这里使用了一个梯度的“Exponential smoothing” \( v_t \)去替代。上面的式子中,\( v_t \) 如果展开有如下表达式:

    $$
    \begin{align}
    v_t & = (1-\beta)\partial{w}_{t} + \beta(1-\beta)\partial{w}_{t-1} + \beta^2(1-\beta)\partial{w}_{t-2} … \\
    & = \sum\limits_{i=0}^{t} \beta^{i}(1-\beta)\partial{w}_{i}
    \end{align}
    $$

    使用“Exponential smoothing” 之后,新的迭代方向 \( v_t \),可以理解为一个前面所有梯度方向的加权平均。离得越近的梯度,权重越高,例如,\( \partial{w}_{t} \)的权重是\( (1-\beta) \);而之前的梯度,则每次乘以一个 \( \beta \)衰减。

    Exponential moving average的“冷启动问题”与修正

    仔细观测上诉的 “Exponential moving average” 公式,可以注意到一个问题,就是其最初的几个点总是会偏小。其原因是,当前值的权重总是为 \( 1- \beta \),而因为是初始的几个值,并没有更前面的数据去“平均”当前值,也就会出现,初始值总是会偏小的问题。

    通常,如果样本量很大的事时候,则可以忽略这个问题,因为初始值偏小的点占比会非常少,可以忽略。如果要一定程度上解决这个问题,也有继续对上述的 “Exponential moving average”做了一些修正,可以考虑对 \( v_t \)的结果值做一个修正:\( v_t := \frac{vt}{1-\beta^t} \)。

    一般的,因为样本的数量总是比较大的,所以我们可以忽略这个问题,而无需做任何修正。

    RMSprop

    在前面的“Gradient Descent with Momentum”中,我们看到为了解决梯度波动较大的问题,使用了 “Exponential moving average” 去尝试将一些比较偏的梯度,拉倒一个较为平均的方向上来。RMSprop的想法也是类似的,这里通过了root mean square的想法进行平均值的计算。具体的,在进行 SGD 时,每次更新梯度,按照如下的方法进行更新:

    $$
    \begin{align}
    s_0 & = 0 \quad \partial{w}_t = \frac{\partial J}{\partial w}|_{(for \, sample \, t)} \\
    s_{t} & = \beta*s_{t-1} + (1-\beta)(\partial{w}_{t})^2 \\
    w & := w – \alpha \frac{\partial w}{\sqrt{s_{t}}}
    \end{align}
    $$

    说明:这里对梯度进行平方时,如果在程序中是一个梯度向量,那么这里“平方”也就是对梯度的每一个分量进行一次平方。

    在“Exponential smoothing”的实现中,是将当前值,使用一个加权平均替代。与“Exponential smoothing”类似的,原本的梯度方向,现在使用如下的方向去替代了:

    $$
    \begin{align}
    s_t & = \frac{\partial{w}_{t}}{\sqrt{(1-\beta)(\partial{w}_{t})^2 + \beta(1-\beta)(\partial{w}_{t-1})^2 + \beta^2(1-\beta)(\partial{w}_{t-2})^2 + \cdots }} \\
    & = \frac{\partial{w}_{t}}{\sqrt{\sum\limits_{i=1}^{t}\beta^i(1-\beta)(\partial{w}_{i})^2}} \\
    \end{align}
    $$

    Adam Gradient Descent

    这可能是实际使用最多的算法,全称是 Adaptive Moment Estimation 。该实现,将 “Momentum” 和 “RMSprop” 做了一定的融合,形成了新的“最佳实践” Adam。在融合上,具体的实现与两个细节点:

    (1) 在 Adam 中均使用了“修正”计算,即 \( \hat{v_t} = \frac{v_t}{1-(\beta_1)^t} \quad \hat{s_t} = \frac{s_t}{1-(\beta_1)^t} \)

    (2) 参数更新公式,使用了两个算法的融合: \( w := w – \alpha \frac{\hat{v_t}}{\sqrt{\hat{s_t}}} \)

    Adam optimization的效果对比

    在 Adam 的论文中对于效果做了非常多的评估,感兴趣的可以参考相关论文。

    这里根据之前完成的训练程序,也进行了优化,实现了Adam算法。在 MNIST 数据集的训练上,我们来看看 Adam 的效果:

    从右图可以看到,Adam(蓝色)明显的提升了迭代效率。依旧一定程度存在 mini-batch(绿色) 的梯度波动的问题。相比于,batch gradient descent (红色)算法,迭代效率大大增加,约在第10次迭代,即在第一个epoch 的第十批样本进行训练时,cost 就下降到了比较低的程度。

    关于 root mean square

    root mean square也叫二次平均值,考虑一组数据:\( {x_1,x_2, \cdots , x_n } \),其RMS则为:

    $$ x_{rms} = \sqrt{\frac{1}{n} \sum_{i=1}^n x_i^2} = \sqrt{\frac{1}{n} (x_1^2 + x_2^2 + \cdots + x_n^2)} $$

    补充说明

    可以看到,所有的这些优化都是面向“最优化”问题的。梯度下降是一个一阶优化(First-order Optimization)的方法,其核心就在与每次迭代时,应该如何去更新响应的参数值,在梯度下降中也就是如何去选择合适的学习率。

    牛顿法是典型的二阶优化(Second-order Optimization),在迭代时使用了二阶导数,所以,通常可以获得更好的迭代效率。但是因为二阶导数的计算复杂度会上升非常多(对应的矩阵可能是所有参数的平方,应该也有人尝试去算过了…)。这也是为什么在这个场景下,依旧是使用一阶优化方法的原因。

    如果想比较好的理解学习率、Momentum、RMSprop、Adam等内容,建议先了解梯度、数值方法、最优化问题等数学方法。

    到这里这个系列算是一个小阶段了,这是一个个人学习的笔记,从数学的梯度概念开始,逐步到神经网络训练的Adam优化算法,也包含部分动手实践的神经网络算法实现。完成的系列包括了:

  • Terraform 可以自动化的创建云端的资源,但是要想实现更高的灵活度,则需要更为灵活的使用Terraform的“Data Sources”能力。例如,在自动化的创建数据库时,云厂商允许创建的版本号是在动态变化的,例如,当前最新的允许的创建的MySQL版本通常是 8.0.40,但通常过了一个季度之后,就变成了 8.0.41。这时,对应的 Terraform 的脚本就需要调整或者传递参数就需要发生变化。而 Terraform 提供的 “Data Sources” 能力则可以很好的解决这个问题。

    在 Oracle 的 Terraform 中可以使用 “Data Source: oci_mysql_mysql_versions” 实现该能力。

    示例

    首先使用 data 命令定义该对象:

    data "oci_mysql_mysql_versions" "gmv" {
        compartment_id = oci_identity_compartment.oic.id
    }

    这里会获取该租户环境下支持的所有MySQL版本。

    然后,再使用 output 命令就可以获取并输出这些版本信息。详细的output命令如下:

    output "mysql_version" {
      value       = data.oci_mysql_mysql_versions.gmv.versions
    }

    详细的输出示例如下:

    mysql_version = tolist([
      {
        "version_family" = "8.0"
        "versions" = tolist([
          {
            "description" = "8.0.36"
            "version" = "8.0.36"
          },
          {
            "description" = "8.0.37"
            "version" = "8.0.37"
          },
          {
            "description" = "8.0.38"
            "version" = "8.0.38"
          },
          {
            "description" = "8.0.39"
            "version" = "8.0.39"
          },
          {
            "description" = "8.0.40"
            "version" = "8.0.40"
          },
          {
            "description" = "8.0.41"
            "version" = "8.0.41"
          },
        ])
      },
      {
        "version_family" = "8.4 - LTS"
        "versions" = tolist([
          {
            "description" = "8.4.0"
            "version" = "8.4.0"
          },
          {
            "description" = "8.4.1"
            "version" = "8.4.1"
          },
          {
            "description" = "8.4.2"
            "version" = "8.4.2"
          },
          {
            "description" = "8.4.3"
            "version" = "8.4.3"
          },
          {
            "description" = "8.4.4"
            "version" = "8.4.4"
          },
        ])
      },
      {
        "version_family" = "9 - Innovation"
        "versions" = tolist([
          {
            "description" = "9.1.0"
            "version" = "9.1.0"
          },
          {
            "description" = "9.1.1"
            "version" = "9.1.1"
          },
          {
            "description" = "9.1.2"
            "version" = "9.1.2"
          },
          {
            "description" = "9.2.0"
            "version" = "9.2.0"
          },
        ])
      },
    ])

    获取特定大版本的各小版本

    可以通过 data资源中新增filter模块以过滤出需要的对象。

    在 Terraform 中,关于 data 资源是否可以使用 filter,以及filter支持的完整度视乎并没有明确的说明。这需要更具不同的供应商的实现。常见的,在data resourcefilter可以支持“列表匹配”、“通配符匹配”或者“正则匹配”。具体的匹配方式,则需要通过文档、或者测试区验证。

    添加带正则匹配的 filter
    data "oci_mysql_mysql_versions" "gmv" {
        #Required
        compartment_id = oci_identity_compartment.oic.id
        filter {
            name = "version_family"
            values = ["8.0.*"]
            regex  = true
        }
    }

    通过 HCL 语言获取最新的版本
    output "latest_versions" {
      value = {
        for db_version in data.oci_mysql_mysql_versions.gmv.versions : db_version.version_family => sort([
          for v in db_version.versions : v.version
        ])[length(db_version.versions) - 1] // 取排序后的最后一个版本
      }
    }

    最后的输出如下:

    latest_versions = {
      "8.0" = "8.0.41"
      "8.4 - LTS" = "8.4.4"
      "9 - Innovation" = "9.2.0"
    }

    参考链接

  • 快速了解 Aurora DSQL

    ·

    上周在 AWS re:Invent大会(类似于阿里云的云栖大会)上推出了新的产品 Aurora DSQL[1] ,在数据库层面提供了多区域、多点一致性写入的能力,兼容 PostgreSQL。并声称,在多语句跨区域的场景下,延迟只有Google Spanner的1/4。

    Aurora DSQL 提供了多可用区、多区域的多点一致性写入的内容。在技术层面,Aurora DSQL 通过把数据库的 log 模块和 block (或者说是cache)模块做了分离,从而更好的实现多点/多区域分布式能力,这与 Google AlloyDB 是比较类似的;此外,在跨区域强一致性实现上,则使用“Amazon Time Sync Service” [3] 来保障多个区域之间事务顺序的一致性。

    在产品层面,分为两个场景,一个是 Aurora DSQL(region内模式)和一个 Aurora DSQL Global 模式(多 region 内模式)。在 Region 内场景下,相比于普通 Aurora PostgreSQL ,Aurora DSQL 在多个可用区内都可以提供强一致的读写接入点,而Aurora PostgreSQL只在一个可用区提供写,其他可用区仅提供只读节点。

    在跨 Region 的场景下,Aurora DSQL 则提供了同步的、跨区域的多点写入能力。这对于业务在全球分布的客户,则可以进一步的降低业务的复杂度。而原来的 Aurora Global Database 仅提供单个 Region 的写入能力,并且,在其他 Region 的读节点需要承受一定的数据访问延迟,这对于很多的在线业务场景可能是无法接受的,或者需要在应用层面做针对性的改造。

    这是 Aurora 发布的10周年,AWS 依旧是创新、技术能力非常强的一家公司。此外,产品是在内测阶段,普通用户还无法体验。

    参考文档

  • 在前述的文章(参考)中,我们实现了带有一个隐藏层的神经网络,并使用该神经网络对手写数字0/1进行识别。本文对该神经网络的识别效果以及相关的超参数的配置做一些分析与优化。

    这里涉及的超参数包括了学习率、迭代次数、隐藏层神经元的个数,这里对这三个参数的不同取值进行了相关测试,并观察训练时间与模型效果。

    不同学习率的模型训练

    学习率应该是这里最为重要参数了。在相同的迭代次数下(这里取500),不同的学习率展现出了非常大的差异。这里从0.001开始、尝试了:0.001、0.005、0.01、0.1、0.5等取值。详细的数据如下:

    可以看到,不同的学习率展现出了训练效率的差异非常大:

    • 在相同的迭代次数(均取500)情况下,学习率增加到0.1之后,预测错误率降低到了0.09%,并且再增加学习率,预测错误率并没有提升
    • 在学习率,从0.001增加到了最后的0.5之后,在进行了相同的迭代次数时,训练的目标函数取值下降一直都较为明显

    学习率如何影响目标函数的收敛速度

    右图展示了学习率取值分别为0.1和0.01时,目标函数的收敛速度趋势图。可以看到:

    • 学习率为 0.1 时,在迭代约40次以前,目的函数的收敛速度非常快,并快速的收敛到了非常低的水平
    • 学习率为0.01时,迭代到100次时,代价依旧非常高

    从这次实现代码也可以看到,学习率对于模型的训练效率有这至关重要的影响。如果学习率选择不合适,则会耗费大量计算资源进行非常慢的训练。那么,如果选择合适的学习率以进行更加高效进行梯度下降迭代,这是一个比较复杂的问题,这里暂时先挖个小坑在这里,待后续再做更多讨论。

    迭代次数 epoch 如何影响模型

    这里选取学习率为0.01,隐藏层10个人工神经元,从而观测随着“迭代次数”效率如何影响:

    可以看到,当迭代不够充分时,目标函数收敛还不够时,模型效果也会比较差。随着迭代次数不断增加,目标函数下降就不再明显了。完整的目标函数收敛趋势如下图:

    隐藏层神经元个数与模型效果

    这里观察隐藏层神经元个数与模型效果趋势图。这里分别测试了1、10、50、100、150、300个神经元时候模型的表现,如下图:

    从测试来看,在这个案例中,随着隐藏层神经元个数的增加并不会提升模型性能的。这可能暗示了,此类任务(图像识别相关)使用前馈神经网络时,其性能可能较差。

    部分识别失败的图片

    在该模型与训练下,部分识别失败率比较高的图片如下:

    9879
    8325
    9634
    3073
    2185

  • 更进一步,本文将从零构建一个浅层的前馈神经网络,实现手写数字0和1的识别。示例手写数字如下图所示。整体实现包括数据预处理、参数初始化、forward/backword propagation、训练与预测。

    问题描述

    使用 Python ,但不使用任何深度学习框架,编写一个浅层的神经网络实现,识别 MNIST 数据集中的手写数字 0 和 1。


    前置知识

    • 熟悉Python 、NumPy的使用
    • 了解神经网络的表示与forward propagation
    • 了解 backward propagation
    • 了解 MNIST 数据集

    数据说明

    这里的手写数字使用了数据集 MNIST,该数据集有着机器学习界的“果蝇”之称。里面包含了大量经过预处理的手写数字,这里选择其中标记为0或1的图片进行训练。并使用对应的测试集,验证训练的神经网络的效果。

    这里取一张 MNIST 中的图片,以及对应的像素数据,以帮助直观理解其中的数据。MNIST 中的图片数据可以理解为一个如下的 28*28 的数组,其所代表的图片也如下图所示:

    更为具体的,数组中的每一个数值代表了图片中的某个像素的灰度值,0代表黑色,255代表白色,这之间的值则代表介于之间的灰色。

    MNIST 原始数据可以在:“THE MNIST DATABASE of handwritten digits”页面获取。也可以直接使用 Python 中模块keras.datasets获取。本程序使用后者方式获取。

    实现概述

    这里的使用了两层神经网络,即一个隐藏层,一个输出层,因为这是二分类问题,所以输出层使用logistic函数,而隐藏层则使用了 ReLU 作为激活函数。隐藏层的神经元是一个动态参数(超参数),在测试时,可以调整为10~300之间不等。输入层则是,MNIST 的图片,也就是一个28*28个数组。

    这里使用了最为基础的Gradient Descent算法,有时候也被称为Batch Gradient Descent。

    数据预处理

    在 Python 中可以使用 keras.datasets模块便捷的获取 MNIST 的数据集。这里编写了一个简单的函数,将该数据集中0和1相关的图片筛选出来,然后将数组的维度,转化为需要的维度。例如,期望的输入数组的维度为\( (n^{[0]},m) \),其中 \( n^{[0]} \) 表示输入的特性(feature)数量,这里也就是 784(即28*28),\( m \)表示样本个数。

    具体代码如下:

    # return only data lable 0 or 1 from MNIST for the binary classification
    def filter_mnist_data(data_X, data_y):
        data_filter = np.where((data_y == 0) | (data_y == 1))
        filtered_data_X, filtered_data_y = data_X[data_filter], data_y[data_filter]
        r_data_X = filtered_data_X.reshape(filtered_data_X.shape[0],filtered_data_X.shape[1]*filtered_data_X.shape[2])
        return (r_data_X, filtered_data_y)
    
    (train_all_X, train_all_y), (test_all_X, test_all_y) = mnist.load_data()
    (train_X,train_y) = filter_mnist_data(train_all_X, train_all_y)
    (test_X ,test_y ) = filter_mnist_data(test_all_X, test_all_y)
    
    X  = train_X.T
    Y  = train_y.reshape(1,train_y.shape[0])

    超参数的配置

    该程序涉及的超参数(hyper-parameter)包括:隐藏层神经元个数 \( n_1 = 10 \)、学习率 \( \alpha = 0.5 \)、迭代次数等。

    # hyper-parameter; read the comments above for structure of the NN
    n0 = X.shape[0]   # number of input features
    n1 = 10           # nerons of the hidden layer
    n2 = 1            # nerons of the output layer
    iteration_count = 500
    learning_rate   = 0.5

    feature scaling和参数初始化

    为了增加训练的速度,这里对输入进行了“标准化”,将所有的数据,归一化为均值为0、方差为1的数据集。具体的步骤如下:

    # feature scaling / Normalization
    mean = np.mean(X,axis = 1,keepdims = True)
    std  = np.std( X,axis = 1,keepdims = True)+0.000000001
    X  = (X-mean)/std

    注意:实现过程中需要注意的是,如果对输入进行了标准化,那么在预测时也需要对预测的输入进行相同的归一化。所以这里的均值和方差,需要记录下来,以备后续使用。这里给方差额外加了一个非常小的数字,一般是没有必要的,这里是为了防止输入数据全部都相同,即方差为0。

    此外,为了增加训练的速度,也参考业界通用做法,对输入也进行了随机初始化,具体如下:

    # initial parameters: W1 W2 b1 b2 size
    np.random.seed(561)
    W1 = np.random.randn(n1,n0)*0.01
    W2 = np.random.randn(n2,n1)*0.01
    b1 = np.zeros([n1,1])
    b2 = np.zeros([n2,1])

    Forward Propagation / Backward Propagation

    Forward Propagation 总是简单的。Backward Propagation 本身的计算与推导是非常复杂的,这里不打算详述(后续再单独介绍)。这里把 Backward Propagation 的计算结果列举如下。对于输出层:

    $$
    \begin{align}
    dW^{[l]} & = \frac{\partial J}{\partial W^{[l]}} \\
    & = \frac{\partial J}{\partial Z^{[l]}} @ (A^{[l-1]})^T
    \\
    dZ^{[l]} & = \frac{\partial J}{\partial Z^{[l]}} \\
    & = \frac{1}{m}(A^{[l]} – Y)
    \end{align}
    $$

    这里的 @ 表示矩阵乘法符号。 上标 T表示矩阵转置。这里一共两层,故这里的 \( l = 2 \) 。

    对于隐藏层:

    $$
    \begin{align*}
    dW^{[k-1]} & = \frac{\partial J}{\partial W^{[k-1]}} \\
    & = \frac{\partial J}{\partial Z^{[k-1]}} @ (A^{[l-2]})^T &
    \\
    dZ^{[k-1]} & = \frac{\partial J}{\partial Z^{[k-1]}} \\
    & = (W^{[k]})^T @ \partial Z^{[k]} \cdot g\prime(Z^{[k-1]}) &
    \end{align*}
    $$

    因为这里只有一个隐藏层和一个输出层,故这里的 \( k = 2 \);相同的,这里的 @ 表示矩阵乘法符号;这里的 \( \cdot \) 表示对应元素相乘(element-wise);函数 \( g \) 表示 ReLU函数。

    具体的代码实现:

        # forward propagation
        A0 = X
    
        Z1 = W1@X + b1  # W1 (n1,n0)  X: (n0,m)
        A1 = np.maximum(Z1,0) # relu
        Z2 = W2@A1 + b2
        A2 = logistic_function(Z2)
    
        dZ2 = (A2-Y)/m
        dW2 = dZ2@A1.T
        db2 = np.sum(dZ2,axis=1,keepdims = True)
    
        dZ1 = W2.T@dZ2*(np.where(Z1 > 0, 1, 0)) # np.where(Z1 > 0, 1, 0) is derivative of relu function
        dW1 = dZ1@A0.T
        db1 = np.sum(dZ1,axis=1,keepdims = True)

    训练迭代

    每次训练迭代都需要根据样本数据,进行一次上述的 Forward Propagation / Backward Propagation ,然后更新新的参数值,以用于下一次迭代,具体的实现如下:

    cost_last = np.inf # very large data,maybe better , what about np.inf
    for i in range(iteration_count):
        ...
        Forward Propagation / Backward Propagation
        ...
        W1 = W1 - learning_rate*dW1
        W2 = W2 - learning_rate*dW2
        b1 = b1 - learning_rate*db1
        b2 = b2 - learning_rate*db2

    使用测试数据集验证训练效果

    在完成训练后,就可以对测试集中的数据进行验证了。需要注意的是,在前面做了“feature scaling”,这里需要对应将输入数据做对应的处理。将预测结果,与标注数据进行对比,就可以确定该模型的效果了。

    # Normalization for test dataset
    X  = (test_X.T - mean)/std
    Y  = test_y.reshape(1,test_y.shape[0])
    
    Y_predict = (logistic_function(W2@np.maximum((W1@X+b1),0)+b2) > 0.5).astype(int)
    
    for index in (np.where(Y != Y_predict)[1]):
        print(f"failed to recognize: {index}")
        # np.set_printoptions(threshold=np.inf)
        # np.set_printoptions(linewidth=np.inf)
        # print(test_X[index].reshape(28,28))
    
    print("total test set:" + str(Y.shape[1]) + ",and err rate:"+str((np.sum(np.square(Y-Y_predict)))/Y.shape[1]))

    完整的代码

    完整的代码可以参考 GitHub 仓库中的 ssnn_bear.py 脚本:参考。这里张贴当前版本如下:

    """
    super simple neural networks (bear version)
      * input x is matrix 784x1 (where 784 = 28*28 which is MNIST image data)
      * 30 neurons for the only hidden layer, as n^{[1]} = 30
      * output layer: one neuron for classification(logistic)
      * using relu for activation function in hidden layer
    
    input layer:
        x: (shape 784x1)
            X: m samples where m = 12665 , X: 784 x 12665
            as : a^{[0]}: 784 x 12665  n^{[0]} = 784
    hidden layer:
        n^{[1]}:  30
        W^{[1]}: (30,784)   as (n^{[1]},n^{[0]})
        Z^{[1]}: (30,12665) as (n^{[1]},m)
        A^{[1]}: (30,12665) as (n^{[1]},m)
    output layer:
        n^{[2]}: 1
        W^{[2]}: (1,30)     as (n^{[2]},n^{[1]})
        Z^{[2]}: (1,12665)  as (n^{[2]},m)
        A^{[2]}: (1,12665)  as (n^{[2]},m)
    
    output:
        y \in [0,1] or  p \in {0,1}
            Y: (1 x m) ndarray
    structure for only one sample:
          x_1   ->   W*X + B   ->  relu  ->
          x_2   ->   W*X + B   ->  relu  ->  \
          ...   ->     ...     ->     .. ->  -> w*x+b -> logistic
          x_784 ->   W*X + B   ->  relu  ->  /
         ------     --------------------       ------------------
           |                |                          |
           V                V                          V
         input         30 neurons                 one neuron
        feature      relu activation             output layer
    
      By numpy with m samples:
        np.logistic(W2@g(W1@X+b1)+b2) as \hat{Y}: (1 x m) ndarray
    
        dimension analysis:
            W2        : (n2,n1)
            g(W1@X+b1): (n1,m)
                W1 : (n1,n0)
                X  : (n0,m)
                b1 : (n1,1)  with broadcasting to (n1,m)
            b2: (n2,1) with broadcasting to (n2,m)
    
    grad and notaion:
        forward propagation : A1 A2 Z1 Z2
        backward propagation: dW1 dW2 db1 db2
    
        more details:
            Z1 = W1@X  + b1
            Z2 = W2@A1 + b2
            A1 = g(Z1)      -- g     for relu
            A2 = \sigma(Z2) -- sigma for logistic
    
            dW2 = ((1/m)*(A2-Y))@A1.T
                dW2 = dZ2@A1.T  where dZ2 = (1/m)*(A2-Y)
                A2.shape:(1,m) Y.shape:(1,m) A1.T.shape:(n1,m)
                so: dW2.shape: (1,n1)
    
            dW1 = (W2.T@((1/m)*(A2-Y))*g_prime(Z1))@A0.T
                dW1 = dZ1@A1.T
                    where
                        dZ1 = W2.T@dZ2 * g_prime(Z1)
                        g_prime is derivative of relu
                    dW2.shape: (n1,n0)
            note: @ for matrix multiply;   * for dot product/element-wise
    
    Challenges
        1. Understanding the MNIST dataset and labels
        2. Understanding gradient caculate and the gradient descent
        3. Understanding logistic regression loss function and the caculation
        3. Knowing feature normalization
    
    about it:
        it's a simple project for human learning how machine learning
        version ant : scalar input/one neuron/one layer/binary classification
        version bear: vector input/30+1 neurons /two layer/binary classification
        by orczhou.com
    """
    
    from keras.datasets import mnist
    import numpy as np
    
    # return only data lable 0 or 1 from MNIST for the binary classification
    def filter_mnist_data(data_X, data_y):
        data_filter = np.where((data_y == 0) | (data_y == 1))
        filtered_data_X, filtered_data_y = data_X[data_filter], data_y[data_filter]
        r_data_X = filtered_data_X.reshape(filtered_data_X.shape[0],filtered_data_X.shape[1]*filtered_data_X.shape[2])
        return (r_data_X, filtered_data_y)
    
    (train_all_X, train_all_y), (test_all_X, test_all_y) = mnist.load_data()
    (train_X,train_y) = filter_mnist_data(train_all_X, train_all_y)
    (test_X ,test_y ) = filter_mnist_data(test_all_X, test_all_y)
    
    X  = train_X.T
    Y  = train_y.reshape(1,train_y.shape[0])
    
    m = X.shape[1]    # number of samples
    
    # hyper-parameter; read the comments above for structure of the NN
    n0 = X.shape[0]   # number of input features
    n1 = 10           # nerons of the hidden layer
    n2 = 1            # nerons of the output layer
    iteration_count = 500
    learning_rate   = 0.5
    
    # feature scaling / Normalization
    mean = np.mean(X,axis = 1,keepdims = True)
    std  = np.std( X,axis = 1,keepdims = True)+0.000000001
    X  = (X-mean)/std
    
    # initial parameters: W1 W2 b1 b2 size
    np.random.seed(561)
    W1 = np.random.randn(n1,n0)*0.01
    W2 = np.random.randn(n2,n1)*0.01
    b1 = np.zeros([n1,1])
    b2 = np.zeros([n2,1])
    
    # logistic function
    def logistic_function(x):
        return 1/(1+np.exp(-x))
    
    about_the_train = '''\
    try to train the model with:
      learning rate: {:f}
      iteration    : {:d}
      neurons in hidden layer: {:d}
    \
    '''
    print(about_the_train.format(learning_rate,iteration_count,n1))
    
    # forward/backward propagation (read the comment above:"grad and notaion")
    cost_last = np.inf # very large data,maybe better , what about np.inf
    for i in range(iteration_count):
        # forward propagation
        A0 = X
    
        Z1 = W1@X + b1  # W1 (n1,n0)  X: (n0,m)
        A1 = np.maximum(Z1,0) # relu
        Z2 = W2@A1 + b2
        A2 = logistic_function(Z2)
    
        dZ2 = (A2-Y)/m
        dW2 = dZ2@A1.T
        db2 = np.sum(dZ2,axis=1,keepdims = True)
    
        dZ1 = W2.T@dZ2*(np.where(Z1 > 0, 1, 0)) # np.where(Z1 > 0, 1, 0) is derivative of relu function
        dW1 = dZ1@A0.T
        db1 = np.sum(dZ1,axis=1,keepdims = True)
    
        cost_current = np.sum(-(Y*(np.log(A2))) - ((1-Y)*(np.log(1-A2))))/m
        if (i+1)%(iteration_count/20) == 0:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
    
        cost_last = cost_current
        W1 = W1 - learning_rate*dW1
        W2 = W2 - learning_rate*dW2
        b1 = b1 - learning_rate*db1
        b2 = b2 - learning_rate*db2
    
    print("Label:")
    print(np.round( Y[0][:20]+0.,0))
    print("Predict:")
    print(np.round(A2[0][:20],0))
    
    # Normalization for test dataset
    X  = (test_X.T - mean)/std
    Y  = test_y.reshape(1,test_y.shape[0])
    
    Y_predict = (logistic_function(W2@np.maximum((W1@X+b1),0)+b2) > 0.5).astype(int)
    
    for index in (np.where(Y != Y_predict)[1]):
        print(f"failed to recognize: {index}")
        # np.set_printoptions(threshold=np.inf)
        # np.set_printoptions(linewidth=np.inf)
        # print(test_X[index].reshape(28,28))
    
    print("total test set:" + str(Y.shape[1]) + ",and err rate:"+str((np.sum(np.square(Y-Y_predict)))/Y.shape[1]))

    简单统计,代码总计180行,其中注释约90行。

    运行结果说明

    运行该程序会有如下输出:

    try to train the model with:
      learning rate: 0.500000
      iteration    : 500
      neurons in hidden layer: 10
    iteration:    25,cost_current:0.009138,cost_last:0.009497,cost reduce:0.000358
    ...
    iteration:   500,cost_current:0.000849,cost_last:0.000850,cost reduce:0.000001
    failed to recognize: 1749
    failed to recognize: 2031
    total test set:2115,and err rate:0.0009456264775413711

    在本次训练后,对于测试集中的 2115 中图片,识别的错误率为 0.094%,即有两张图片未能够正确识别。这两张图片的编号分别是1749和2031,对应的图像如下:

    最后

    该神经网络在训练后,对于2115张测试集中的图片,仅有2张识别失败。通过该程序,可以看到,神经网络的强大与神奇。

    因为神经网络的 Backward Propagation 的推导非常复杂,本文中直接给出了相关公式,后续将独立介绍该部分内容。如果有任何问题可以留言讨论或公众号后台给作者留言。

    补充记录

    在这个该程序实现中,还遇到了一些其他的问题,记录如下。

    一些报错

    “RuntimeWarning: overflow encountered in exp: return 1/(1+np.exp(-x))”

    当训练或者预测计算过程中,如果出现部分-x值比较大,那么就会出现np.exp(-x)溢出的问题。不过,还好该问题在该场景下并不会影响计算结果,因为这种情况下,1/(1+np.exp(-x))依旧会返回0,而这也是预期中的。

    难以分析

    在这个识别程序中,最后发现有部分图片是难以识别的,而且无论如何调整参数,部分图片还是难以识别。最好的处理方法当然是增加训练数据集,但现实中可能并不总是能够做到。

    而关于模型应该如何优化,这似乎是难以分析,无从入手。

  • 在开始之前,也没想到99行代码就够了,原以为里面有个“梯度下降”,代码行数应该是数百级别吧,实际完成后,发现加上注释(约35%)才99行。


    “极简”神经网络的结构

    先来看“极简”有多简:

    • 一个输入层、一个输出层,中间没有隐藏层
    • 输入层的样本数据就一个特性(feature),总计6个样本
    • 这是二分类问题,所以输出层就一个输出值
       x ->   w*x + b   ->   logistic function  -> output
            ----------------------------------
                 |                    |
                 V                    V
             one neuron     activation function

    在后续的实现中,我们构造了六个样本用于该神经网络的训练。

    鉴于这个“极简”神经网络,没有任何隐藏层,所以,这也是一个典型的“logistic regression”问题。

    前置知识

    你需要了解如下的前置知识,以很好的理解该神经网络的实现与训练:

    • 了解神经网络的基础:浅层神经网络
    • 了解 梯度下降算法,了解基本最优化算法概念,了解链式法则
    • 了解 logistic function 的基本特性
    • 了解 Python 和 NumPy 的基本使用

    问题描述与符号约定

    用于训练的样本数据有\( (x,\hat{y}) \): \( (1,0)、(2,0)、(3,0)、(4,0)、(5,1)、(6,1) \)。

    一个具体的样本,在下面的公式中通常使用 \( (x^{(j)}, y^{(j)}) \)表示, 其中,\( j = 1…m \)。

    \( \hat{y} \)则表示根据参数计算出的预测值,更为具体的 \( \hat{y}^{(j)} \)表示 \(x = x^{(j)} \)时的预测值。

    构建目标函数

    从样本数据可以看到,这是一个二分类问题,可以使用logistic function作为输出层的激活函数,如果输出值大于0.5,则预测为1,否则预测为0

    对于任何一个样本,就可以如下函数作为logistic function的损失函数\( L \):

    $$
    L(y,\hat{y}) = – (yln(\hat{y}) + (1-y)ln(1-\hat{y}))
    $$

    所以,全局的目标函数就是:

    $$
    \begin{aligned}
    J(w,b) & = \frac{1}{m} \sum\limits_{j=0}^{m} L(y^{(j)},\hat{y}^{(j)}) \\
    & = \frac{1}{m} \sum\limits_{j=0}^{m} – (y^{(j)}ln(\hat{y}^{(j)}) + (1-y^{(j)})ln(1-\hat{y}^{(j)}))
    \end{aligned}
    $$

    其中 \( m \)表示总样本数量,这里取值是6。在这个极简的神经网络中 \( \hat{y} \)有如下计算表达式:

    $$
    \hat{y} = \frac{1}{1+e^{-(wx+b)}}
    $$

    最终,该神经网络的参数求解(也就是“训练”)过程,就是求解如下的极值问题:

    $$
    (w,b) = \min_{w, b} J(w, b) = \min_{w,b} \frac{1}{m} \sum\limits_{j=0}^{m} L(y^{(j)},\hat{y}^{(j)})
    $$

    目标函数计算的具体代码如下:

    def cost_function(w_p,b_p,x_p,y_p):
        c = 0
        for i in range(m):
            y = function_f(x_p[i],w_p,b_p)
            c += -y_p[i]*math.log(y) - (1-y_p[i])*math.log(1-y)
        return c

    梯度计算

    前面介绍了很多梯度的内容,这里不再详述。在这个具体的问题中,需要求解的梯度为:

    $$
    (\frac{\partial J}{\partial w},\frac{\partial J}{\partial b})
    $$

    在这里,简单展示该梯度的计算,主要需要使用的是链式法则和基本的求导/微分运算。

    首先,为了便于计算,这里记:

    $$
    \begin{aligned}
    \hat{y} & = \frac{1}{1+e^{-z}} \\
    z & = w*x + b
    \end{aligned}
    $$

    所以,根据链式法则容易有:

    $$
    \frac{\partial L}{\partial w} = \frac{\partial L}{\partial \hat{y}} * \frac{\partial \hat{y}}{\partial z} * \frac{\partial z}{\partial w} \\
    \frac{\partial L}{\partial b} = \frac{\partial L}{\partial \hat{y}} * \frac{\partial \hat{y}}{\partial z} * \frac{\partial z}{\partial b}
    $$

    这其中,\( \frac{\partial L}{\partial \hat{y}} \) 和 \( \frac{\partial \hat{y}}{\partial z} \)略有一些计算量,\( \frac{\partial z}{\partial w} \) 和\( \frac{\partial z}{\partial b} \)比较简单,具体的:

    $$
    \begin{aligned}
    \frac{\partial L}{\partial \hat{y}} * \frac{\partial \hat{y}}{\partial z} & = \hat{y} – y \\
    \frac{\partial z}{\partial w} & = x \\
    \frac{\partial z}{\partial b} & = 1
    \end{aligned}
    $$

    所以,最终的梯度计算公式如下:

    $$
    \begin{aligned}
    \frac{\partial J}{\partial w} & = \frac{1}{m} \sum\limits_{j=1}^{m} (\hat{y}^{(j)} – y^{(j)})*x^{(j)} \\
    \frac{\partial J}{\partial b} & = \frac{1}{m} \sum\limits_{j=1}^{m} (\hat{y}^{(j)} – y^{(j)})
    \end{aligned}
    $$

    在实际的计算中,先通过正向传播(Forward Propagation)计算出\( \hat{y}^{(j)} \),然后在计算出梯度。此外,可以使用NumPyndarray简化表达,同时增加计算的并行性。这里为便于理解,全部都使用标量计算,在文章的最后也提供了NumPy的对应实现。

    正向传播计算如下:

    # function_f: 
    # x   : scalar
    # w   : scalar
    # b   : scalar
    def function_f(x,w,b):  
        return 1/(1+math.exp(-(x*w+b)))

    梯度(反向传播)计算如下:

    # Gradient caculate 
    # x_p: x_train
    # y_p: y_train
    # w_p: current w
    # b_p: current b
    def gradient_caculate(x_p,y_p,w_p,b_p):
        gradient_w,gradient_b = (0.,0.)
        for i in range(m):
            gradient_w += x_p[i]*(function_f(x_p[i],w_p,b_p)-y_p[i])
            gradient_b += function_f(x_p[i],w_p,b_p)-y_p[i]
        return gradient_w,gradient_b

    梯度下降迭代

    这里设置迭代次数为50000次,学习率设置为0.01,当迭代目标函数变化值小于0.000001时也停止迭代(这并不是必须的)。具体的:

    iteration_count = 50000
    learning_rate = 0.01
    cost_reduce_threshold = 0.000001

    于是又如下梯度下降迭代过程的代码:

    cost_last = 0
    for i in range(iteration_count):
        grad_w,grad_b = gradient_caculate(x_train,y_train,w,b)
        w = w - learning_rate*grad_w
        b = b - learning_rate*grad_b
        cost_current = cost_function(w,b,x_train,y_train)
        if i >= iteration_count/2 and cost_last - cost_current<= cost_reduce_threshold:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
            break
        if (i+1)%(iteration_count/10) == 0:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
        cost_last = cost_current

    预测

    完成训练后,则可以对输入值进行预测。代码如下:

    print("after the training, parameter w = {:f} and b = {:f}".format(w,b))
    
    for i in range(m):
        y = function_f(x_train[i],w,b)
        p  = 0
        if y>= 0.5: p  = 1
        print("sample: x[{:d}]:{:d},y[{:d}]:{:d}; the prediction is {:d} with probability:{:4f}".format(i,x_train[i],i,y_train[i],p,y))

    上述代码会产生如下的输出:

    after the training, parameter w = 5.056985 and b = -22.644516
    sample: x[0]:0,y[0]:0; the prediction is 0 with probability:0.000000
    sample: x[1]:1,y[1]:0; the prediction is 0 with probability:0.000000
    sample: x[2]:2,y[2]:0; the prediction is 0 with probability:0.000004
    sample: x[3]:3,y[3]:0; the prediction is 0 with probability:0.000568
    sample: x[4]:4,y[4]:0; the prediction is 0 with probability:0.081917
    sample: x[5]:5,y[5]:1; the prediction is 1 with probability:0.933417

    可以看到,在完成训练后的这个极简神经网络能够较为准确的预测样本中的数据。

    完整的代码

    完成在的代码可以在 GitHub 上查看:https://github.com/orczhou/ssnn/ 。包括了三个段程序:

    • ssnn_ant.py : 最为基础的上述神经网络的实现
    • ssnn_ant_np.py : 使用numpy对上述实现进行向量化
    • ssnn_ant_tf.py : 使用 TensorFlow 框架实现上述程序

    这里也简单列出相关程序如下(最新代码可以参考上述 GitHub 仓库):

    ssnn_ant.py
    """
    super simple neural networks 
      * only one neuron in only the one hidden layer
      * input x is scalar (one-dimension)
      * using logistic function as the activation function
    
    input layer:
        x: scalar 
    parameters: 
        w: scalar
        b: scalar
    output:
        y \in [0,1] or  p \in {0,1}
    structure:
             
       x ->   w*x + b   ->   logistic function  -> output
            -----------      -----------------
                 |                    |
                 V                    V
             one neuron     activation function
    
    about it:
        it's a simple project for human learning how machine learning 
        by orczhou.com
    """
    import numpy as np
    import math
    
    # function_f: 
    # x   : scalar
    # w   : scalar
    # b   : scalar
    def function_f(x,w,b):  
        return 1/(1+math.exp(-(x*w+b)))
    
    # initial w,b
    w,b = (0,0)
    
    # samples
    x_train = np.array([0,1,2,3,4,5])
    y_train = np.array([0,0,0,0,0,1])
    #y_train = np.array([0,0,0,1,1,1])
    
    # m for sample counts
    m = x_train.shape[0]
    
    iteration_count = 50000
    learning_rate   = 0.01
    cost_reduce_threshold = 0.000001
    
    # Gradient caculate 
    # x_p: x_train
    # y_p: y_train
    # w_p: current w
    # b_p: current b
    def gradient_caculate(x_p,y_p,w_p,b_p):
        gradient_w,gradient_b = (0.,0.)
        for i in range(m):
            gradient_w += x_p[i]*(function_f(x_p[i],w_p,b_p)-y_p[i])
            gradient_b += function_f(x_p[i],w_p,b_p)-y_p[i]
        return gradient_w,gradient_b
    
    def cost_function(w_p,b_p,x_p,y_p):
        c = 0
        for i in range(m):
            y = function_f(x_p[i],w_p,b_p)
            c += -y_p[i]*math.log(y) - (1-y_p[i])*math.log(1-y)
        return c
    
    about_the_train = '''\
    try to train the model with:
      learning rate: {:f}
      max iteration : {:d}
      cost reduction threshold: {:f}
    \
    '''
    print(about_the_train.format(learning_rate,iteration_count,cost_reduce_threshold))
    
    # start training
    cost_last = 0
    for i in range(iteration_count):
        grad_w,grad_b = gradient_caculate(x_train,y_train,w,b)
        w = w - learning_rate*grad_w
        b = b - learning_rate*grad_b
        cost_current = cost_function(w,b,x_train,y_train)
        if i >= iteration_count/2 and cost_last - cost_current<= cost_reduce_threshold:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
            break
        if (i+1)%(iteration_count/10) == 0:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
        cost_last = cost_current
    
    print("after the training, parameter w = {:f} and b = {:f}".format(w,b))
    
    for i in range(m):
        y = function_f(x_train[i],w,b)
        p  = 0
        if y>= 0.5: p  = 1
        print("sample: x[{:d}]:{:d},y[{:d}]:{:d}; the prediction is {:d} with probability:{:4f}".format(i,x_train[i],i,y_train[i],p,y))

    使用NumPy向量化 ssnn_ant_np.py
    """
    super simple neural networks(using numpy,snn.py not using numpy)
      * only one neuron in only the one hidden layer
      * input x is scalar (0-dimension)
      * using logistic function as the activation function
    
    input layer:
        x: scalar
    parameters:
        w: scalar
        b: scalar
    output:
        y \in [0,1] or  p \in {0,1}
    structure:
    
       x ->   w*x + b   ->   logistic function  -> output
            -----------      -----------------
                 |                    |
                 V                    V
             one neuron     activation function
    
    about it:
        it's a simple project for human learning how machine learning
        by orczhou.com
    """
    import numpy as np
    import math
    
    # function_f:
    # x   : scalar or ndarray
    # w   : scalar
    # b   : scalar
    def function_f(x,w,b):
        return 1/(1+np.exp(-(x*w+b)))
    
    # initial w,b
    w,b = (0,0)
    
    # samples
    x_train = np.array([0,1,2,3,4,5])
    y_train = np.array([0,0,0,0,0,1])
    #y_train = np.array([0,0,0,1,1,1])
    
    # m for sample counts
    m = x_train.shape[0]
    
    iteration_count = 50000
    learning_rate   = 0.01
    cost_reduce_threshold = 0.000001
    
    # Gradient caculate
    # w_p: current w
    # b_p: current b
    def gradient_caculate(w_p,b_p):
        gradient_w = np.sum((function_f(x_train,w_p,b_p) - y_train)*x_train)
        gradient_b = np.sum(function_f(x_train,w_p,b_p) - y_train)
        return gradient_w,gradient_b
    
    def cost_function(w_p,b_p,x_p,y_p):
        hat_y = function_f(x_p,w_p,b_p)
        c = np.sum(-y_p*np.log(hat_y) - (1-y_p)*np.log(1-hat_y))
        return c/m
    
    about_the_train = '''\
    try to train the model with:
      learning rate: {:f}
      max iteration : {:d}
      cost reduction threshold: {:f}
    \
    '''
    print(about_the_train.format(learning_rate,iteration_count,cost_reduce_threshold))
    
    # start training
    cost_last = 0
    for i in range(iteration_count):
        grad_w,grad_b = gradient_caculate(w,b)
        w = w - learning_rate*grad_w
        b = b - learning_rate*grad_b
        cost_current = cost_function(w,b,x_train,y_train)
        if i >= iteration_count/2 and cost_last - cost_current<= cost_reduce_threshold:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
            break
        if (i+1)%(iteration_count/10) == 0:
            print("iteration: {:5d},cost_current:{:f},cost_last:{:f},cost reduce:{:f}".format( i+1,cost_current,cost_last,cost_last-cost_current))
        cost_last = cost_current
    
    print("after the training, parameter w = {:f} and b = {:f}".format(w,b))
    
    for i in range(m):
        y = function_f(x_train[i],w,b)
        p  = 0
        if y>= 0.5: p  = 1
        print("sample: x[{:d}]:{:d},y[{:d}]:{:d}; the prediction is {:d} with probability:{:4f}".format(i,x_train[i],i,y_train[i],p,y))

    使用TensorFlow实现该功能

    这里也是使用 TensorFlow 对上述问题中的数据进行训练并预测。详细代码和 TensorFlow 输出参考小结“TensorFlow 代码”和“TensorFlow的输出”。

    这里对其训练结果做简要的分析。在输出中,可以看到训练后的参数分别是:\( w = 1.374991 \quad b = -5.9958787 \),那么对应的预测表达式为:

    $$ \frac{1}{1+e^{-(w*x+b)}} $$

    代入 \( x = 1 \),其计算结果为:\( np.float64(0.009748092866213252) \),这与 TensorFlow 输出的 \( [0.00974809] \) 是一致的,这也验证了训练程序其实现与理解的事完全一致的。

    TensorFlow 代码 ssnn_ant_tf.py
    import tensorflow as tf
    import numpy as np
    
    tf.random.set_seed(1)
    X_train = np.array([[1], [2], [3], [4], [5],[6]], dtype=float)
    y_train = np.array([[0], [0], [0], [0], [1],[1]], dtype=float)
    
    model = tf.keras.Sequential([
        tf.keras.layers.Input(shape=(1,)),
        tf.keras.layers.Dense(units=1, activation='sigmoid')
    ])
    
    # model.compile(optimizer='adam', loss='binary_crossentropy', metrics=['accuracy'])
    model.compile(optimizer=tf.keras.optimizers.SGD(learning_rate=0.1), loss='binary_crossentropy', metrics=['accuracy'])
    
    model.fit(X_train, y_train, epochs=1000, verbose=0)
    model.summary()
    
    model.evaluate(X_train,  y_train, verbose=2)
    
    predictions = model.predict(X_train)
    print("Predictions:", predictions)
    
    for layer in model.layers:
        weights, biases = layer.get_weights()
        print("weights::", weights)
        print("biases:", biases)
    TensorFlow的输出
    Model: "sequential"
    ┏━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┳━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┳━━━━━━━━━━━━━━━━━┓
    ┃ Layer (type)                         ┃ Output Shape                ┃         Param # ┃
    ┡━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━╇━━━━━━━━━━━━━━━━━━━━━━━━━━━━━╇━━━━━━━━━━━━━━━━━┩
    │ dense (Dense)                        │ (None, 1)                   │               2 │
    └──────────────────────────────────────┴─────────────────────────────┴─────────────────┘
     Total params: 4 (20.00 B)
     Trainable params: 2 (8.00 B)
     Non-trainable params: 0 (0.00 B)
     Optimizer params: 2 (12.00 B)
    1/1 - 0s - 32ms/step - accuracy: 1.0000 - loss: 0.1856
    1/1 ━━━━━━━━━━━━━━━━━━━━ 0s 10ms/step
    Predictions: 
     [[0.00974809]
     [0.03747462]
     [0.13343701]
     [0.37850127]
     [0.7066308 ]
     [0.90500087]]
    weights:: [[1.374991]]
    biases: [-5.9958787]

    补充说明

    一个一般的前馈神经网络(FNN)通常至少需要一个隐藏层,并且在隐藏层有多个神经元。一个具有多个神经元的多层网络的结构和训练,其复杂度会高很多,后续会再做介绍。本文实现代码虽然只有99行代码,去掉注释只有70行代码,但麻雀虽小、五脏俱全,包含了梯度下降的实现、链式法则的应用等,如果理解了该示例,则可以很好帮助开发者打好基础,以便更好的理解更为复杂的神经网络的训练。

    此外,在计算中省略了一些求导计算,其中略微有一些复杂度的是 \( \frac{\partial L}{\partial \hat{y}} * \frac{\partial \hat{y}}{\partial z} \),感兴趣的可以自行补全。