五、线性回归建模及检验★

在第四章中，已经完成回归分析之线性回归模型构建。接下来是模型诊断，可以结合第四章线性回归一起学习。

1. 当Y服从正态分布时，模型表现指标会好。如$R^2$(精确度)、β(稳定性)。
2. 模型要求：
   精确度$R^2$变化不大且稳定性β长期保持在一个水平上。
3. 要保持模型精确和稳定，就要进行模型检验，使模型满足假设。

1 多元回归假设★★★

多元回归假设：

• Y的平均值能够准确地被由X组成的线性函数建模出来的。(回归均值)
• 解释变量X和随机扰动项$ϵ$不存在线性关系。(回归公式：$Y=\beta X+ϵ$)
• 解释变量X之间不存在线性关系(或强相关)
• 假设随机误差项$ϵ$是一个均值为0的正态分布
• 假设随机误差项$ϵ$的方差恒为${\sigma }^2$
• 误差是独立的

注：

线性关系：解释变量和被解释变量之间是线性关系

解释1和解释4：$ϵ$（不能被X所解释的其他因素）和Y的分布是完全一样的，只是Y的均值是βX，所以$ϵ$是均值为0 的正态分布。

解释2：X是拿到的数据对Y的解释。$ϵ$是没有拿到的数据对Y的解释。所以解释变量和扰动项之间不存在线性关系。
多元线性回归：$Y={\beta }_1{X}_1+{\beta }_2{X}_2+ϵ$【扰动项：影响Y但是没有被纳X中的】

例子：孩子的收入【Y】 = 孩子教育水平【X】 + 父亲职业
父亲职业不影响孩子教育水平
=>父亲职业纳入扰动项：孩子的收入【Y】 = 孩子教育水平【X】 + 父亲职业【扰动项】
父亲职业影响孩子教育水平
=>父亲职业纳入解释变量 : 孩子收入=$\beta_1$孩子教育+$\beta_2$父亲收入
=>若未纳入解释变量，则变为：孩子收入=(${\beta }_1+{\beta }_2$)孩子教育 +父亲职业【扰动项】，导致孩子教育对孩子收入的影响虚高。
所以解释变量和扰动项之间不存在线性关系，如果存在线性关系，就将跟X相关的进行多的变量纳入模型中。这也是多元回归存在的必要性。

解释3：
一元线性回归 : $Y={\beta }_1{X}_1+{\beta }_0+ϵ$
一元线性回归进行估计时，估计的是回归系数的均值(见第四章)和回归系数的标准差。
$$S_{\beta }=\frac{S_{ϵ}}{|X|}$$
$S_{\beta }$是回归系数的标准差，$S_{ϵ}$是扰动项(残差)的标准差。残差 $e=y-\hat{y}$。
|X|中如果X中两列线性相关，那么行列式为0，由于分母不能为0，所以解释变量(X)之间不存在线性相关。

解释5：
下图是Y的预测值和残差的关系图。

$\hat{y}=\beta x$,$\beta$是常量，所以$\hat{y}$与x的分布是一致的。由图可知，随着x的增大，扰动项的方差是不变的。

解释6：Y是客观存在的、独立的，产生的误差是独立的。

2 模型检验之残差检验

模型检验目的：检验模型和假设是否一致，若不一致调整到一致即可。

• 检验的根本：画图，对数据有一个感性认识。

• 残差检验标准：
  1.残差是否离群
  2.残差是否和某个X有曲线关系
  3.残差是否和某个X有关

2.1 查看残差图

• 残差独立同分布且服从正态分布。
• 图形横坐标是Y的预测值，纵坐标是扰动项(残差)。
• 点到线的竖直距离是残差。
• 翻转的正态分布图形是Y的预测值，$\hat{y}=\beta x$

正常残差图：随着X的增加，残差是随机分布的。
X为非线性：随着X的变化，残差呈现抛物线变化，需要将$X^2$放进模型。
异方差：随着X的增大，残差呈现扇面式的分布。常出现在横截面数据(样本在同一个时间采集到的)。
自相关：随着X的变化，残差呈现周期式波动。常出现在时间序列数据上。

2.2 解决方案★★

• X为非线性：加入X的高阶形式(如$X^2$)
• 异方差：横截面数据经常出现异方差现象，修正的方法如：加权最小二乘法、稳健回归，最简单的是对Y取自然对数。
• 自相关：分析时间序列和空间数据时易出现，复杂方法是使用时间序列或空间计量方法分析，简单方法是加入Y的一阶滞后项进行回归。

注：
在描述性统计分析中，连续变量只会有四种可能：正态分布、泊松分布、伽马分布、对数正态分布；其中对数正态分布又是右偏最严重的函数，对于对数正态分布来说，对函数取对数之后就是正态分布。

2.3 案例：线性回归诊断——模型检验

2.3.1 导入数据

import matplotlib.pyplot as plt
import os
import numpy as np
import pandas as pd
import statsmodels.api as sm
from statsmodels.formula.api import ols

os.chdir(r"E:\7linearmodel")
raw = pd.read_csv(r'creditcard_exp.csv', skipinitialspace=True)
raw.head()

- 数据基本情况
# 有支出数据的情况
exp = raw[raw['avg_exp'].notnull()].copy().iloc[:, 2:].drop('age2',axis=1)
# 没有支出数据的情况，数据待预测
exp_new = raw[raw['avg_exp'].isnull()].copy().iloc[:, 2:].drop('age2',axis=1)

exp.describe(include='all')

2.3.2 相关性分析

- 1.散点图：收入与支出关系
exp.plot('Income', 'avg_exp', kind='scatter')
plt.show()

- 2.pearson相关系数
exp[['Income', 'avg_exp', 'Age', 'dist_home_val']].corr(method='pearson')

2.3.3 简单线性回归

- 1.建立线性回归模型
lm_s = ols('avg_exp ~ Income', data=exp).fit()
lm_s.summary()

- 2.线性回归模型预测、残差值
lm_s.predict(exp)
lm_s.resid
pd.DataFrame([lm_s.predict(exp),lm_s.resid],index=['predict','resid']).T.head()

- 3.待预测数据集中得到预测值
lm_s.predict(exp_new)[:5]

2.3.4 多元线性回归

1. 建立多元线性回归模型

- 建立多元线性回归模型
lm_m = ols('avg_exp ~ Age + Income + dist_home_val + dist_avg_income',
           data=exp).fit()
lm_m.summary()

2. 模型筛选，找出最好模型【见四、线性回归模型构建】

2.3.5 线性回归诊断——模型检验(模型调优)

2.3.5.1 检验残差分析

- 1.建立线性回归模型
ana1 = ols('avg_exp ~ Income',exp).fit()
ana1.summary()

- 2.散点图：计算预测值与残差值并散点图分析
# 回归模型
ana1 = lm_s
# 保留预测值
exp['Pred'] = ana1.predict(exp)
# 保留残差
exp['resid'] = ana1.resid
# 收入和残差进行分析，作收入的预测值
exp.plot('Income','resid',kind = 'scatter')   
plt.show()

模型解释：
图形呈扇形(喇叭状)，所以是异方差，其解决方法是取对数。

2.3.5.2 模型调优(检验)

• 根据残差分析的散点图进行调优

- 3.建立新的线性回归模型，计算预测值、残差值
exp['avg_exp_ln'] = np.log(exp['avg_exp']) #对Y取对数
ana2 = ols('avg_exp_ln ~ Income',exp).fit()
exp['Pred'] = ana2.predict(exp)
exp['resid'] = ana2.resid
ana2.summary()

- 4.新的散点图分析
exp.plot('Income','resid',kind='scatter')

模型解释：
根据结果$R^2$比未取对数时下降了，所以模型效果变差了。需要再次调优。

2.3.5.3 模型调优

因为根据散点图判断残差分析是异方差，所以需要对Y、X都取对数再次检验。

- 5.建立新的线性回归模型，判断残差分析
# 对X取对数
exp['Income_ln'] = np.log(exp['Income'])
# 回归模型
ana3 = ols('avg_exp_ln ~ Income_ln',exp).fit()
# 预测Y值
exp['Pred'] = ana3.predict(exp)
# 计算残差
exp['resid'] = ana3.resid
# 散点图
exp.plot('Income_ln','resid',kind='scatter')

- 6.模型效果
ana3.summary()

模型解释：
$R^2$增加，模型效果变好；一个人的年收入每增加一个百分点，支出增加0.89个百分点。

2.4 强影响点

强影响点作用：构建关于强影响点的统计量去更好识别强影响点，如：学生化残差。

学生化残差 = 残差/标准误差，将样本残差作学生标准化。
强影响点的识别：
样本量为三位数，|SR|>2为强影响点；
样本量为四位数，|SR|>3为强影响点。

2.4.1 案例：强影响点

- 1.学生化残差:标准化残差
exp['resid_t'] = (exp['resid'] - exp['resid'].mean())/exp['resid'].std()
exp['resid_t'].sort_values()

- 2.选择标准化残差小于等于2的样本，利用数据进行建模
# 样本量为三位数，|SR|>2为异常值
exp2 = exp[abs(exp['resid_t'])<=2].copy()  
ana4 = ols('avg_exp_ln ~ Income_ln',exp2).fit()
exp2['Pred'] = ana4.predict(exp2)
exp2['resid'] = ana4.resid
exp2.plot('Income','resid',kind='scatter')
plt.show()

去除异常值，异常值只去除一次

- 3.回归模型
- R^2上升，模型效果变好了
ana4.summary()

3 共线性问题——X与X之间的关系

• 若X之间有共线性问题，模型是不稳定的。方差膨胀因子>10表示某变量的多重共线性严重。

• 方差膨胀因子公式：
  $$V{I}_i=\frac{1}{1-R_i^2}$$
  例：模型$Y={\beta }_1{X}_1+{\beta }_2{X}_2$检验共线性：

  1. 利用X1和X2做回归即 $\overline{X_1}=\beta X_2$,得到$R^2$
  2. 将$R^2$带入方差膨胀因子公式，方差膨胀因子>10则存在多重共线性。

- 1.建立回归模型
exp2['dist_home_val_ln'] = np.log(exp2['dist_home_val'])
exp2['dist_avg_income_ln'] = np.log(exp2['dist_avg_income'])
ana5 = ols('avg_exp_ln ~ Age + Income_ln + dist_home_val_ln + dist_avg_income_ln',exp2).fit()

ana5.summary()

- 2. 方差膨胀因子(VIF)
def vif(df,col_i):
	 '''
    将每个X与剩下的X做线性回归，得到R^2
    返回VIF = 1. / (1. - r2)
    '''
    cols = list(df.columns)
    cols.remove(col_i)
    cols_noti = cols
    formula = col_i + '~' + '+'.join(cols_noti)
    r2 = ols(formula, df).fit().rsquared
    return 1. / (1. - r2)

exog = exp2[['Income_ln', 'dist_home_val_ln',
             'dist_avg_income_ln']]

for i in exog.columns:
    print(i, '\t', vif(df=exog, col_i=i))

模型解释：
步骤一： 方差膨胀因子>10的表示某变量多重共线性严重。可看出Income、avg_income存在多重共线性问题。
步骤二：对模型进行调整，删除不显著的变量。图中存在共线性问题且不显著的是Income_ln，所以删除Income_ln重新建立回归模型。

ana6 = ols('avg_exp_ln ~ Age + dist_home_val_ln + dist_avg_income_ln',exp2).fit()

ana6.summary()

注：根据一个标准(如VIF、AIC等)，一般是不能将所有不显著的变量都去除掉的，需要结合使用。

4 线性回归建模流程★

1. 初始分析：分析X、Y的分布【单变量：Y是否右偏、X和Y有无异常值】
2. 变量选择：观察X和Y之间的关系，保留显著的X并放入模型【检验X、Y是否显著】
3. 验证模型假设：X和Y的函数关系【看散点图】
4. 多重共线性和强影响点：检验共线性和影响点的问题

所有后续模型的检验都少于线性模型，因为只有线性模型做的特别精确。

5 线性回归假设归纳★★★

多元线性回归的六个假设：
假设一：模型设置，选择何种回归方法、如何选变量、变量以何种形式放入模型【看散点图】
假设二：解释变量和扰动项不相关【理论或常识判断，无法检验】
假设三：解释变量之间不能强相关【方差膨胀因子】
假设四：扰动项独立同分布【异方差检验、DW检验】
假设五：扰动项服从正态分析【看直方图是否是正态分布】
注：只要Y不是正态分布，残差必然是异方差