异方差性实验.doc

1.2 实验二 异方差性及其性质1.2.1 实验目的我们已经知道，在经典条件下，线性模型回归参数的OLS估计是具有最小方差的线性无偏估计量。随机误差项的异方差性，是线性回归模型中常见的不满足经典条件的情形。与满足经典条件的情形相比，当模型中出现异方差性时，模型参数的普通最小二乘（OLS）估计的统计性质将发生什么样的变化？如何理解和把握这些变化？如何纠正模型估计因为异方差性而产生的问题？通过本实验，可以帮助学生理解异方差性本身的概念、存在异方差性时模型参数的OLS估计量的性质、加权最小二乘法等。1.2.2 实验背景与理论基础1. 异方差性本实验以二元线性回归模型为例进行说明。线性回归模型，假设模型满足除“同方差性”之外的所有经典假设：（1），或表示为，从而有；（3），随机误差无序列相关； （4）解释变量是确定性变量，与随机误差项不相关：，（5）自变量之间不存在精确（完全）的线性关系。矩阵X是列满秩的：。（要求样本容量）（6）随机误差的正态性：，。2. 异方差性条件下OLS估计量的统计性质（1）的无偏性：模型回归参数的OLS估计量为:可以证明，即使在异方差性条件下，上述估计量依然满足无偏性：（2）的方差及协方差：在模型满足经典条件时，OLS估计量的方差协方差矩阵为，但是在异方差性条件下，不存在独立于X的随机误差项方差，因此不再存在这一简单公式。另一方面，在计量分析实践中，即使在线性回归模型的经典条件下，随机误差项的方差本身也不是可直接观察的，实践中我们用对其进行估计（大多数统计分析软件正是如此处理的），也即用矩阵去估计OLS估计量的方差协方差矩阵，并在此基础上对模型进行各种检验。在线性回归模型的经典条件下，这种估计将是无偏的（参见本章实验一）。重要的问题是，在异方差性条件下，如果无视异方差性的存在，仍用去估计OLS估计量的方差协方差矩阵，这种估计是否仍具有无偏性？建立在这种估计之上的各种模型检验是否依然有效？3. 加权最小二乘法修正异方差性的常用方法是加权最小二乘法，它是广义最小二乘法中的一种。具体的方法是：如果模型中的标准差为，在原模型中乘以，模型变为：将此模型看成是对的线性回归模型，此模型将具有同方差性，由于原模型满足除同方差性外的所有经典条件，因此上述模型将满足所有线性模型的经典条件，从此模型中利用最小二乘法估计参数将获得具有最小方差的线性无偏估计量，这就是加权最小二乘法及其原理。1.2.3 实验原理本实验仍然通过一个虚构的二元线性回归模型来展开。与本章实验一一样，我们首先设定一定二元线性回归模型的回归参数，取定解释变量的样本值。由于是存在异方差性的模型，不能再设定随机误差项的方差，但是我们可以设定随机误差项的方差与解释变量值之间的函数关系。这样，在总体上我们已经完全掌握了模型。接着我们使用Matlab进行模拟随机抽样，对于得到的每一个模拟随机样本，我们分别使用普通最小二乘法和加权最小二乘法得到模型参数的两种不同的估计量，分别记为和。反复以上模拟抽样和估计，我们将分别得到每个模型参数的普通最小二乘估计量和加权最小二乘估计量的样本。通过这两个样本，我们可以探讨普通最小二乘估计量和加权最小二乘估计量的统计性质，分析两者之间的共同性质和区别。1.2.4 实验过程和步骤1. 程序设计以下将实验过程通过编制一个简单的Matlab程序来进行。程序分为以下几个部分：（1）第一步，设置模型基本参数，解释变量的样本值，这一步与实验一的相应步骤是类似的。Matlab程序段如下：clearn=20; beta0=10;beta1=5;beta2=-3; x1=15*rand(n,1)+1;x2=10*rand(n,1)+1;e=ones(n,1);X=e,x1,x2;（2）第二步，反复抽取样本，进行普通最小二乘估计和加权最小二乘估计，并将估计结果保存在相应的向量中。Matlab程序段如下：b0=;b1=;b2=;sigma=;c0=;c1=;c2=;XX=X./x1,x1,x1;times=5000;for j=1:times uu=normrnd(0,se,n,1); u=2*x1.*uu;Y=beta0+beta1*x1+beta2*x2+u; b,bint,r=regress(Y,X); b0=b0;b(1);b1=b1;b(2);b2=b2;b(3);sigma=sigma,sum(r.2)/(n-3); YY=Y./x1; c,bint,r=regress(YY,XX); c0=c0;c(1);c1=c1;c(2);c2=c2;c(3); end代码解释： “b0=;b1=;b2=;sigma1=;”生成4个维数可变的动态向量，准备分别存放每次抽样所产生的普通最小二乘估计量及，其中为普通最小二乘回归的残差。 “c0=;c1=;c2=;sigma2=;”生成4个维数可变的动态向量，准备分别存放每次抽样所产生的加权最小二乘估计量估计量及，其中为加权最小二乘回归的残差。“XX=X./x1,x1,x1;”表示将矩阵X的每一列向量的元素对应除以列向量X1的元素，也即得到矩阵生成这一矩阵的目的将在下文中揭示，实际上是为下文中的加权最小二乘估计做准备。“times=5000;”设定反复抽样和回归的次数，你可以根据需要设定成另外的整数；“uu=normrnd(0,1,n,1); ”随机生成分布N(0,1)的简单随机样本，构成n维列向量， “u=2*X1.*uu”表示生成一个n维列向量u，其元素是列向量X1和uu的对应元素的乘积再乘以2，即这表示随机误差相互间不相关，但其标准差为。“Y=beta0+beta1*x1+beta2*x2+u;”利用以上生成的向量生成被解释变量Y的一个模拟样本：这是一个具有异方差性的二元线性回归模型的模拟样本，随机误差的标准差与解释变量x1的关系是。 “b,bint,r=regress(Y,X);”将 Y 对 X 进行普通最小二乘回归，获取估计结果参数，其中b为回归系数点估计向量，r为残差列向量。 “b0=b0;b(1);” 将的值逐个存入数表b0，使b0于循环结束时成为times维列向量。“b1=b1;b(2);”将的值逐个存入数表b1，使b1于循环结束时成为times维列向量。“b2=b2;b(3);”将的值逐列存入数表b2，使b2于循环结束时成为times维列向量。“sigma=sigma,sum(r.2)/(n-3);”将由回归残差计算所得的值逐列存入数表sigma1，使sigma于循环结束时成为times维列向量。 “YY=Y./x1;”将向量Y的分量分别除以向量x1的对应分量，得到列向量YY，即“c,bint,r=regress(YY,XX);”将 YY 对 XX （参见前文中的XX的构造）进行普通最小二乘回归，获取估计结果参数。实际上，此时YY 对 XX 的普通最小二乘回归正是对线性模型的普通最小二乘回归，由于的标准差，因此上述模型又等价于对此模型的普通最小二乘估计正是对原模型的加权最小二乘估计。因此上述命令中的c为即为原模型回归参数的加权最小二乘估计向量，r为残差列向量。 “c0=c0;c(1);” 将的值逐个存入数表c0，使c0于循环结束时成为times维列向量；“c1=c1;c(2);” 将的值逐个存入数表c1，使c1于循环结束时成为times维列向量；“c2=c2;c(3);” 将的值逐列存入数表c2，使c2于循环结束时成为times维列向量。2. 输出实验结果运行上述程序后，模型回归参数的普通最小二乘估计量的样本存放在向量b0，b1，b2中；模型回归参数的加权最小二乘估计量的样本存放在向量c0，c1，c2中，同时，两种估计法所产生的残差所计算的和的所有数据分别存放在向量sigma1和sigma2中。这些向量的值均被暂时保存于内存中，我们可以用相应的Matlab命令输出我们所希望获得的各种结果。这里需要特别指出的是，以下结果只是某一次实验（上述程序某一次运行）的结果，该程序的每一次运行都将使解释变量的值，以及被解释变量的随机样本发生变化，因此你实验时得到的结果与以下结果将会略有出入。（1）解释变量的值利用命令x1,x2可直接输出程序中生成的解释变量的值（经整理成为下表）。解释变量的值并没有在对被解释变量的反复抽样中发生变化。x1x210.12812.21051.23645.50751.24538.15883.85119.92849.80383.7311.86373.54776.51359.65610.47183.323511.76459.048711.3910.0842.26123.31897.81533.39317.62741.49756.29881.78383.30417.408211.13472.908911.48829.438711.91262.7398.17582.70799.322610.943 （2）的无偏性即使在异方差性条件下，依然分别是的无偏估计。在上述程序运行结束后，使用以下MATLAB命令：b_theo=beta0,beta1,beta2mean_b_ols=mean(b0),mean(b1),mean(b2)mean_b_wols=mean(c0),mean(c1),mean(c2)其中第一个命令输出模型被事先设定的回归参数值，记为b_theo；第二个命令输出回归参数的普通最小二乘估计量的样本均值，记为mean_b_ols；第三个命令输出回归参数的加权最小二乘估计量的样本均值，记为mean_b_wols。运行上述命令，将分别得到：b_theo = 10 5 -3mean_b_ols = 9.8762 5.0029 -2.9843mean_b_wols = 9.9659 4.9989 -2.9954实验结果显示，回归参数的普通最小二乘估计量和加权最小二乘估计量均是无偏估计量。（3）普通最小二乘估计量的方差协方差矩阵在模型存在异方差性时，如果无视异方差性，仍然用去估计随机误差项的方差，进而用矩阵去估计的方差协方差矩阵，这种估计具有无偏性吗？为此，我们使用以下Matlab命令：fangcha_est_ols=mean(sigma);cov_b_est_ols=fangcha_est_ols*(inv(X*X)cov_b_sample_ols=cov(b0,b1,b2)代码解释：“fangcha_est_ols=mean(sigma1);”计算普通最小二乘估计法产生的的均值，赋予变量fangcha_est_ols。“cov_b_est_ols=fangcha_est_ols*(inv(X*X)”计算并输出矩阵的均值，记为“cov_b_est_ols”。“cov_b_sample_ols=cov(b0,b1,b2)”计算并输出向量b0, b1, b2 的样本方差协方差矩阵，也即的样本方差协方差矩阵，记为“cov_b_sample_ols”。输出结果为：cov_b_est_ols = 27.0496 -1.8265 -1.8278 -1.8265 0.2447 0.0037 -1.8278 0.0037 0.3235cov_b_sample_ols = 18.1221 -0.9746 -2.2516 -0.9746 0.1969 0.0329 -2.2516 0.0329 0.3660这只是一次实验的结果。反复运行之前的程序和上述命令，我们将发现， 的方差的最小二乘估计量相对于本身的方差而言明显偏大，同样的结果出现在对的方差的估计上。相反，的方差的最小二乘估计量相对于本身的方差而言却明显偏小。这说明，当模型存在异方差性时，尽管回归参数的普通最小二乘估计量是无偏估计量，但是利用普通最小二乘法估计回归参数的方差，将出现系统性高估或低估，从而使建立在这种估计之上的各种检验（如t-检验和F-检验）无效。实际上，从上述结果可以看出，用矩阵去估计的方差协方差矩阵除了高估或低估估计量的方差外，对估计量之间的协方差的估计也是有偏的。两者元素连“成比例”都达不到，也就是说，不可能表示成任何常数与矩阵的乘积。（4）与的方差比较。加权最小二乘估计量实际上来自一个满足经典条件的线性模型的OLS估计，因此是具有最小方差的线性无偏估计量。相反，由于原模型存在异方差性，所以尽管其普通最小二乘估计量仍是线性无偏估计量，但是却不具有最小方差性。前面的（3）中我们已经得到了的样本方差协方差矩阵：cov_b_sample_ols = 18.1221 -0.9746 -2.2516 -0.9746 0.1969 0.0329 -2.2516 0.0329 0.3660再使用以下Matlab命令：cov_b_sample_wols=cov(c0,c1,c2)可以得到加权最小二乘估计量的样本方差协方差矩阵：cov_b_sample_wols = 5.0167 -0.3118 -0.6360 -0.3118 0.1058 0.0130 -0.6360 0.0130 0.1002实验结果显示，的样本方差5.0167, 0.1058, 0.1002，明显分别小于的样本方差18.1221, 0.1969, 0.3660，这至少说明了不再是具有最小方差的线性无偏估计量。实际上，加权最小二乘估计量是一个特殊的消除了异方差性的模型的普通最小二乘估计量，因此它们是具有最小方差的线性无偏估计量。（5）和的样本分布直方图1）的样本分布直方图Matlab命令：histfit(b0)输出图形：图1.5（根据前文数据，此分布的均值为 9.8762，方差为18.1221），即使在异方差性条件下，普通最小二乘估计量依然服从正态分布。图1.5 的样本b0的频数直方图（带正态拟合）2）的样本分布直方图Matlab命令：histfit(c0)输出图形：图1.6（根据前文数据，此分布的均值为 9.9659，方差为5.0167），加权最小二乘估计量服从正态分布。图1.6 的样本c0的频数直方图（带正态拟合）图1.7 的样本b1的频数直方图（带正态拟合）3）的样本分布直方图Matlab命令：histfit(b1)输出图形：图1.7（根据前文数据，此分布的均值为 5.0029，方差为0.1969），普通最小二乘估计量服从正态分布。4）的样本分布直方图Matlab命令：histfit(c1)输出图形：图1.8（根据前文数据，此分布的均值为 4.9989，方差为0.1058），加权最小二乘估计量服从正态分布。图1