Theano学习四：单层简单线性模型识别手写字

深度学习是机器学习的一个研究领域,本文作为介绍theano深度学习框架的一部分,会介绍深度学习的一些概念,如损失函数,交叉熵,反向传播,梯度下降等.

1.1MNIST数据集

美国国家标准与技术研究所改进(MNIST)数据集是一个用于训练和测试分类模型的非常著名的手写体数字数据集{0,1,2,3,4,5,6,7,8,9}.
分类模型是指在给定输入下预测所观察类的概率的模型.
训练是一种参数学习任务,以使得模型适合于数据,并可以对任何输入图像预测正确的标签.对于该训练任务,MNIST数据集中包含了60000幅图像,并对每个示例具有目标标签(0~9之间的数字).
为验证训练是否有效并决定何时停止训练,通常将训练数据集分为来能个数据集:80%~90%的图像用于训练,剩下10%~20%的图像没有用于训练算法,而是验证该模型是否可推广到未观察的数据.
另外,还有一个在训练过程中从未经历过算法的独立数据集,称为测试集,是由MNIST数据集中的10000幅图像组成的.
在MNIST数据集中,每个示例的输入数据是28*28的标准化单色图像以及由0~9之间简单整数表示的标签.现在,展示其中的一些图像:
一.首先,下载数据集的封装版本,以便更容易在Python中加载:

wget http://www.iro.umontreal.ca/~lisa/deep/data/mnist/mnist.pkl.gz

二.然后将数据加载到一个Python会话中:

import pickle, gzip
import matplotlib
import numpy
import matplotlib.pyplot as plt

with gzip.open('/home/wc/test/mnist.pkl.gz', 'rb') as f:
    train_set, valid_set, test_set = pickle.load(f, encoding='latin1')

print(len(train_set[0]), len(valid_set[0]), len(test_set[0]))#50000 10000 10000
print(train_set[1].shape)       #就是一个数字
print(train_set[0].shape[0], train_set[0].shape[1])
plt.rcParams['figure.figsize'] = (20, 20)
plt.rcParams['image.cmap'] = 'gray'

for i in range(9):
    plt.subplot(1,10,i+1)
    plt.imshow(train_set[0][i].reshape(28,28))
    plt.axis('off')
    plt.title(str(train_set[1][i]))
plt.show()

运行结果显示了数据集前9个样本及其对应的标签(实际值,即预期的分类算法正确答案).

为了避免对GPU的过多传输,且由于完整数据集足够小以适应GPU的内存,通常将完整训练集置于共享变量中:

train_set_x = theano.shared(numpy.asarray(train_set[0],  dtype=theano.config.floatX))
train_set_y = theano.shared(numpy.asarray(train_set[1],  dtype='int32'))

1.2 训练程序架构
训练程序架构通常包括以下步骤:

• 1)脚本环境设置:如封装包的导入,GPU的使用等
• 2)数据加载:在训练过程中访问数据的数据加载类通常是以随机顺序加载,以避免出现太多同类中的相似样本,但是在先是简单样本而后是复杂样本的课程学习情况下,有时会以精确顺序加载数据.
• 3)数据预处理:是指进行一系列转换,如在图像上交换维度,增加模糊或噪声.增加一些数据增强的转换是很常见的,如随机剪裁,比例,亮度或对比度抖动,可得到更多的样本,并减少数据过拟合的风险.如果模型中自由参数的个数对于训练数据集大小非常重要,那么该模型可能会从所提供的样本中学习.此外,如果数据集太小,且对相同数据迭代执行过多,则该模型可能变得只针对训练样本,而不能很好地推广到新的未知样本.
• 4)模型构建:根据持续变量(共享变量)中的参数来定义模型结构,以便在训练期间更新其值以适应训练数据.
• 5)训练:在完整数据集上进行整体训练或每个样本依次单个训练都有相应的不同算法.通常按批次进行训练可得到最佳收敛,一个批次是从几十到几百个样本所组成的一个小样本子集.
  采用批处理的另一个原因是为了提高GPU的寻来呢速度,这是因为单个数据传输的成本太高,且GPU内存不足以保存整个数据集.GPU是一个并行架构,因此处理一批样本通常比逐一处理样本跟快.在一定程度上,同事观察更多样本会更快收敛[即墙上时间(in wall-time)-进程运行时间].即使GPU内存足够大足以保存整个数据集,也是如此:批大小递减通常会处理得越快,即处理较小的批要快于处理整个数据集.值得注意的是,对于目前最先进的CPU也是如此,只是优化的批大小通常较小.
• 6)在训练过程中,经过一定次数的迭代,通常需要使用一部分训练数据或未经过学习的验证数据集进行验证.在验证数据集上计算训练损失.虽然算法目的是在给定训练数据下减少损失,但并不能保证对未知数据的泛化.验证数据是用于估计泛话性能的未知数据.在训练数据不具有代表性,或是由于异常和未正确采样,或模型过拟合训练数据时,可能会缺乏泛华.
  验证数据一切正常,且当验证损失不在减少时停止训练,即使训练损失可能会继续下降:不值得进一步训练,且可能会导致过拟合.
• 7)保存模型参数并显示结果,如用于收敛分析的最佳训练/验证损失值,训练损失曲线等.

在分类情况下,需计算训练过程的准确率(正确分类的把百分比)或错误率(错误分类的百分比)以及损失.在训练结束时,混淆矩阵有助于评估分类器性能.

1.3分类损失函数

损失函数是一个训练过程中最小化的目标函数,以获得最优模型.现有多种不同的损失函数.
在分类问题中,目标是从k类中预测正确的类,由于需要测量各类的真实概率分布q和预测概率分布p之差,通常采用交叉熵:

$$\sum _{i=1}^n\sum _{c=1}^k{q}_i(c)log(p_i(c))$$
式子中,i是数据集中的样本索引;n是数据集中的样本个数;k是类的个数.
尽管各类的实际概率q(c|xi)未知,但实际上可根据经验分布进行简单近似,即随机抽取数据集中的一个样本.同理,任何预测概率p的交叉熵都可由经验交叉熵来近似:
$$cross-entropy=\frac{1}{N}\sum _{i=1}^{N}log(p_i)$$
式子中,p_i是模型估计的样本x_i正确分类的概率.
精度和交叉熵都是同一目标,但是测量内容不同.精度是量测预测分类的正确程度,而交叉熵是测量概率之差.交叉熵减少表明预测正确分类的概率更高,但精度可能保持不变或下降.
精度是离散的且不可微的,而交叉熵损失则是可用于训练模型的一个可微函数.

1.4 单层线性模型

最简单的模型是线性模型,其中对于每个分类c,输出是输入的线性组合:

且输出是无界的.
要得到和为1的概率分布pi,线性模型的输出可通过一个softmax函数:
s o f t m a x : { o c } c → { e c 0 ∑ ξ e ξ } c softmax:\lbrace o_c \rbrace_c\to \lbrace \frac{e_c^0}{\sum_{\xi}e^\xi} \rbrace_c softmax:{oc​}c​→{∑ξ​eξec0​​}c​
其中,对于一个输入x,分类c的估计概率可由向量表示:

p=softmax(Wx+b)

python实现中,具体如下:

batch_size = 600
n_in=28 * 28
n_out=10

x = T.matrix('x')
y = T.ivector('y')
W = theano.shared(
            value=numpy.zeros(
                (n_in, n_out),#784*10
                dtype=theano.config.floatX
            ),
            name='W',
            borrow=True
        )
b = theano.shared(
    value=numpy.zeros(
        (n_out,),
        dtype=theano.config.floatX
    ),
    name='b',
    borrow=True
)
model = T.nnet.softmax(T.dot(x, W) + b)

给定输入下的预测由最大可能分类(最大概率)确定:
y_pred = T.argmax(model, axis=1)
在单层线性模型中,信息从输入转移到输出:这称为前馈网络.给定输入下计算输出的过程称为前向传播.
该层是完全连接的,因为所有输出oc都是以相应系数的所有输入之和:

1.5 成本函数和误差

由模型给出预测概率的成本函数如下:

cost = -T.Mean(T.log(model)[T.arange(y.shape[0]), y])

误差是指与实际分类不同的预测分类数除以总的个数,这可由均值表示:

error = T.mean(T.neq(y_pred, y))

相反,精度是对应于正确分类个数除以预测总数.误差和精度之和为1.
对于其他类型的问题,其他损失函数及其实现如下表:

分类交叉熵 上述的等效实现方法	T.nnet.categorical_crossentropy(model, y_true).mean()
二值交叉熵 适用于输出只能取两种值{0,1}的情况 通常在预测概率的sigmoid激活函数之后使用	T.nnet.binary_crossentropy(model, y_true).mean()
均方误差 回归问题的L2范数	T.sqr(model-y_true).mean()
绝对平均误差 回归问题的L1范数	T.abs_(model-y_true).mean()
L1平滑 称为回归问题的离群损失	T.switch( T.lt(T.abs_(model-y_true), 1./sigma), 0.5sigmaT.sqr(model-y_true), T.abs_(model-y_true)-0.5/sigma). sum(axis=1).mean();
方铰链损失 特别用于非监督式问题	T.sqr(T.maximum(1.-y_true*model, 0)).mean()
铰链损失	T.maximum(1.-y_true*model, 0.).mean()

1.6反向传播算法和随机梯度下降

反向传播或误差反向传播是自适应连接权重的监督式学习算法中的最常用发发.
考虑误差或成本是权重W和b的函数时,可用梯度下降法来逼近成本函数的局部极小值,这主要是沿负误差梯度方向改变权重来实现的:
$$W\leftarrow W-\lambda \frac{\partial E}{\partial w}$$
式中,λ是学习率,一个定义下降速度的正常数.
每次前馈运行后更新变量的编译函数如下:

g_W = T.grad(cost=cost, wrt=W)
g_b = T.grad(cost=cost, wrt=b)

learning_rate=0.13
index = T.lscalar()

train_model = theano.function(
    inputs=[index],
    outputs=[cost,error],
    updates=[(W, W - learning_rate * g_W),(b, b - learning_rate * g_b)],
    givens={
        x: train_set_x[index * batch_size: (index + 1) * batch_size],
        y: train_set_y[index * batch_size: (index + 1) * batch_size]
    }
)

输入变量是批的索引,这是因为所有数据集都是在共享变量中一次传输给CPU的.训练主要是每个样本迭代赋予模型(迭代)并重复执行多次(周期):

n_epochs = 100 #轮数
n_train_batches = train_set[0].shape[0] // batch_size #训练集个数除以batch_size得出分多少批 50000//600=83 

n_iters = n_epochs * n_train_batches    #总批次
train_loss = numpy.zeros(n_iters)
train_error = numpy.zeros(n_iters)
for epoch in range(n_epochs):
    for minibatch_index in range(n_train_batches):
        iteration = minibatch_index + n_train_batches * epoch
        train_loss[iteration], train_error[iteration] = train_model(minibatch_index)

        if (epoch* train_set[0].shape[0] + minibatch_index) %print_every == 0 :
            print('epoch {}, minibatch {}/{}, validation error {:02.2f} %, validation loss {}'.format(
                epoch,
                minibatch_index + 1,
                n_train_batches,
                valid_error[val_iteration] * 100,
                valid_loss[val_iteration]
            ))

尽管这只是给出了一个小批量的损失和误差,但对于整个数据集的平均值也不错.
在第一次迭代中误差率快速下降,然后逐渐减慢.
在GPU为GeForce GTX 2070的台式计算机上执行时间为3.6s.而在Intel i7 9700 CPU上,执行时间为7.90s.
经过一段时间后,模型收敛到5.3%~5.5%的误差率,且迭代次数更多,误差率会进一步减小,但同时也会导致过拟合.过拟合是指模型能够很好地拟合训练数据但对于未知数据却无法得到同样的误差率.
在这种情况下,模型是过于简单而只能拟合这些数据.
太简单的模型学习性能不佳.深度学习的原理是增加更多的层,即增加深度,由此构建更深度的网络来获得更好的精度.

最后贴出整段代码：

#from __future__ import print_function
import numpy
from theano import theano
import theano.tensor as T
import pickle, gzip
import timeit

data_dir = "/home/wc/test/"

print("Using device", theano.config.device)

print("Loading data")
with gzip.open(data_dir + "mnist.pkl.gz", 'rb') as f:
    train_set, valid_set, test_set = pickle.load(f, encoding='bytes')#训练集50000 验证集10000 测试集10000
    
train_set_x = theano.shared(numpy.asarray(train_set[0],  dtype=theano.config.floatX))   #图片集 50000张图片
train_set_y = theano.shared(numpy.asarray(train_set[1],  dtype='int32'))                        #结果集 50000张图片对应的结果

print("Building model")

batch_size = 600    #每个batch 训练600张图片
n_in=28 * 28        #784, 每隔手写图片有784个像素点,每个像素点存储的是一个0~1之间的数,每个图片是一个1*784的张量
n_out=10                #输出 10, 结果有10种可能

x = T.matrix('x')
y = T.ivector('y')
W = theano.shared(
            value=numpy.zeros(
                (n_in, n_out),#784*10
                dtype=theano.config.floatX
            ),
            name='W',
            borrow=True#在CPU上跑时,false是深拷贝, true是浅拷贝; 在GPU上跑都是深拷贝
        )
b = theano.shared(
    value=numpy.zeros(
        (n_out,),   #注意默认一行,即1*10
        dtype=theano.config.floatX
    ),
    name='b',
    borrow=True
)
model = T.nnet.softmax(T.dot(x, W) + b) #计算出来的结果是一个1*10的向量,每个分量的值在0~1之间,10个分量加起来的和刚好为1
#x应该不会只喂入一个数据，应该是batch个数据
y_pred = T.argmax(model, axis=1)    #model最大值列下标索引，也就是输入图片预测的值
error = T.mean(T.neq(y_pred, y))    #一个向量和一个标量怎么能比??

cost = -T.mean(T.log(model)[T.arange(y.shape[0]), y])
g_W = T.grad(cost=cost, wrt=W)
g_b = T.grad(cost=cost, wrt=b)

learning_rate=0.13
index = T.lscalar()

train_model = theano.function(
    inputs=[index],
    outputs=[cost,error],
    updates=[(W, W - learning_rate * g_W),(b, b - learning_rate * g_b)],
    givens={#一个一个小批量地喂数据
        x: train_set_x[index * batch_size: (index + 1) * batch_size],
        y: train_set_y[index * batch_size: (index + 1) * batch_size]
    }
)

validate_model = theano.function(#验证模型,验证训练模型训练出来的参数W,b
    inputs=[x,y],
    outputs=[cost,error]
)

print("Training")

n_epochs = 100 #轮数
n_train_batches = train_set[0].shape[0] // batch_size #训练集个数除以batch_size得出分多少批 50000//600=83 

n_iters = n_epochs * n_train_batches    #总批次
train_loss = numpy.zeros(n_iters)
train_error = numpy.zeros(n_iters)

validation_interval = 100   #一共训练83*100个interval,没训练100个interval 验证一下参数
n_valid_batches = valid_set[0].shape[0] // batch_size#验证集个数除以batch_size  结果是分多少批验证 10000//60 == 166
nums = (int) (n_iters/validation_interval)  #总批次除以验证间隔 得到验证次数 100*83/100 = 83
valid_loss = numpy.zeros(nums)
valid_error = numpy.zeros(nums)

start_time = timeit.default_timer()
for epoch in range(n_epochs):#0~99
    for minibatch_index in range(n_train_batches):  #0~82
        iteration = minibatch_index + n_train_batches * epoch   #batch编号,从0开始
        train_loss[iteration], train_error[iteration] = train_model(minibatch_index)    #训练

        if (iteration+1) % validation_interval == 0 :   #每100个minibatch验证一次
            val_iteration = iteration // validation_interval    #第几次验证
            valid_loss[val_iteration], valid_error[val_iteration] = numpy.mean([
                    validate_model( #验证
                        valid_set[0][i * batch_size: (i + 1) * batch_size],
                        numpy.asarray(valid_set[1][i * batch_size: (i + 1) * batch_size], dtype="int32")
                        )
                        for i in range(n_valid_batches) #每次验证的时候,将验证集所有数据验证一轮
                     ],axis=0)

            print('epoch {}, minibatch {}/{}, validation error {:02.2f} %, validation loss {}'.format(
                epoch,
                minibatch_index + 1,
                n_train_batches,
                valid_error[val_iteration] * 100,
                valid_loss[val_iteration]
            ))
end_time = timeit.default_timer()
print( end_time -start_time )

numpy.save("simple_train_loss", train_loss)
numpy.save("simple_valid_loss", valid_loss)