偏差-方差均衡
由上一讲可知，普通最小二乘算出的权重\(w_{\rm LS}=(X^\mathsf{T}X)^{-1}X^\mathsf{T}y\)，为无偏差估计但是有比较高的方差，即 \[ \mathbb{E}[w_{\rm LS}] = w,\ \ {\rm var}[w_{\rm LS}] = \sigma^2 (X^\mathsf{T}X)^{-1} \]
而岭回归算出的权重\(w_{\rm RR} = (\lambda I + X^\mathsf{T}X)^{-1}X^\mathsf{T}y\)，用前面的方法，可以算出 \[ \mathbb{E}[w_{\rm RR}] = (\lambda I + X^\mathsf{T}X)^{-1}X^\mathsf{T}Xw,\ \ {\rm var}[w_{\rm RR}] = \sigma^2 Z(X^\mathsf{T}X)^{-1}Z^\mathsf{T} \\ Z = (I+\lambda(X^\mathsf{T}X)^{-1})^{-1} ( = (X^\mathsf{T}X+\lambda I)^{-1}X^\mathsf{T}X) \] 可以看出来是一个有偏估计，但是方差更小。那么哪种方法更好？实际上，我们关心所得到的权重的通用性，即\(w\)对于新的数据\((x_0, y_0)\)表现如何（这里我们只知道\(x_0\)，还不知道\(y_0\)）。我们可以知道使用普通最小二乘预估的\(\hat{y}_0 = x_0^\mathsf{T} w_{\rm LS}\)，使用岭回归预估的\(\hat{y}_0 = x_0^\mathsf{T}w_{\rm RR}\)
对于建模得出的\(\hat{w}\)（这里\(\hat{w}\)可以是\(w_{\rm LS}\)或者\(w_{\rm RR}\)），可以使用平方误差和期望的定义来计算误差的期望，即 \[ \mathbb{E}[(y_0 - x_0^\mathsf{T}\hat{w})^2|X, x_0] = \int_{\mathbb{R}}\int_{\mathbb{R}^n} (y_0 - x_0^\mathsf{T}\hat{w})^2 p(y|X,w) p(y_0|x_0, w)dydy_0 \]
也就是说，假设我知道\(X\)，\(x_0\)，同时假设有一个真实的\(w\)（尽管不知道具体是什么），而且结果\(y \sim N(Xw, \sigma^2I)\)，要使用\(\hat{w} = w_{\rm LS}\ {\rm or}\ w_{\rm RR}\)去近似\(w\)，然后用\(y_0 \approx x_0^\mathsf{T}\hat{w}\)去预测\(y_0 \sim N(x_0^\mathsf{T}w, \sigma^2)\)，计算预测的期望平方误差
假设是计算以\(x_0\)和\(X\)的条件误差，可以套用公式 \[ \mathbb{E}[(y_0 - x_0^\mathsf{T}\hat{w})^2] = \mathbb{E}[y_0^2] - 2\mathbb{E}[y_0]x_0^\mathsf{T}\mathbb{E}[\hat{w}] + x_0^\mathsf{T}\mathbb{E}[\hat{w}\hat{w}^\mathsf{T}]x_0 \]
由于\(y_0\)和\(\hat{w}\)是独立的，因此\(\mathbb{E}[y_0\hat{w}] = \mathbb{E}[y_0]\mathbb{E}[\hat{w}]\)。且由前面的推导结论，有 \[\begin{align*} \mathbb{E}[\hat{w}\hat{w}^\mathsf{T}] &= {\rm var}[\hat{w}] + \mathbb{E}[\hat{w}]\mathbb{E}[\hat{w}]^\mathsf{T} \\ \mathbb{E}[y_0^2] &= \sigma^2 + (x_0^\mathsf{T}w)^2 \end{align*}\]
代入上面的式子，有 \[\begin{align*} \mathbb{E}[(y_0 - x_0^\mathsf{T}\hat{w})^2] &= \sigma^2 + (x_0^\mathsf{T}w)^2 - 2(x_0^\mathsf{T}w)(x_0^\mathsf{T}\mathbb{E}[\hat{w}]) + (x_0^\mathsf{T}\mathbb{E}[\hat{w}])^2 + x_0^\mathsf{T}{\rm var}[\hat{w}]x_0 \\ &= \sigma^2 + x_0^\mathsf{T}(w-\mathbb{E}[\hat{w}])(w-\mathbb{E}[\hat{w}])^\mathsf{T}x_0 + x_0^\mathsf{T}{\rm var}[\hat{w}]x_0 \end{align*}\]
所得到的泛化误差可以分为三个部分
- \(\sigma^2\)：测量误差，无法控制
- 中间项：模型的偏差：得到的模型与真实模型之间的平均距离
- 末项：得到的模型对训练数据有多敏感
推而广之，可以认为真实模型都可以表示成\(y=f(x;w)+\epsilon\)的形式，其中\(\mathbb{E}(\epsilon) = 0, {\rm var}(\epsilon) = \sigma^2\)。我们通过最小化损失函数\(\hat{f} = {\rm arg}\min_f \mathcal{L}_f\)来近似\(f\)，然后将\(\hat{f}\)用于新的数据，即$y_0 (x_0) _0 $。可知
\[\begin{align*} \mathbb{E}[(y_0 - \hat{f}_0)^2] &= \mathbb{E}[y_0^2] - 2\mathbb{E}[y_0\hat{f}_0] + \mathbb{E}[\hat{f}_0^2] \\ &= \sigma^2 + f_0^2 - 2f_0\mathbb{E}[\hat{f}_0] + \mathbb{E}[\hat{f}_0]^2 + {\rm var}[\hat{f}_0] \\ &= \sigma^2 + (f_0 - \mathbb{E}[\hat{f}_0])^2 + {\rm var}[\hat{f}_0] \end{align*}\]
第一项是数据本身的噪声，第二项是偏差的平方，第三项是方差
选择参数一个很好的方法是交叉验证，过程大略为
- 将数据分成大致相等的\(K\)份
- 在\(K-1\)份数据上学习模型，在剩下的第\(K\)份上进行预测
- 重复\(K\)次，每次留出一组数据
- 用预测的平均情况评价性能
贝叶斯定律
考虑岭回归的目标函数 \[ \mathcal{L} = \sum_{i=1}^n (y_i - x_i^\mathsf{T}w)^2 + \lambda w^\mathsf{T}w \]
这里惩罚项对比较大的\(w\)进行了惩罚，其隐含的意思是，我们对“什么样的\(w\)是好的\(w\)”有一个先验知识。贝叶斯定律可以将这种知识形式化。我们有
- \(P(A, B) = P(A|B)P(B) = P(B|A)P(A)\)
- \(P(A) = \sum_bP(A, B=b)\)
- \(P(B) = \sum_aP(A=a, B)\)
综上所述，对离散随机变量有 \[ P(A|B) = \frac{P(B|A)P(A)}{P(B)} = \frac{P(B|A)P(A)}{\sum_a P(B|A=a)P(A=a)} \\ P(B|A) = \frac{P(A|B)P(B)}{P(A)} = \frac{P(A|B)P(B)}{\sum_b P(A|B=b)P(B=b)} \]
对于\(P(B|A) = P(A|B)P(B) / P(A)\)来讲，\(P(B|A)\)称为后验，\(P(A|B)\)称为似然，\(P(B)\)称为先验，\(P(A)\)称为边际（marginal）
离散的情况可以推广到连续的情况。假设\(\theta\)是模型参数，为连续值；\(X\)为要处理的数据，根据贝叶斯定律，有
\[ P(\theta|X) = \frac{P(X|\theta)P(\theta)}{\int P(X|\theta)P(\theta)d\theta} = \frac{P(X|\theta)P(\theta)}{P(X)} \]
这里\(P(X|\theta)\)是似然，可以从模型定义得到；\(P(\theta)\)是定义的先验分布
可以看一个掷硬币的例子。假设我们手中有一枚硬币，掷出正面的概率是\(\pi\)，现在掷了\(n\)次，得到一个序列（假设正面记为1反面记为0）。假设每次投掷是独立的，也就是说有 \[ p(x_1, \ldots , x_n|\pi) = \prod_{i=1}^np(x_i|\pi) = \prod_{i=1}^n \pi^{x_i}(1-\pi)^{1-x_i} \]
假设\(\pi\)的先验概率满足Beta分布 \[ p(\pi) = Beta(\pi|a,b) = \frac{\Gamma(a+b)}{\Gamma(a)\Gamma(b)}\pi^{a-1}(1-\pi)^{b-1} \]
如何求给定\(x_1, \ldots, x_n\)条件下\(\pi\)的后验分布？
可以先套入贝叶斯定律 \[ p(\pi | x_1, \ldots, x_n) = \frac{p(x_1, \ldots, x_n|\pi)p(\pi)}{\int_0^1 p(x_1, \ldots, x_n|\pi)p(\pi)d\pi} \]
由于分母所做的就是将分子归一化，与\(\pi\)无关，因此可以忽略，即可以写作\(p(\pi|x) \propto p(x|\pi)p(\pi)\)。代入有 \[\begin{align*} p(\pi|x_1, \ldots, x_n) &\propto \left[\prod_{i=1}^n \pi^{x_i}(1-\pi)^{1-x_i}\right]\left[\frac{\Gamma(a+b)}{\Gamma(a)\Gamma(b)}\pi^{a-1}(1-\pi)^{b-1}\right] \\ & \propto \pi^{\sum_{i=1}^n x_i + a - 1}(1-\pi)^{\sum_{i=1}^n (1-x_i) + b - 1} \end{align*}\]
即\(p(\pi|x_1, \ldots, x_n) = Beta(\sum_{i=1}^n x_i + a, \sum_{i=1}^n (1-x_i) + b)\)
最大后验（Maximum A Priori -- MAP）
之前最小二乘的\(w\)可以通过MLE得到，那么岭回归也可以使用概率解释。沿用之前的似然模型，\(y\sim N(Xw, \sigma^2I)\)，不同的是这里我们还要为\(w\)提出一个先验模型。假设这个模型满足高斯分布，\(w\sim N(0, \lambda^{-1}I)\)，那么 \[ p(w) = \left(\frac{\lambda}{2\pi}\right)^{\frac{d}{2}} {\rm e}^{-\frac{\lambda}{2}w^\mathsf{T}w} \]
可以用最大后验估计来求在后验分布下最可能的\(w\)，即 \[\begin{align*} w_{\rm MAP} &= {\rm arg}\max_w \ln p(w|y, X) \\ &= {\rm arg}\max_w \ln \frac{p(y|w, X)p(w)}{p(y|X)} \\ &= {\rm arg}\max_w \ln p(y|w, X) + \ln p(w) - \ln p(y|X) \end{align*}\]
这里第一项是原最小二乘法的目标函数（最大似然），第二项是先验，第三项可以不考虑，因为没有涉及到\(w\)，而且在很多模型里我们也不知道这个项会是什么。代入得
\[\begin{align*} w_{\rm MAP} &= {\rm arg}\max_w \ln p(y|w, X) + \ln p(w) \\ &= {\rm arg}\max_w -\frac{1}{2\sigma^2}(y-Xw)^\mathsf{T}(y-Xw) - \frac{\lambda}{2}w^\mathsf{T}w + C \end{align*}\]
令这个函数为目标函数\(\mathcal{L}\)，则可以通过令其对\(w\)的偏导数为0求解最优的\(w\) \[ \nabla_w \mathcal{L} = \frac{1}{\sigma^2} X^\mathsf{T}y - \frac{1}{\sigma^2}X^\mathsf{T}Xw - \lambda w = 0 \] 解得\(w_{\rm MAP} = (\lambda \sigma^2I +X^\mathsf{T}X)^{-1}X^\mathsf{T}y\)。因此可以看出 \[ w_{\rm MAP} = w_{\rm RR} \\ w_{\rm ML} = w_{\rm LS} \]
