PRML第五章读书笔记——Neural Networks 二次下降、Hessian矩阵的近似求解和精确求解、正切传播/Tikhonov正则化/软权值共享、混合密度网络、贝叶斯神经网络

（神经网络发展很快，书中的思想停留在2006年左右那个时代。所以这里只对其中有价值的部分进行摘录。对于一般的知识略去，有更好的介绍书籍）

5.1 Feed-forward Network Functions

P231 单隐层神经网络的拟合能力

这里只用了三个隐层单元，tanh激活，就已经能拟合的这么好了。感受一下。

分类能力

5.2 Network Training

P235 逻辑回归中的交叉熵和平方损失

在逻辑回归中，用sum-of-squares做loss其实也还凑合（注意和感知机做区分，感知机因为没有sigmoid，所以直接用平方损失会受离群点影响）。不过Simard et al.(2003)指出，对分类问题，用交叉熵会比平方损失更快，也能提高泛化性

P239 梯度下降和二次下降的效率对比

如果参数量为$W$，

• 二次下降：求Hessian矩阵需要$\mathcal{O}(W^2)$的复杂度，目标优化中的$H^{-1}g$一项又至少需要$\mathcal{O}(W^3)$的复杂度
• 梯度下降：求梯度$\mathcal{O}(W)$复杂度，找极小值需要$\mathcal{O}(W)$的复杂度（待证明），所以一共只需要$\mathcal{O}(W^2)$的复杂度

P241 随机梯度下降相比于梯度下降的优势

1. 容易逃出局部最优解
2. 更高效。想象数据集复制成原先两倍，梯度下降要全部过一遍，而随机梯度下降则不受影响。

5.3 Error Backpropagation

P246 梯度的数值解近似

方法1为
$$\frac{\partial E_n}{\partial w_{ji}}=\frac{E_n(w_{ji}+ϵ)-E_n(w_{ji})}{ϵ}+\mathcal{O}(ϵ)$$
方法2为
$$\frac{\partial E_n}{\partial w_{ji}}=\frac{E_n(w_{ji}+ϵ)-E_n(w_{ji}-ϵ)}{2ϵ}+\mathcal{O}({ϵ}^2)$$
方法2的无穷小项阶数更小，可以用泰勒展开进行证明。这好神奇啊。
不过方法2的计算量是方法1的两倍（在神经网络已知$E_n(w_{ji})$的前提下）

P247 Jacobian矩阵

神经网络中可以计算jacobian矩阵，计算的时候可以前向计算，也可以反向计算。

5.4 The Hessian Matrix

计算Hessian矩阵时，通常考虑所有参数，包括weights和bias，一起算一个大的Hessian矩阵$\text{H}$.

P250 对角近似

• 考虑到有时候只是用Hessian矩阵的逆，所以更愿意估计一个对角阵$\text{H}$。存在一种方法能让计算复杂度退化到$\mathcal{O}(W)$（Becker and Le Cun, 1989; Le Cun et al., 1990），具体计算过程略去，翻书……
• 如果损失函数$E={\sum }_n{E}_n$，是$n$个样本相加，那么$\text{H}$也可以每个样本单独算，最后加起来

P251 外积近似

考虑最小二乘为损失函数的回归问题，
$$E=\frac{1}{2}\sum _{n=1}^N(y_n-t_n{)}^2$$
$$\text{H}={\nabla }^{2}E=\sum _{n=1}^N\nabla y_n^T\nabla y_n+\sum _{n=1}^N(y_n-t_n){\nabla }^2{y}_n$$
如果$y_n$和$t_n$很接近，第二项很小，忽略。（或者假定第二项中$y_n-t_n$与 ${\nabla }^2{y}_n$无关，则求和后，因为$y_n-t_n$的误差期望是0，苏所以也能忽略）
剩下
$$\text{H}\simeq \sum _{n=1}^N{\text{b}}_n{\text{b}}_n^T$$
其中${\text{b}}_n=\nabla a_n$，$a_n$是没有激活的网络最后一层（logit值）
对于逻辑回归，则有类似结论
$$\text{H}\simeq \sum _{n=1}^N{y}_n(1-y_n){\text{b}}_n{\text{b}}_n^T$$
这种方式的计算，比较高效，$\mathcal{O}(W^2)$的复杂度发生在${\text{b}}_n$的矩阵乘法那里

用这种方式可以近似计算${\text{H}}^{-1}$。考虑数据是序列进来的，
$${\text{H}}_{L+1}={\text{H}}_L+{\text{b}}_{L+1}{\text{b}}_{L+1}^T$$
利用woodbury等式可以得到
$${\text{H}}_{L+1}^{-1}={\text{H}}_L^{-1}-\frac{{\text{H}}_L^{-1}{\text{b}}_{L+1}{\text{b}}_{L+1}^T{\text{H}}_L^{-1}}{1+{\text{b}}_{L+1}^T{\text{H}}_L^{-1}{\text{b}}_{L+1}}$$
初始的${\text{H}}_0$可以选$\alpha \text{I}$，其中$\alpha$比较小。这样计算的是$\text{H}+\alpha \text{I}$的逆，但结果并不会对$\alpha$过于敏感

P252 有限差

数值方法计算
∂ 2 E ∂ w j i ∂ w l k = 1 4 ϵ 2 { E ( w j i + ϵ , w l k + ϵ ) − E ( w j i + ϵ , w l k − ϵ ) − E ( w j i − ϵ , w l k + ϵ ) + E ( w j i − ϵ , w l k − ϵ ) } + O ( ϵ 2 ) \frac{\partial^2 E}{\partial w_{ji}\partial w_{lk}}=\frac{1}{4\epsilon^2}\{ E(w_{ji}+\epsilon, w_{lk}+\epsilon) - E(w_{ji}+\epsilon, w_{lk}-\epsilon) -E(w_{ji}-\epsilon, w_{lk}+\epsilon) +E(w_{ji}-\epsilon, w_{lk}-\epsilon) \} + \mathcal{O}(\epsilon^2) ∂wji​∂wlk​∂2E​=4ϵ21​{E(wji​+ϵ,wlk​+ϵ)−E(wji​+ϵ,wlk​−ϵ)−E(wji​−ϵ,wlk​+ϵ)+E(wji​−ϵ,wlk​−ϵ)}+O(ϵ2)
这里的无穷小项是$\mathcal{O}({ϵ}^2)$，而不是$\mathcal{O}(ϵ)$.
这样的复杂度是$\mathcal{O}(W^3)$，太高了。一种有效的方法是先求梯度，再用数值方法
$$\frac{{\partial }^{2}E}{\partial w_{ji}\partial w_{lk}}=\frac{1}{2ϵ}\left\{\frac{\partial E}{\partial w_{ji}}(w_{lk}+ϵ)-\frac{\partial E}{\partial w_{ji}}(w_{lk}-ϵ)\right\}+\mathcal{O}({ϵ}^2)$$
复杂度为$\mathcal{O}(W^2)$

P253 精确计算

Hessian可以精确计算，不难，但是繁琐，要细心，想清楚计算图的依赖关系。具体公式翻书

P254 Hessian矩阵的快乘法

很多时候，我们需要的是Hessian矩阵的有关乘法，而不是其本身
如果要计算${\text{v}}^T\text{H}={\text{v}}^T\nabla (\nabla E)$，记 R { ⋅ } \mathcal R\{\cdot\} R{⋅}表示算子${\text{v}}^T\nabla$
这样可以对这个算子采用通常的微积分规则，并且已知
R { w } = v \mathcal R\{\textbf w\} = \textbf v R{w}=v
直接前向-反向传播就行了！真的神奇。（细节翻书）

5.5 Regularization in Neural Networks

P257 L2正则化的尺度

• 如果缩放输入和输出尺度，那么可以通过改变参数的尺度得到同样的映射。但是在有L2约束时，缩放后的映射不会是和原先的映射一样。所以需要给每一层都采用不同的L2权重系数。
  当尺度变化时，L2权重系数也要改变。
  以两层参数为例，先验通常为
  $$p(\text{w}|{\alpha }_1,{\alpha }_2)\propto \exp \left(-\frac{{\alpha }_1}{2}\sum _{w\in {\mathcal{W}}_1}{w}^2-\frac{{\alpha }_2}{2}\sum _{w\in {\mathcal{W}}_2}{w}^2\right)$$

• 注意这是一个反常先验improper. 不能标准化。因为偏置项$b$没有约束。
  反常先验会导致在贝叶斯框架下进行模型比较困难。因为模型证据为0。常用的做法是在偏置项上单独施加一个先验（例如高斯，虽然这样破坏了平移不变性）

P259 早停

早停有两种解释

1. 模型在训练过程中，有效自由度的数量会越来越高，导致模型越来越复杂。早停是一种限制模型复杂度的做法
2. 早停和L2权重衰减有相似性，且可以定量证明两者之间的联系，即$\tau \eta$（$\tau$是迭代次数，$\eta$是学习率）和正则化权重系数$\lambda$倒数有联系。这一点，Goodfellow的Deep Learning中也有写。

P261 不变性

如果希望神经网络对于输入的某种变换具有不变性（例如尺度不变性、平移不变性等）。
通常有四种方法：

1. 多采带变换的样本
2. 加入正则化，对变换进行惩罚（也即正切传播的作用）
3. 通过预处理，消除变化
4. 把不变性的性质整合到神经网络，如CNN

P264 正切传播

为了保持不变性，考虑变换为
$$\text{s}(\text{x},\xi )$$
其中$\text{s}(\text{x},0)=\text{x}$
则在$\text{x}$的流形$\mathcal{M}$上，样本移动方向为
$${\tau }_n={\frac{\partial \text{s}(\text{x},\xi )}{\partial \xi }|}_{\xi =0}$$
神经网络需要加入正则项
$$\Omega =\frac{1}{2}\sum _n\sum _k{\left({\frac{\partial y_{nk}}{\partial \xi }|}_{\xi =0}\right)}^2=\frac{1}{2}\sum _n\sum _k{\left(\sum _{i=1}^D\frac{\partial y_{nk}}{\partial x_{ni}}\frac{\partial x_{ni}}{\partial \xi }\right)}^2=\frac{1}{2}\sum _n\sum _k{\left(\sum _{i=1}^D{J}_{nki}{\tau }_{ni}\right)}^2$$
其中$J_{nki}$是Jacobian矩阵里的元素

• 如果有多个变换，对每个变换都具有不变性，那么对变换的组合就会有（局部）不变性

P265 Tikhonov正则化

$$\Omega =\frac{1}{2}\int \Vert \nabla y(\text{x}){\Vert }^{2}p(\text{x})d\text{x}$$
这个正则化项对输入的高斯噪声有一定鲁棒性，可以认为这是一种不变性。
实际上，这个正则化项可以通过误差函数
E ~ = 1 2 ∭ { y ( s ( x , ξ ) ) − t } 2 p ( t ∣ x ) p ( x ) p ( ξ ) d x d t d ξ \tilde E = \frac{1}{2}\iiint \{ y(\textbf s(\textbf x,\xi))-t \}^2p(t|\textbf x)p(\textbf x)p(\xi)d\textbf xdtd\xi E~=21​∭{y(s(x,ξ))−t}2p(t∣x)p(x)p(ξ)dxdtdξ
推导出，其中对于高斯噪声，$\text{s}(\text{x},\xi )=\text{x}+\xi$，$p(\xi )\sim \mathcal{N}(0,\lambda )$
可以推出
$$\tilde{E}=E+\lambda \Omega$$
具体过程略去

P269 软权值共享

为了减少参数，可以认为参数是从一堆$M$个混合高斯中来的，也即对每个参数$w_i$，有
$$p(w_i)=\sum _{j=1}^M{\pi }_j\mathcal{N}(w_i|{\mu }_j,{\sigma }_j^2)$$
正则化的先验为
$$\Omega (\text{w})=-\sum _i\ln \left(\sum _{j=1}^M{\pi }_j\mathcal{N}(w_i|{\mu }_j,{\sigma }_j^2)\right)$$
总损失函数为
$$\tilde{E}=E+\lambda \Omega$$
这里$\pi ,\mu ,\sigma$需要和参数联合优化，所以不能用EM算法。仍然采用梯度下降。先令
$${\gamma }_j(w)=\frac{{\pi }_j\mathcal{N}(w|{\mu }_j,{\sigma }_j^2)}{{\sum }_k{\pi }_k\mathcal{N}(w|{\mu }_k,{\sigma }_k^2)}$$
可以求得
$$\frac{\partial \tilde{E}}{\partial w_i}=\frac{\partial E}{\partial w_i}+\lambda \sum _j{\gamma }_j(w_i)\frac{(w_i-{\mu }_j)}{{\sigma }_j^2}$$
其直观含义是让$w$尽量靠近每个类的均值
$$\frac{\partial \tilde{E}}{\partial {\mu }_j}=\lambda \sum _j{\gamma }_j(w_i)\frac{({\mu }_j-w_i)}{{\sigma }_j^2}$$
其直观含义是让每个均值尽量靠近类中心
$$\frac{\partial \tilde{E}}{\partial {\sigma }_j}=\lambda \sum _i{\gamma }_j(w_i)\left(\frac{1}{{\sigma }_j}-\frac{(w_i-{\mu }_j{)}^2}{{\sigma }_j^3}\right)$$
为了防止$\sigma$为负，或为0，得到病态解，常用${\sigma }_j^2=\exp ({\zeta }_j)$来重参数化
对于$\pi$，因为限制了和为1，所以用$\eta$进行softmax重参数化，
$${\pi }_j=\frac{\exp ({\eta }_j)}{{\sum }_{k=1}^M\exp ({\eta }_k)}$$
求得
∂ E ~ ∂ η j = ∑ i { π j − γ j ( w i ) } \frac{\partial \tilde E}{\partial \eta_j}=\sum_i \{ \pi_j - \gamma_j(w_i) \} ∂ηj​∂E~​=i∑​{πj​−γj​(wi​)}
可以看到对$\eta$的梯度希望它和参数表现上的$\gamma$一致

5.6 Mixture Density Networks

如果输出是有几个分布混合的，例如在回归问题中，并不是单纯的高斯分布（比如输入房屋面积，输出房价，但是样本有北京和武汉两个地方，这对应两个不同分布）。这时候，如果强行用高斯分布拟合，结果对两边都会很差。
书中的例子是$t=x+0.3\sin (2\pi x)+noise$，然后把$x,t$对调

这时可以建模成混合分布，$$p(\text{t}|\text{x})=\sum _{k=1}^K{\pi }_k(\text{x})\mathcal{N}(\text{t}|{\mu }_k(\text{x},{\sigma }_k^2(\text{x}))$$
这是一个异方差heteroscedastic的例子，因为噪声来自$\text{x}$
如果$t$的维度是$L$，打算用$K$个混合高斯分布建模，此时神经网络的输出设计为：

1. $K$个分布先验概率${\pi }_k(\text{x})$，可以用无约束的${\eta }_k(\text{x})$代替，并做softmax重参数化
2. $K$个分布各自的方差${\sigma }_k^2(\text{x})$，可以用无约束的${\zeta }_k(\text{x})$代替，并用指数函数$\sigma =\exp (\zeta )$重参数化
3. $K$个类各自的均值${\mu }_k$，注意一共有$K\times L$个参数
   综上，一共需要$(K+2)L$个参数

对应损失函数为

关于${\mathbf{\mu }}_k,{\eta }_k,{\zeta }_k$的偏导数和软权值共享中的形式非常像，这里不抄了

训练完成后，推断时，可以求得均值
$$\mathbb{E}[\text{t}|\text{x}]=\sum _{k=1}^K{\pi }_k(\text{x}){\mu }_k(\text{x})$$
方差
$$s^2(\text{x})=\mathbb{E}[\Vert \text{t}-\mathbb{E}[\text{t}|\text{x}]{\Vert }^2|\text{x}]=\sum _{k=1}^K{\pi }_k(\text{x})\left\{{\sigma }_k^2(\text{x})+{\Vert {\mu }_k(\text{x})-\sum _{l=1}^K{\pi }_l(\text{x}){\mu }_l(\text{x})\Vert }^2\right\}$$

这个求混合分布的想法挺有意思的，注意推断时，均值$\mathbb{E}[\text{t}|\text{x}]$可能并不能很好地反映分布，更好的办法是找条件众数conditional mode，也即概率最大的$p(\text{t}|\text{x})$对应的$\text{t}$，但是这没有解析解，只能通过数值方法优化。一种简单的代替方法是直接找先验概率最大的部分，如上图(d)所示

5.7 Bayesian Neural Networks

贝叶斯方法中，为了进行预测，需要对参数进行积分。这在神经网络中非常困难。变分推断在假定后验是高斯的情况下进行。不过最完整的方法是基于拉普拉斯近似。
这里假设

1. 参数后验以众数为中心的高斯
2. 协方差很小
   （其实就是把后验建模成高斯，同时解决输出和输入之间不是线性关系的问题，当协方差很小时，后者就可以用一阶泰勒展开近似。这样问题就和前两章的回归和分类很像了）

P278 后验参数分布与贝叶斯回归网络

$$\begin{array}{rl}p(\text{w}|\alpha )& =\mathcal{N}(\text{w}|\text{0},{\alpha }^{-1}\text{I})\\ p(t|\text{x},\text{w},\beta )& =\mathcal{N}(t|y(\text{x},\text{w}),{\beta }^{-1})\end{array}$$
对于数据集$\mathcal{D}$，
$$p(\text{w}|\mathcal{D},\alpha ,\beta )\propto p(\text{w}|\alpha )p(\mathcal{D}|\text{w},\beta )=p(\text{w}|\alpha )\prod _{n=1}^N\mathcal{N}(t_n|y({\text{x}}_n,\text{w}),{\beta }^{-1})$$
因为$y({\text{x}}_n,\text{w})$非线性，后验不是高斯
但可以采用拉普拉斯近似找高斯形式的近似后验
最大化
ln ⁡ p ( w ∣ D ) = − α 2 w T w − β 2 ∑ n = 1 N { y ( x n , w ) − t n } 2 + const \ln p(\textbf w|\mathcal D)=-\frac{\alpha}{2} \textbf w^T\textbf w-\frac{\beta}{2}\sum_{n=1}^N\{ y(\textbf x_n,\textbf w)-t_n \}^2+ \text{const} lnp(w∣D)=−2α​wTw−2β​n=1∑N​{y(xn​,w)−tn​}2+const
中心为概率最高的参数点，即${\text{w}}_{MAP}$
该点二阶导为
$$\text{A}=-{\nabla }^2\ln p(\text{w}|\mathcal{D},\alpha ,\beta )=\alpha \text{I}+\beta \text{H}$$
其中$\text{H}$是平方和误差关于$\text{w}$的Hessian矩阵
所以，近似的后验为
$$q(\text{w}|\mathcal{D})=\mathcal{N}(\text{w}|{\text{w}}_{MAP},{\text{A}}^{-1})$$
预测分布为$$p(t|\text{x},\mathcal{D})=\int p(t|\text{x},\text{w})q(\text{w}|\mathcal{D})d\text{w}$$
该积分仍然无法理论处理，因为$y$是
$\text{w}$的非线性函数（也就是说这不是线性高斯模型）。但是如果认为$q$方差比较小，那么可以用泰勒展开得到
$$y(\text{x},\text{w})\simeq y(\text{x},{\text{w}}_{MAP})+{\text{g}}^T(\text{w}-{\text{w}}_{MAP})$$
其中$g={\nabla }_{\text{w}}y(\text{x},\text{w}){|}_{\text{w}={\text{w}}_{MAP}}$，这就变成了一个线性高斯模型
结果为$$p(t|\text{x},\mathcal{D},\alpha ,\beta )=\mathcal{N}(t|y(\text{x},{\text{w}}_{MAP}),{\sigma }^2(\text{x}))$$
其中${\sigma }^2(\text{x})={\beta }^{-1}+{\text{g}}^T{\text{A}}^{-1}\text{g}$

P280 超参数优化

可以对$\alpha ,\beta$进行优化，方法类似第三章的最大化证据函数，这里略去。
不过和第三章对比，不同的是这里$\alpha$的改变，会引起$\text{H}$的改变，进而影响$\lambda$。这里把这个影响忽略掉。
注意这里，可以根据$p(\mathcal{D})$去选模型，也即调整隐层中节点的个数

P281 贝叶斯分类网络

做法与回归类似，与贝叶斯逻辑回归也类似，还是要用拉普拉斯近似。
注意在优化超参数$\alpha$时，要计算计算$\ln p(\mathcal{D}|\alpha )$，方法是要用模型证据中的归一化系数估计方法（详见第4章）

感觉这样选出来的$\alpha$确实能比较好的防止过拟合

参考文献：
[1] Christopher M. Bishop. Pattern Recognition and Machine Learning. 2006