数据分析模型 第八章

一. 模型选择标准(Model selection criteria)

我们回顾下监督学习
当我们有$p+1$个变量和n组数据(n行)，如果我们要预测其中的一个变量，我们要利用剩下的部分变量(<=p个)来预测我们的这个变量。

如果我们的预测变量是分类数据，那么我们利用分类有关的模型。

如果我们的预测变量是数字型数据，那么我们利用回归相关的模型。

在第7章和第6章，我们有提到过，剩下的p个变量我们不一定全部都会用上，因为有的变量很重要，有的变量不重要。举个例子，我们预测房价，但是，假如我们有个变量是关于客户的血压，这个变量很明显和我们的预测房价这个变量多多少少关系不大，那我们则会忽略这个变量，如果我们有个变量是关于客户的薪资，这个变量我们自然会把它考虑在内。至于如何挑选有关变量请看第7章逻辑回归模型(分类模型)的变量选择，模型选择，和第6章的模型选择，利用模型的信息准则和正向逻辑或者逆向逻辑来选择相关的变量。

如果我们忽视了某个重要的变量，那么我们的模型很容非饱和，也就是我们的预测(对于训练数据内部的预测)会有很大的偏差。

如果我们把全部的变量都考虑在内，那么我们的模型很容易过饱和。这是因为我们的模型也“学习”到了数据里的噪音和随机变化。那么对于我们未知数据的预测(对于测试数据的预测)会产生极其不准确的影响。

非饱和，过饱和例子
假如我们这里有n组数据(n行)，我们有一个x变量和一个预测y变量。我们现在要用一个模型来诠释我们的数据，即x变量和y变量的关系。根据第6章所述，我们用一个简单粗暴的模型即多项式回归模型，它很简单，很灵活，可以表达非线性的关系。
我们先随便写个多项式模型：
$$y={\beta }_0+{\beta }_{1}x+{\beta }_2{x}^2+{\beta }_3{x}^3+...+{\beta }_p{x}^p$$
上述多项式模型有p+1项，x最高为p次，没错$x^2$意味平方，$x^p$意味着p次方。小弟这里只是为了举个非饱和和过饱和的例子，我们可以把$(x,x^2,x^3,...,x^p)$看成有p个x变量即$(x_1,x_2,...,x_p)$。

假设我们有50个数据，那么先把它的正确曲线标画出来，其实这50个数据就是由该曲线加上误差$\epsilon$产生的,如下图:

上图的黄色曲线为我们的正确回归曲线。

如果我们的多项式回归模型仅仅有两项即$(x,x^2)$,那么我们的图像虚线为:

从这张图像中我们可以看到，我们的虚线过于简单了，我们的虚线无法诠释大部分蓝色数据的点。这即为非饱和。

如果我们的多项式回归模型有20项$(x,x^2,..,x^{2}0)$,那么我们的图像虚线为:

从上图可以看到，我们的多项式回归模型过于复杂了，导致我们的虚线甚至都学习到了蓝色数据点中的噪音，变得更加弯弯曲曲。虽然对于训练数据来说，会有很好的表现即偏差较小，但会有很大的方差即对于预测测试数据的表现很差。这即为过饱和。

如果我们的多项式回归模型有5项$(x,x^2,.,x^5)$,那么我们的图像虚线为:

上述这张图，很明显就自然很多了，这样的多项回归模型所展示的虚线和我们的真实曲线就接近很多，由于测量的随机性和噪音导致这两条曲线无法完全一样，但足够接近。

均方预测误差(Mean squared prediction error)
我们回顾下之前的知识点:
如果是回归模型,那么:
$$E[Y|x_1,...,x_p]=f(x_1,...,x_p)$$
如果是二元分类模型,那么：
$$log(\frac{P(Y=1|x_1,...,x_p)}{1-P(Y=1|x_1,..,x_p)})=f(x_1,...,x_p)$$
$f(\cdot )$定义某个方程，即为真实的关于$y,(x_1,x_2,...,x_p)$关系方程。我们很明显是不知道这个$f(\cdot )$。我们会从总体数据中抽取样本这里我们写作$y=(y_1,y_2,...,y_n)$作为我们的训练数据。我们利用我们的训练数据去估计$f(\cdot )$得到我们的模型写作$\bar{M}(y)\equiv M$。之后，我们再从总体数据中取一组新的数据$y^{\prime }=(y_1^{\prime },y_2^{\prime },...,y_m^{\prime })$作为我们的测试数据。那么，我们可以利用我们之前得到的模型$\bar{M}(y)\equiv M$来做预测，得:
$${\bar{y}}_i(\bar{M}(y))=\bar{y}(x_{i,1}^{\prime },x_{i,2}^{\prime },...,x_{i,p}^{\prime },\bar{M}(y))$$
别被这公式的写法吓住了，这其实就是在说利用测试数据的$(x_{i,1}^{\prime },x_{i,2}^{\prime },...,x_{i,p}^{\prime })$和我们通过训练数据估计出的模型来预测我们的测试数据的y值($y_i^{\prime }$).

还记得我们第三章内容估计的均方误差(mean squared error,简写MSE)么，没错，一样的操作，那么我们可以写成:
$$MSPE(\bar{M}(y))=\frac{1}{m}\sum _{i=1}^m(y_i^{\prime }-{\bar{y}}_i(\bar{M}(y)){)}^2$$
这就是基于我们测试数据的均方预测误差. 那么MSPE越小，说明我们的预测越准确，因此我们需要找一个模型具有较小的MSPE

偏差和方差
根据我们的MSPE公式，我们会发现，MSPE这个公式很依赖于我们的估计的模型，而我们的估计的模型很依赖于我们的训练样本。还记得我们第二章和第三章讲置信区间时所说，我们一般会多次取样，每次取样本数为n，也就是说，MSPE很依赖于我们抽取一次样本该样本的具体数值，那么，为了解决这个我们MSPE公式对于抽一次样后样本具体数值的依赖，我们利用期望即均值对MSPE公式进行优化，那么得EMSPE（expected MSPE）:
$$EMSPE=E[MSPE(\bar{M}(y))]=偏{差}^2+方差$$
这里的期望意味着，多次取样后的均值。也就是，多次取样估计模型，然后分别计算MSPE，最后取均值得MSPE。这里的偏差意味着，我的预测与测试数据的偏差。方差为关于测试数据预测的变化程度大小。

那么我们会发现，如果要得到较小的EMSPE，我们需要使偏差和方差都要比较小一些。

再次回到我们的非饱和，过饱和，利用多项式回归模型的例子中去。如果我们多次取样，根据我们第二章和第三章的知识，我们可以得到95%的置信区间关于我们对于训练数据的预测。

当我们的多项式回归模型仅仅有两项即$(x,x^2)$时：

虚线为95%的置信区间预测，中间红色虚线为我们根据多次取样后分别估计参数得到的模型，进而分别预测出的值取均值。绿色的实线为真实的曲线。那么上图的 EMSPE=3.33, 偏差的平方为3.27，方差为0.058. 根据我们之前的结论我们知道这是非饱和的。

当我们的多项式回归模型有20项$(x,x^2,..,x^{2}0)$时：

EMSPE=0.4860,偏差的平方为0，方差为0.4860.我们可以看到在曲线两端的95%区间很明显范围变大了，这说明，我们的预测数值会有较大的变化。

当我们的多项式回归模型有5项$(x,x^2,.,x^5)$时：

EMSPE=0.147,方差为0.147，偏差的平方为0
我们可以看到在曲线两端的95%区间很明显范围还是比较窄的，这说明，我们的预测数值不会有较大的变化。

所以，总的来说:
1.模型越简单(变量越小)会有较大的偏差，较低的方差。
2.模型越复杂(变量越多)会有较低的偏差，较高的方差。
3.方差会随着样本数目n的增加而减少，但偏差不会。
4.方差会随着变量数目的增加而增加。
5.忽视一个变量会减少方差，但如果那个变量很重要则会增加偏差。

所以对于模型的变量选择就需要我们对于偏差和方差之间做出衡量，从而选择一个适当的模型复杂度。

假设检验和多次检验问题(hypothesis testing and multiple testing problems)
我们再次回到了小弟经常说的一个问题，如何判断该变量是否重要。我们会用到一种方法：利用假设检验
$$H_0:{\beta }_j=0$$
$$H_1:{\beta }_j\ne 0$$
如果p值较小，那么我们拒绝原假设。

但是，在运用中，这个方法有个潜在问题，当我们仅仅测试一个变量，如果它的p值小于0.05我们拒绝原假设$\beta =0$，也就是说该变量与我们的预测y值无关。因为我们仅仅测试了这一个变量，那么我们也可以认为5%的概率该变量与我们的y值预测有关。那么问题来了，如果我们有p个变量，那么我们就会做p次检验，从概率的角度上来看那么会有px5%个变量会和y值的预测有关，即有px5%个变量是重要的。如果我们的p很大，假如我们的p=10000，那么会有10000*5%=500个变量很重要。即使我们这10000个变量其实跟y值无关，但我们从概率的角度来看依然有500个很多的变量跟y值有关。也就是说，其实并没有那么多的变量跟y值有关。这个问题被称为数据疏浚(data dredging)意味着当我们在数据中寻找更多有意义的信息的时候，会导致本来无关系的变量之间变的有关系甚至变的关系重要。但单纯由于p值导致的数据疏浚我们也可以把该问题称为"p-hacking".但归根结底我们也知道它本身是因为大量多次检验导致的。

关于这个问题，有很多方法可以解决它。最常用的方法是邦费罗尼矫正(Bonferroni Procedure).关于它的由来和更深的原理，小弟也不知道，我们这里只学习他如何解决由于p值导致的数据疏浚问题。

也很简单，举个例子就明白了，如果我们要做p次检验。那么邦费罗尼这个大哥告诉我们，只有当：
$$P值<\frac{\alpha }{p}$$
时，我们才要拒绝原假设。这里的$\alpha =0.05$。

似然和信息准则
如何判断该变量是否重要。我们会用到另外一种方法：似然和信息准则
小弟在第6章已经简单分享了该方面知识，现在带大家重新回顾一下。
我们经常用负log似然来估计模型的参数，从而使得我们的负log似然的值为最小值：
$$L(y|\bar{\theta })=-\sum _{i=1}^{n}log(p(y_i|\bar{\theta }))$$
然后在给负log似然加上不同的惩罚方程形成我们的信息准则，我们可以通过选取不同的变量搭配使我们的信息准则为最小值，这样该不同的变量搭配即为我们认为的重要变量。

我们先了解下写法:
$M$意味着模型，$\bar{M}$意味着通过数据估计出模型参数后的模型，即该模型可以诠释我们的数据。
$L(y|\bar{M})$ 为最小负log似然值，也就是我们把我们估计后参数带入负log似然公式中而已。
$k_M$意味着变量的个数

那么有这么几种常用的信息准则

1.赤池信息准则(AIC)
$$AIC(\bar{M})=L(y|\bar{M})+k_M$$

2.库尔贝克信息准则(Kullback-Leibler information criterion. KIC)
库尔贝克这个大哥跟费雪大哥的名气是一个级别的，我们后续会学更多关于他的理论知识。
$$KIC(\bar{M})=L(y|\bar{M})+\frac{3}{2}{k}_M$$
当n很大时，AIC和KIC都很容易使我们最后的模型过饱和。但很不幸，到目前为止，没有一个特别严谨的方法可以中和过饱和，当我们训练的数据n过多使我们的参数更加精准，很容易过饱和，减少偏差。但是在运用时，KIC比AIC 表现更好点，因为KIC 可以减少一些过饱和的程度。

3.贝叶斯信息准则(BIC)
$$BIC(\bar{M})=L(y|\bar{M})+\frac{k_M}{2}logn$$
当我们的训练数据n增加，BIC 过饱和的概率会趋近于0.相反，当n很小，BIC会使得我们的模型非过饱和。所以我们一般会利用AIC 和BIC 来中和我们的模型，一会算算看AIC 一会算算看BIC，这样时我们的模型既不过饱和，也不非过饱和，当我们算的AIC和BIC之间差>=3，那么我们认为这个模型正好。

4.风险膨胀信息准则(Risk inflation information criterion. RIC)
$$RIC(\bar{M})=L(y|\bar{M})+k_{M}logp$$
$p$意味着我们认为可能的变量个数。RIC这个方法本质和多次检验方法类似。当我们常识过多的变量个数时，我们要尽量避免过饱和。当p值很大，通常我们会用RIC.

我们在利用我们之前多项式模型的例子，来画张图，如下:

该图的y轴为我们的信息准则数值，x轴代表多项式的项数。从该图中我们会发现，如果用AIC的话，我们的多项式模型的项达到7项，相比其他信息准则，更加过饱和些。

交叉验证(Cross validation. CV)
交叉验证是一个很普遍的方法，它可以控制我们模型的复杂度。

根据小弟在本章最早提到，我们选择我们的模型通过减小预测值和测试数据里测试数据里的真值($y_1^{\prime },y_2^{\prime },..,y_m^{\prime }$)进行比较，从而选择适合的模型。
$$MSPE(\bar{M}(y))=\frac{1}{m}\sum _{i=1}^m(y_i^{\prime }-{\bar{y}}_i(\bar{M}(y)){)}^2$$
但计算MSPE会有一个小小的问题，那就抽取测试数据，小弟在本章讲MSPE的时候说的很轻松抽取一组样本数据作为测试数据。但这组测试数据不仅要提前准备好，而且和训练数据不能有相同。不仅如此，训练数据和测试数据的数据分布也要相同。这样的话如果我们从总体中要抽取训练数据，那么无疑困难是很大的。

在交叉验证这个方法中，我们可以抽取样本，并将样本分为两部分，第一部分为训练数据，用来估计出参数，得到我们的模型。第二部分为测试数据用来看我们的模型对于未来数据预测的表现，从而选取适当的模型。

交叉验证法逻辑:
1.取样数据y，并分为两个部分：
-训练集$y_{train}$
-测验集$y_{test}$
2.利用$y_{train}$估计出参数，得到我们的模型$M$
3.计算$MSPE$,根据$y_{test}$
4.反复这个流程多次(尽可能的多)，计算$MSPE$均值，即$EMSPE$,只不过这里我们不爱叫它$EMSPE$而叫它“CV error”。那么我们变可以知道我们的模型对于未来预测数据的表现情况。CV error越小，模型对于预测未来数据的表现越好。另外注意点是，在第4步中，我们尽可能多来几次，计算MSPE均值，为了避免随机抽样后导致MSPE的变化，从而使我们的CV error变的更低。

另外这里有两种交叉验证法CV的变种:
第一种: K-fold CV (KF CV)
将样本均分为K份，训练集为K-1份，剩下的一份作为测试集。然后一样反复流程，最后求均值MSPE。

第二种: Leave-one-out CV(LOO CV)
如果我们的样本数据有n个，n-1个数据作为训练集，剩下的那一个数据作为测试集。然后一样反复流程，最后求均值MSPE。

看到这里，可能大家会发现一个问题，那就是小弟在第6章选择模型内容中所提到的，如果我们有p个变量，那么一共有$2^p$种组合，对没错，都尝试一次，根据多次假设检验也好，根据信息准则也罢，根据CV也行，这些计算量无论哪一种方法大家可想而知是很大的。令一个问题那就是，如果我们的训练样本的数据有少许变化，那么我们最终的结果会有很大不同，这也是这些传统选择模型方法的不足，即为统计的不稳定性(Statistical Instability)。

二. 惩罚回归(Penalized Regression)

到目前为止，我们学到了回归模型中的简单线性模型，和多元线性模型，以及多项式回归模型，和关于分类模型里的二元逻辑分类模型。现在小弟再简单的介绍下惩罚回归模型里的一些著名的模型。

惩罚回归其实是,针对于我们估计参数的惩罚，再将我们惩罚后的估计出的参数放入我们原有的回归模型中去罢了。对于估计参数的惩罚一般会写成:
( β ˉ 0 , β ˉ ) = a r g m i n β 0 , β { R S S ( β 0 , β ) + λ ∑ j = 1 p g ( β j ) } (\bar \beta_0,\bar \beta)=arg min_{\beta_0,\beta}\{RSS(\beta_0,\beta)+\lambda \sum_{j=1}^{p}g(\beta_j)\} (βˉ​0​,βˉ​)=argminβ0​,β​{RSS(β0​,β)+λj=1∑p​g(βj​)}

上述公式里的$RSS(\cdot )$为拟合优度，例如我们在线性模型里用的最小二乘法，在二元逻辑回归里用的似然，其实就是我们正常情况下估计参数的方法。然后加上惩罚方程$g(\cdot )$

现在我们来理解下这个惩罚方程，一般情况下$|{\beta }_j|$越大，我们对它的惩罚越大，为什么？因为如果一个参数的值越大那么它对我们的模型和预测会有很大的影响，如果参数值很小那么对于预测变量会有较小的影响。但惩罚方程它想做的就是鼓励较小的参数值，惩罚较大的参数值，这样做的话会使减少模型的过饱和几率。对没错，可能有的同学发现了，那么其实它无疑是增加了偏差，从而减小了方差，所以惩罚回归牺牲了偏差，从而增加了对于未来预测数据的稳定性。这里我们会注意到，惩罚一般是不会惩罚${\beta }_0$的，因为${\beta }_0$的数值其实由我们的预测变量y的单位所决定的，例如如果我们的y是米，那么${\beta }_0$可能以1来计算。如果y是厘米，那么${\beta }_0$可能以100来计算。

说到单位，现在我们出现了个很头疼的问题那就是变量的范围。例如有的变量数据在(0-100)的范围内，但有的数据在(0-1)范围内。那么我们需要先对数据进行标准化即数据标准化（data standardise），使变量满足:
$$\sum _{i=1}^n{x}_{i,j}=0$$
$$\sum _{i=1}^n{x}_{i,j}^2=n$$
使我们的变量数据均满足上述两个公式，其实也很容易做到，只要减去均值除以标准差即可。但为什么非要以该形式进行标准化呢？因为这样有利于我们之后对于数据的运算，例如估计参数，例如寻找重要的变量等等。

回归正题
( β ˉ 0 , β ˉ ) = a r g m i n β 0 , β { R S S ( β 0 , β ) + λ ∑ j = 1 p g ( β j ) } (\bar \beta_0,\bar \beta)=arg min_{\beta_0,\beta}\{RSS(\beta_0,\beta)+\lambda \sum_{j=1}^{p}g(\beta_j)\} (βˉ​0​,βˉ​)=argminβ0​,β​{RSS(β0​,β)+λj=1∑p​g(βj​)}
上述公式的$\lambda$即为惩罚的力度，它最为人所知的名称叫为超参数(hyperparameter).
例如我们的线性回归，如果$\lambda =0$，那我们的惩罚回归就是最小二乘法。惩罚回归有个好处就是在选择模型当中，它解决了统计的不稳定性。训练数据的小变化，使我们的最后结果也是很小的变化。

看到这里，想必大家对惩罚回归有了初步的了解，小弟在此简单的介绍两个用的最多的惩罚回归模型。

1.岭回归(Ridge regression)
( β ˉ 0 , β ˉ ) = a r g m i n β 0 , β { R S S ( β 0 , β ) + λ ∑ j = 1 p β j 2 } (\bar \beta_0,\bar \beta)=arg min_{\beta_0,\beta}\{RSS(\beta_0,\beta)+\lambda \sum_{j=1}^{p}\beta_j^2\} (βˉ​0​,βˉ​)=argminβ0​,β​{RSS(β0​,β)+λj=1∑p​βj2​}
岭回归的优点: 第一对于大量的样本数目n，和大量的变量，它都算的很快。第二当很多变量之间有很大关联时，它的方差会很低。但无论$\lambda$为多少时，即$\lambda <\infty$,对于参数的估计不会等于0，但会收敛于0

2.拉索回归(lasso regression)
( β ˉ 0 , β ˉ ) = a r g m i n β 0 , β { R S S ( β 0 , β ) + λ ∑ j = 1 p ∣ β j ∣ } (\bar \beta_0,\bar \beta)=arg min_{\beta_0,\beta}\{RSS(\beta_0,\beta)+\lambda \sum_{j=1}^{p}|\beta_j |\} (βˉ​0​,βˉ​)=argminβ0​,β​{RSS(β0​,β)+λj=1∑p​∣βj​∣}
相比岭回归，对于参数的估计是会等于0的。

那么如何选择适合的惩罚力度$\lambda$呢？
最常用的方法是利用CV来选择$\lambda$，根据预测的值和测验数据的真值差的多不多即CV error。
1.我们定我们的$\lambda$属于0(无惩罚)到一个很大的值(较大惩罚)
2.对于每个$\lambda$我们都尝试一次，利用CV来计算估计预测误差。
3.选择一个较小CV error的$\lambda$
4.利用该$\lambda$来估计我们的模型参数
如果用MATLAB的同学，我们会用到lasso()和lassoglm()这两个公式。如果是用R的同学，我们需要用glmnet

惩罚回归的偏差和方差
根据：
$$EMSPE=E[MSPE(\bar{M}(y))]=偏{差}^2+方差$$
我们的惩罚回归会降低方差，即提升对于未来数据预测的稳定性，但降低了偏差准确度。因为惩罚回归会将估计的参数逐渐向0收敛，惩罚力度越大，参数估计会趋向于0.对于本身参数值较小的，惩罚力度会小，对于参数值较大的，惩罚力度会大。根据这个特点估计出来的参数也被称为收缩估计(shrinkage estimate), 但拉索的惩罚可以将参数估计为0，所以根据拉索的惩罚回归估计出的参数也被称为稀疏估计(Sparse estimate).

最后来两张图大家就明白了

岭回归:

拉索回归:

x轴即为我们的$\lambda$惩罚力度，y轴为参数估计的具体数值。
无论是岭回归还是拉索回归，我们可以看到蓝线(y轴0下面的蓝线)几乎没有受到什么惩罚，因为它是${\beta }_0$，其次我们可以看到，随着惩罚力度的增加，拉索估计出的参数可以为0，但岭回归是逐渐向0收敛。另外我们也可以看到无论是哪种惩罚回归，参数的绝对值越大，惩罚力度越大，弧线更狠，但对于初始绝对值很小的参数随着惩罚力度增加，它受到的惩罚不是很大，弧线很平缓。

三. 结语

有公式推导错误的或理论有谬误的，请大家指出，我好及时更正，感谢。