python手写神经网络之BatchNormalization实现

也是来源于《深度学习入门——基于Python的理论与实现》附加代码，书中只是给了BN的对比结果，展示了BN的效果，没有再赘述实现（可能因为有点复杂），所以这里研究一下BN的代码。

之前我曾经使用过TensorFlow的BN，它提供了两三种接口，透明程度和使用方法不相同，有的是透明到你可以自定义参数并传给BN层，然后训练参数，也有只定义一个层，全自动使用的，但是都没有自己纯手写一个python实现更透彻。而且TensorFlow有很多高级封装和其他功能的整合，乱花渐欲迷人眼。不如来一个纯版的实现更好。

关于理论基础，这里就不赘述了，先将x归一化，再通过学习到的参数缩放平移x，达到分布归一与分布还原拟合的trade-off，公式如下：

可以看到，正向传播有两个计算参数mu和sigma，然后完成其余归一化计算，最后再通过两个可学习参数beta和gamma进行分布恢复。

反向传播，反向传播有两个关注点，计算gamma和beta的梯度，以便进行更新；计算对输入（x）的梯度，以便完成梯度的传导工作。

如果是工程版的实现，根据Tensorflow接口的功能，还应该分训练和测试阶段，训练阶段应该进行一个mu（mean）和sigma（var）的总平均值的更新。这个好像没悬疑，必须要有的，没有就不完整（失去了降噪和降低过拟合的功能）

实现框架，类似本书其他网络层，实现一个BN类封装，包含初始化、正向传播、反向传播等函数。

其实难点就是把公式逐步拆解，再把所有的导数公式都找出来，这一步最好在纸上做。给每个中间变量都想好名字。
因为考虑不一定周到，比如维度（convNet默认（N,C,H,W）），时间关系，我就不重新造轮子了，主要是先自己列个大纲和伪代码，然后解读实例代码，看自己有哪些遗漏或者想不到的，他细节都是怎么实现的，算一个学习过程，为了直观，解读以注释形式给出。
这个方法的收获还是有很多的，反向传播的逐个计算可能是最难的，或者叫繁琐，我认为本质上和sigmoid或者affine层或者softmax_with_loss是一样的，如果有前边的基础，这个只是有空再练手，其实并不影响你理解，所以我跳过了，没在纸面上推。

下面是代码解读

• 前向传播的外层封装，卷积网络被改造成N*D的模式，最后恢复（好像用了channel_first模式）：

    def forward(self, x, train_flg=True):#封装：修改卷积模式到FC模式，
        # BN应该是有一个channel一个平均值的模式？这里channel被平铺了，每个feature“像素点”一个平均值。
        # 所以这可能不算最好的conv实例，权当是个FC罢,就是遇到Conv强行转FC的一个形式
        self.input_shape = x.shape
        if x.ndim != 2:
            N, C, H, W = x.shape#用了channel_first
            x = x.reshape(N, -1)

        out = self.__forward(x, train_flg)

        return out.reshape(*self.input_shape)#恢复原样

• 前向传播内层实现：

    def __forward(self, x, train_flg):
        if self.running_mean is None:#全局平均值的初始化，作者的参考代码全都是根据第一次前向传播确定形状的模式
            N, D = x.shape
            #这里也可以看出，是根据每个feature一个平均值，所以好像ConvNet不友好？但是BN层又不该改变维度，做了平均值最终也是要恢复的，注意，主要是平均值数值本身的不同，一个channel一个平均值理论上更平均一些罢了。
            self.running_mean = np.zeros(D)
            self.running_var = np.zeros(D)

        if train_flg:#训练时（在纸上做公式的分解）
            mu = x.mean(axis=0)#计算N个xi的平均值，N坍塌，D维持不变
            xc = x - mu#x零中心化（N,D）
            var = np.mean(xc ** 2, axis=0)#分子：计算标准差（简化了，因为xc=x-mu是零中心，所以xi-mu就成了xc，然后np.mean）
            std = np.sqrt(var + 10e-7)#分母，标准差，防零
            xn = xc / std#这是标准化后的输出，恢复分布的操作out= self.gamma*xn + self.beta提到了if外边，为了和else共享代码

            self.batch_size = x.shape[0]
            self.xc = xc#这些都被记录了下来，为了反向传播方便使用
            self.xn = xn
            self.std = std
            # 这两个是记录全局平均，我也好奇全局平均怎么记录，总不能一直累加，最后测试时记住数量除一次？更不能每次简单对冲进去？因为batch不一样！
            # 就算除以batch再去更新，那么权重也不太对，之前batch=10000被算成1份，现在batch=10也被算成一份，就冲淡了10000份？所以看到，他是用了momentum来做滑动平均的更新。
            # 这些细节不去实现真的很不容易想到！！！！
            self.running_mean = self.momentum * self.running_mean + (1 - self.momentum) * mu
            self.running_var = self.momentum * self.running_var + (1 - self.momentum) * var
        else:#这个简单：测试时候直接用全局平均值套公式计算一下就行了
            xc = x - self.running_mean
            xn = xc / ((np.sqrt(self.running_var + 10e-7)))

        out = self.gamma * xn + self.beta#共享操作：恢复分布
        return out

• 反向传播外层

    def backward(self, dout):#同样的，封装了一个ConvNet维度处理的过程
        if dout.ndim != 2:
            N, C, H, W = dout.shape
            dout = dout.reshape(N, -1)

        dx = self.__backward(dout)

        dx = dx.reshape(*self.input_shape)
        return dx

• 反向传播内层

    def __backward(self, dout):
        #beta和gamma是最外层的两个操作，反向传播最简单，由xn*gamma+beta=out可知
        dbeta = dout.sum(axis=0)#细节：x的坍缩，一个batch合并成一个值。
        dgamma = np.sum(self.xn * dout, axis=0)#同上
        #下面的不墨迹了，其实就是前向传播的过程，逐步反推回去
        dxn = self.gamma * dout
        dxc = dxn / self.std
        dstd = -np.sum((dxn * self.xc) / (self.std * self.std), axis=0)
        dvar = 0.5 * dstd / self.std
        dxc += (2.0 / self.batch_size) * self.xc * dvar
        dmu = np.sum(dxc, axis=0)
        dx = dxc - dmu / self.batch_size#用来传递给前一层

        self.dgamma = dgamma#本层需要的结果，保存，用于更新参数
        self.dbeta = dbeta#本层需要的结果，保存，用于更新参数

        return dx

• 最后简单带过init

    def __init__(self, gamma, beta, momentum=0.9, running_mean=None, running_var=None):
        self.gamma = gamma
        self.beta = beta
        self.momentum = momentum#用来做全局平均的衰减
        self.input_shape = None  # Conv层的情况下为4维，全连接层的情况下为2维

        # 测试时使用的全局平均值和方差
        self.running_mean = running_mean
        self.running_var = running_var

        # backward时使用的中间数据
        self.batch_size = None
        self.xc = None
        self.std = None
        self.dgamma = None
        self.dbeta = None

• 这样，外边加上class BatchNormalization就完成了。

• 最后就是使用（需要在外部定义params传给BN类进行初始化，这个类似TF给出的一种半透明的使用方法，可以看到gamma默认1，beta默认0，对xn的分布不产生影响，之后具体什么值就要学习了）

        # 初始化权重
        self.__init_weight(weight_init_std)

        # 生成层
        activation_layer = {'sigmoid': Sigmoid, 'relu': Relu}
        self.layers = OrderedDict()
        for idx in range(1, self.hidden_layer_num+1):
            self.layers['Affine' + str(idx)] = Affine(self.params['W' + str(idx)],
                                                      self.params['b' + str(idx)])
            if self.use_batchnorm:
                self.params['gamma' + str(idx)] = np.ones(hidden_size_list[idx-1])
                self.params['beta' + str(idx)] = np.zeros(hidden_size_list[idx-1])
                self.layers['BatchNorm' + str(idx)] = BatchNormalization(self.params['gamma' + str(idx)], self.params['beta' + str(idx)])
                
            self.layers['Activation_function' + str(idx)] = activation_layer[activation]()
            
            if self.use_dropout:
                self.layers['Dropout' + str(idx)] = Dropout(dropout_ration)

        idx = self.hidden_layer_num + 1
        self.layers['Affine' + str(idx)] = Affine(self.params['W' + str(idx)], self.params['b' + str(idx)])

        self.last_layer = SoftmaxWithLoss()

• 求梯度

    def numerical_gradient(self, X, T):
        loss_W = lambda W: self.loss(X, T, train_flg=True)

        grads = {}
        for idx in range(1, self.hidden_layer_num+2):
            grads['W' + str(idx)] = numerical_gradient(loss_W, self.params['W' + str(idx)])
            grads['b' + str(idx)] = numerical_gradient(loss_W, self.params['b' + str(idx)])
            
            if self.use_batchnorm and idx != self.hidden_layer_num+1:
                grads['gamma' + str(idx)] = numerical_gradient(loss_W, self.params['gamma' + str(idx)])
                grads['beta' + str(idx)] = numerical_gradient(loss_W, self.params['beta' + str(idx)])

        return grads

• 最后就是和其他一样的

w-=lr*grad

• 测试时注意一下参数要写成False：

def predict(self, x, train_flg=False):

小结：

BN层分两个阶段
公式的分解稍微繁琐

ConvNet的优化问题（私以为，这个代码不代表全部情况）