用Python底层编写进行计量经济分析（三）：异方差（原因、结果、检验：White检验、补救：广义线性回归）

淡淡的烟草味﹌ 2023-07-24 02:22 19阅读 0赞

系列文章：

1.  [多元线性回归和显著性检验（参数估计、T检验、F检验、拟合优度）][T_F]
2.  [多重共线性（导致结果、检验——方差膨胀因子、补救措施——岭回归）][Link 1]
3.  [异方差（导致结果、检验——White、补救措施——广义线性回归）][White]
4.  [自相关（导致结果、检验——D-W、补救措施——广义线性回归）][D-W]

写的比较仓促，代码中如有错误欢迎指正！

多元线性回归的基本假定：

1.  模型符合线性模式
2.   X X X满秩（无多重共线）
3.  零均值价值： E ( ε i ∣ X i ) = 0 E(ε\_i|X\_i)=0 E(εi∣Xi)=0（自变量外生）
4.  同方差： V a r ( ε i ∣ X i ) = σ Var(ε\_i|X\_i)=σ Var(εi∣Xi)=σ
5.  无自相关： c o v ( ε i , ε j ) = 0 cov(ε\_i, ε\_j)=0 cov(εi,εj)=0
6.  球形扰动：ε\_i是正态分布

若果模型违反了相应的假设就会犯对应的错误，我们在计量经济分析中的检验就是检验出是否可能犯了某一类错误，若果极有可能犯了一种错误时，我们应该怎么修正它，才能保证分析的结果是有效的。

#### 一、异方差性 ####

**和前面一样，举个例子。如果，我们一个理论上的多元线性回归模型为： Y i = β 0 + β 1 ∗ x 1 i + β 2 ∗ x 2 i + β 3 ∗ x 3 i + u i Y\_i=β\_0+β\_1\*x\_\{1i\}+β\_2\*x\_\{2i\}+β\_3\*x\_\{3i\}+u\_i Yi=β0\+β1∗x1i\+β2∗x2i\+β3∗x3i\+ui  
但我们再建模时不小心遗漏了一个变量 x 3 x\_\{3\} x3，此时我们实际建立的模型变为 Y i = β 0 + β 1 ∗ x 1 i + β 2 ∗ x 2 i + u i Y\_i=β\_0+β\_1\*x\_\{1i\}+β\_2\*x\_\{2i\}+u\_i Yi=β0\+β1∗x1i\+β2∗x2i\+ui  
如果 x 3 x\_\{3\} x3本身具有一定的趋势，这样遗漏一个变量会对我们的估计产生后果？**

**1.1.异方差**

在多元线性回归中，有一条同方差的假定： V a r ( ε i ∣ X i ) = σ Var(ε\_i|X\_i)=σ Var(εi∣Xi)=σ，我们所说的异方差就是违反了同方差的假定，即出现了 V a r ( ε i ∣ X i ) = σ i Var(ε\_i|X\_i)=σ\_i Var(εi∣Xi)=σi，当异方差出现时，我们需要对我们的模型进行一定的修正，才能保证参数检验是有效的。

**1.2.异方差的原因**

V a r ( ε i ∣ X i ) = σ i Var(ε\_i|X\_i)=σ\_i Var(εi∣Xi)=σi说明了扰动项带有了一定的趋势。这往往是由下面几个因素造成的：

1.  模型中缺少了解释变量  
    我们看上面最初例举的例子，我们不小心忽略了 x 3 x\_\{3\} x3， x 3 x\_\{3\} x3本身具有一定趋势性，**此时 x 3 x\_\{3\} x3的变化就会隐藏在 u i u\_i ui中，使得 u i u\_i ui方差不再是一个常数**
2.  样本观测误差  
    随着统计口径的不断改变，或者统计技术的不断进步，我们的被解释变量的观测精度也就是观测误差会发生变化，这些误差的变化会表现在 u i u\_i ui上，导致出现异方差

**1.3.异方差的后果**

1.  参数估计仍然是无偏的  
    因为 u u u虽然自身有一定趋势，但是和 x x x不相关，此时参数的估计仍然是无偏的（具体推导请查书，大体思想是是 β β β估计值等于 β + f ( x u ) β+f(xu) β\+f(xu), x u = 0 时 f ( x u ) = 0 xu=0时f(xu)=0 xu=0时f(xu)=0，因此估计仍是无偏的）。具体几何解释可以理解为， u u u和 x x x仍然是正交的， y y y在 x x x空间的投影不会发生变化。
2.  参数估计方差不再是最小的  
    之前最小二乘的估计具有无偏和方差最小的性质。 β β β估计值等于 β + f ( x u ) β+f(xu) β\+f(xu)，当出现异方差时， β β β的方差估计会出现 E ( u i 2 ) E(u\_i^2) E(ui2)项，由于 u i u\_i ui的方差不再， β β β的方差估计不再是最小的（非有效的）。具体的非有效程度取决于异方差的程度。**此时我们的t检验值会被高估，显著性被放大，此时的t检验是无效的！！！**

**1.4.异方差的检验（White检验）**  
有多多种检验异方差的方法，比如Glodfeld-Quandt检验、Breusch-Pagan检验、White检验等。本文对使用最多的White检验进行介绍和python实现。  
先看一下，异方差 V a r ( ε i ∣ X i ) = σ i Var(ε\_i|X\_i)=σ\_i Var(εi∣Xi)=σi的零假设和备择假设  
H0: σ i σ\_i σi= σ σ σ  
H1:非H0  
可以看到H0中等式的个数和样本相同，我们对这种H0的检验是不可能实现的。因此需要利用其他方法对异方差检验。比如White检验（公式太难打，网上直接搜了一个wihte检验的课件）：  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 1]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 2]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 3]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 4]  
很多人刚开始看到这种检验方法和我一样是懵逼的。为什么 x x x和交叉项对 e i e\_i ei的回归可以判断是否存在异方差？？？？回到第一节，最小二乘得到的 β β β的估计量 b = ( X ′ X ) − 1 X ′ y = β + ( X ′ X ) − 1 X ′ ε b=(X'X)^\{-1\}X'y=β+(X'X)^\{-1\}X'ε b=(X′X)−1X′y=β\+(X′X)−1X′ε，其方差的一直估计量为 s 2 ( X ′ X ) − 1 s^2(X'X)^\{-1\} s2(X′X)−1，存在异方差时，由于估计不再是有效的，则上面结论就不再成立。根据此现象便设计了Wihte检验。**但White检验检验存在几个问题：1.由于还好和交叉项回归，自由度缺失严重；2.检验是非建设性的，不能定位到那里导致了异方差；3.无法确定检验的功效**

**1.5.异方差的补救——广义最小二乘**  
异方差的补救方式有多种，第一种是将模型变为对数形式，降低数据的变化趋势，但也只能减小异方差并不能消除异方差，同时还改变了参数意义（**采取对数变化后，变量前系数表示弹性**）。第二种是加权最小二乘，对较小的 e i 2 e\_i^2 ei2赋予较大的权重，较大的 e i 2 e\_i^2 ei2赋予较小的权重，在进行估计的方法称为加权最小二乘。还有一种方法，也是最常用的方法：**广义最小二乘**。**后面我们可以看到广义最小二乘不仅可以补救异方差，还可以补救异方差问题，现实问题中，广义最小二乘非常常用。我们设 E ( ε ε ′ ) = σ 2 Ω E(εε')=σ^2Ω E(εε′)=σ2Ω，不存在异方差和自相关时 Ω = I Ω=I Ω=I。**  
**使用时，可分为 Ω Ω Ω已知和 Ω Ω Ω未知的情况。在已知的时候使用广义最小二乘（GLS），在 Ω Ω Ω未知的时，使用可行广义最小二乘（FGLS）**（公式太难打，直接粘贴网上广义最小二乘的标准步骤）：  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 5]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 6]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 7]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 8]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 9]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 10]  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 11]  
在 Ω Ω Ω未知的时候，需要找到 Ω Ω Ω的一致估计量，由上面可知，我们需要求解一个 ρ ρ ρ得到 Ω Ω Ω的一致估。什么是 ρ ρ ρ呢，就是 e t 和 e t − 1 e\_\{t\}和e\_\{t-1\} et和et−1的相关系数估计，与下一节说的检验自相关的DW统计量关系为： D W = 2 ( 1 − ρ ) DW=2(1-ρ) DW=2(1−ρ)

#### 二、异方差检验的python实现 ####

因为后面的自相关的补救也会用到广义最小二乘，因此广义最小二乘的python实现放到下一节。此处对异方差的White检验进行实现（前面的估计检验部分直接沿用前几节）：

class mul_linear_model(object):
      #用于计算数据的相关系数矩阵
      def data_relation(self):
          pass
      #用于计算每个变量的VIF
      def VIF_cul(self, ind):
          pass
      #所有变量的分行差膨胀因子检验
      def VIF_test(self):
          pass
      #计算拟合的条件数
      def  condition_num(self):
      	  pass
      #岭回归
      def Ridge_fit(self, k = 0.1):
      	  pass
      #用于计算数据的相关系数矩阵
      def data_relation(self):
          ...
      #用于计算每个变量的VIF    
      def VIF_cul(self, ind):
          ...
      #所有变量的分行差膨胀因子检验
      def VIF_test(self):
          ...
      #计算拟合的条件数
      def condition_num(self):
          ...
      #岭回归
      def Ridge_fit(self, k = 0.1, output=True):
         ...
      #异方差怀特检验
        def White_test(self, alpha):
            #算出y_hat
            y_hat = np.array(np.dot(self.X,self.b).T)[0]
            #计算偏差
            e_square = (y_hat-self.y)*(y_hat-self.y)
            import matplotlib.pyplot as plt
            plt.plot(y_hat-self.y)
            plt.title('ei')
            #计算white检验的自变量white_x,包括每一个x,和x的平方，以及两两之间的交叉项
            white_x = self.data_x.copy()
            x_columns = white_x.columns
            #计算平方项与交叉项，第一列是截距项，无需计算
            for i in range(1, len(x_columns)):
                for j in range(1, i+1):
                    white_x[x_columns[i]+'*'+x_columns[j]] = white_x[x_columns[i]]*white_x[x_columns[j]]
            
            #调用自身类，进行回归求R2,white_x已经包含截距项，无需再添加截距项
            white_model = mul_linear_model(e_square, white_x, intercept=False)
            #求回归参数
            white_model.fit(output=False)
            #求拟合优度
            white_R_2 = white_model.R_square()
            #计算统计量自由度
            white_df = white_model.K-1
            #计算white统计量
            white = white_model.N*white_R_2
            #进行假设检验
            if white > st.chi2.ppf(1-alpha, df = white_df):
                print("white statistic is: {0}, refuse H0".format(white))
            else:
                print("white statistic is: {0}, can't refuse H0".format(white))

我们仍使用前面的房价例子（**其实改例子放在这里不是十分恰当，我们大多数出现异方差问题的回归都是时序性的，而这里当房价例子是横截面数据，放着这里只是为了试验使用**）进行试验：

linear_model = mul_linear_model(house_data['value'],house_data[['HouseAge', 'AveRooms', 'AveBedrms', 'Population']])
    linear_model.fit()
    print('\n')
    linear_model.T_test(0.05)
    print('\n')
    print("R_2 estimate is:", linear_model.R_square(), '\n')
    linear_model.F_test(0.05)
    print('\n')
        
     linear_model.White_test(0.05)
    print('\n')

![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 12]  
检验结果如上图，可以看到White统计量为19.22，小于于卡方检验临界值，我们不能拒绝原假设，认为没有表现出明显的异方差问题。通过下图我们也可以看到，残差图虽然在一定范围内波动，但在前期存在一定的上升趋势，即可能存在自相关问题，我们还需要对自相关进行经一部的简单，[后面一节][D-W]会介绍自相关问题的检验和补救措施：广义线性回归  
![在这里插入图片描述][watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 13]

下一节将对广义线性回归进行实现，来解决异方差问题。

[T_F]: https://blog.csdn.net/u011517132/article/details/105286113
[Link 1]: https://blog.csdn.net/u011517132/article/details/105323135
[White]: https://blog.csdn.net/u011517132/article/details/105392242
[D-W]: https://blog.csdn.net/u011517132/article/details/105581371
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70]: /images/20230528/d5a36ebc886441d381ca9ecac0df23cd.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 1]: /images/20230528/a44e7b030fef4425bfa1ca78b95a0558.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 2]: /images/20230528/b8e7b93e4a144da387f35fcc41807d2a.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 3]: /images/20230528/6e54b2ae8b76461ca37deecc7b4849ab.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 4]: /images/20230528/bc051d8d3acb43c9ad8dfa616a20eb65.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 5]: /images/20230528/6e230bd550944785bd38ab996cc98013.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 6]: /images/20230528/ca24a75fc8f24d3ab0c0ea83fd417688.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 7]: /images/20230528/e8dde852955349239457c90af61fa11f.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 8]: /images/20230528/bccbfe5113c94eaaafd8d77276302c44.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 9]: /images/20230528/e98b287aa5404e83921519ff83dbc769.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 10]: /images/20230528/472489de9c0440b28d41d74bcb9e4086.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 11]: /images/20230528/0a14ab53da1048e2949b058a421089cd.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 12]: /images/20230528/74613a2df8f949a1accf6b8201ec2bb0.png
[watermark_type_ZmFuZ3poZW5naGVpdGk_shadow_10_text_aHR0cHM6Ly9ibG9nLmNzZG4ubmV0L3UwMTE1MTcxMzI_size_16_color_FFFFFF_t_70 13]: /images/20230528/c44f126430464131a9b36e637f13b932.png