Stata命令大全面板数据计量分析与软件实现

1、Stata命令大全面板数据计量分析与软件实现说明：以下do文件相当一部分内容来自于中山大学连玉君STAT徽程,感谢他的贡献。本人做了一定的修改与筛选。* 面板数据模型* 1.静态面板模型：FE和RE* 2.模型选择：FEvsPOLS,REvsPOLS,FEvsRE（pols混合最小二乘估计）* 3.异方差、序列相关和截面相关检验* 4.动态面板模型（DID-GMM,SYS-GMM* 5.面板随机前沿模型* 6.面板协整分析（FMOLS,DOLS*说明：1-5均用STATAfc件实现，6用GAUSS件实现。*生产效率分析（尤其指TFP）:数据包络分析（DEA与随机前沿分析（SFA*说明：DE和D

2、EAP2.依件实现，SF/ttFrontier4.1实现，尤其后者，侧重于比较C-D与Translog生产函数，一步法与两步法的区别。常应用于地区经济差异、FDI溢出效应（SpilloversEffect）、工业行业效率状况等。*空间计量分析：SLM（型与SEM1型* 说明：STATArMatlab结合使用。常应用于空间溢出效应（R&D、财政分权、地方政府公共行为等。* 一、常用的数据处理与作图* 指定面板格式xtsetidyear（id为截面名称，year为时间名称）xtdes/添加标签或更改变量名labelvarh"人力资本"数据特征*/xtsumlogyh/*

3、数据统计特征*/sumlogyh/*数据统计特征*/renamehhum* 排序sortidyear/*是以STATAf板数据格式出现*/sortyearid/*是以DEA&式出现*/* 删除个别年份或省份dropifyear<1992dropifid=2/*注意用=*/* 如何得到连续year或id编号(当完成上述操作时，year或id就不连续，为形成panel格式，需要用egen命令)egenyear_new=group(year)xtsetidyear_new* *保留变量或保留观测值keepinv/*删除变量*/* *日keepifyear=2000* *排序sortid

4、year/*是以STATA3板数据格式出现sortyearid/*是以DEA&式出现* *长数据和宽数据的转换* 长>>>宽数据reshapewidelogy,i(id)j(year)* 宽>>出数据reshapelogy,i(id)j(year)追加数据(用于面板数据和时间序列)xtsetidyear*或者xtdestsappend,add(5)/表示在每个省份再追加5年，用于面板数据/tsset* 或者tsdes.tsappend,add(8)/表示追加8年，用于时间序列/* 方差分解，比如三个变量Y,X,Z都是面板格式的数据，且满足Y=X+Z求方差v

5、ar(Y),协方差Cov(X,Y)和Cov(Z,Y)bysortyear:corrYXZ,cov* *生产虚拟变量* 生成年份虚拟变量tabyear,gen(yr)* 生成省份虚拟变量tabid,gen(dum)* *生成滞后项和差分项xtsetidyeargenylag=l.y/*产生一阶滞后项)，同样可产生二阶滞后项*/genylag2=L2.ygendy=D.y/*产生差分项*/* 求出各省2000年以前的openinv的平均增长率collapse(mean)openinvifyear<2000,by(id)变量排序，当变量太多，按规律排列。可用命令aorder或者orderfdi

6、openinsti*二、静态面板模型* 简介* 面板数据的结构（兼具截面资料和时间序列资料的特征）useproduct.dta,clearbrowsextsetidyearxtdes*固定效应模型* 实质上就是在传统的线性回归模型中加入N-1个虚拟变量,* 使得每个截面都有自己的截距项，* 截距项的不同反映了个体的某些不随时间改变的特征* 例如：lny=a_i+b1散点图+线性拟合直线+置信区间twoway(scatterlogyh)(lfitlogyh)(lfitcilogyh)lnK+b2*lnL+e_it* 考虑中国29个省份的C-Dft产函数*画图* 散点图+线性拟合直线twoway(

7、scatterlogyh)(lfitlogyh)* 散点图+二次拟合曲线twoway(scatterlogyh)(qfitlogyh)*按不同个体画出散点图和拟合线，可以以做出fevsre的初判断stata的估计方法解析*目的：如果截面的个数非常多，那么采用虚拟变量的方式运算量过大*因此，要寻求合理的方式去除掉个体效应*因为，我们关注的是x的系数，而非每个截面的截距项*处理方法：*twoway(scatterlogyhifid<4)(Ifitlogyhifid<4)(Ifitlogyhifid=1)(lfitlogyhifid=2)(lfitlogyhifid=3)* 按不同个体画

8、散点图,sobeautiful!*graphtwowayscatterlogyhifid=1|scatterlogyhifid=2,msymbol(Sh)|scatterlogyhifid=3,msymbol(T)|scatterlogyhifid=4,msymbol(d)|,legend(position(11)ring(0)label(1"北京”)label(2"天津")label(3"河北")label(4"山西")* *每个省份logy与h的散点图，并将各个图形合并twowayscatterlogyh,by(id)y

9、label(,format(%3.0f)xlabel(,format(%3.0f)* 每个个体的时间趋势图*xtlinehifid<11,overlaylegend(on)*一个例子：中国29个省份的C-D生产函数的估计tabid,gen(dum)list* 回归分析reglogylogklogldum*,eststorem_olsxtreglogylogklogl,feeststorem_feesttablem_olsm_fe,b(%6.3f)star(0.10.050.01)* Wald检验testlogk=logl=0testlogk=logl* y_it=u_i+x_it解读xt

10、reg,fe的估计结果b+e_it(1)* ym_i=u_i+xm_i*b+em_i(2)组内平均* ym=um+xm*b+em(3)样本平均* (1)-(2),可得：estimator*/* (y_it-ym_i)=(x_it-xm_i)*b+(e_it-em_i)(4)/*within* (4)+(3),方得：一一一一* (y_it-ym_i+ym)=um+(x_it-xm_i+xm)*b+(e_it-em_i+em)* 可重而表示为：* Y_it=a_0+X_it*b+E_it* 对该模型执行OLSM计，即可得到b的无偏估计量*stata后台操作，揭开fe估计的神秘面纱！egeny_me

11、anw=mean(logy),by(id)/*个体内部平均*/egeny_mean=mean(logy)/*样本平均*/egenk_meanw=mean(logk),by(id)egenk_mean=mean(logk)egenl_meanw=mean(logl),by(id)egenl_mean=mean(logl)gendyw=logy-y_meanwgendkw=logk-k_meanwgendlw=logl-l_meanwregdywdkwdlw,noconseststorem_statagendy=logy-y_meanw+y_meangendk=logk-k_meanw+k_mea

12、ngendl=logl-l_meanw+l_meanregdydkdleststorem_stataesttablem_*,b(%6.3f)star(0.10.050.01)xtreglogyhinvgovopen,fe* -RA2* y_it=a_0+x_it*b_o+e_it(1)pooledOLS* y_it=u_i+x_it*b_w+e_it(2)withinestimator* ym_i=a_0+xm_i*b_b+em_i(3)betweenestimator* ->R-sq:within模型(2)对应的R2,是一个真正意义上的R2* ->R-sq:betweencorr

13、xm_i*b_w,ym_iA2* ->R-sq:overallcorrx_it*b_w,y_itA2* -F(4,373)=855.93检验除常数项外其他解释变量的联合显著性* -corr(u_i,Xb)=-0.2347* -sigma_u,sigma_e,rho* rho=sigma_uA2/(sigma_uA2+sigma_eA2)dise(sigma_u)A2/(e(sigma_u)A2+e(sigma_e)A2)*个体效应是否显著？*F(28,373)=338.86H0:al=a2=a3=a4=a29*Prob>F=0.0000表明，固定效应高度显著*-如何得到调整后的R2

14、,即adj-R2?ereturnlistreglogyhinvgovopendum*-拟合值和残差* y_it=u_i+x_it*b+e_it* predictnewvar,option/*xbxb,fittedvalues;thedefaultstdpcalculatestandarderrorofthefittedvaluesueu_i+e_it,thecombinedresidualxbuxb+u_i,predictionincludingeffectuu_i,thefixed-orrandom-errorcomponentee_it,theoverallerrorcomponent*/

15、xtreglogylogklogl,fepredicty_hatpredicta,upredictres,epredictcres,uegenares=a+reslistarescresin1/10随机效应模型* y_it=x_it*b+(a_i+u_it)* =x_it*b+v_it* 基本应想：将疝机干扰项分成两种* 一种是不随时间改变的，即个体效应a_iu_it* 另一种是随时间改变的，即通常意义上的不扰项* 估计方法：FGLS* Var(v_it)=sigma_aA2+sigma_uA2* Cov(v_it,v_is)=sigma_aA2* Cov(v_it,v_js)=0* 禾1J用

16、PooledOLS,WithinEstimator,BetweenEstimator* 可以估计出sigma_aA2和sigma_uA2,进而采用GL峨FGLS* Re估计量是Fe古正最和Bef古计量的加权平均* yr_it=y_it-theta*ym_i* xr_it=x_it-theta*xm_i* theta=1-sigma_u/sqrt(T*sigma_aA2+sigma_uA2)* 解读xtreg,re的估计结果useproduct.dta,clearxtreglogylogklogl,re*-R2* ->R-sq:withincorr(x_it-xm_i)时间效应*单向固定效

17、应模型*y_it=u_i+x_it*b+e_it*双向而定效应模型*y_it=u_i+f_t+x_it*b+e_itb_r,y_it-ym_iA2* ->R-sq:betweencorrxm_i*b_r,ym_iA2* ->R-sq:overallcorrx_it*b_r,y_itA2* 上述R2®不是真正意分上向R2,因为Ref真型采用的是GLSfc计。* rho=sigma_uA2/(sigma_uA2+sigma_eA2)dise(sigma_u)A2/(e(sigma_u)A2+e(sigma_e)A2)* corr(u_i,X)=0(assumed)* 这无机效

18、应模型的一个最重要，也限制该模型应用的一个重要假设* 然而，采用固定效应模型，我们可以粗略估计出corr(u_i,X)xtregmarketinveststock,fe* Waldchi2(2)=10962.50Prob>chi2=0.0000* 时间效应、模型的筛选和常见问题* 目录* 7.2.1时间效应（双向固定（随机）效应模型）* 7.2.2模型的筛选* 7.2.3面板数据常见问题* 7.2.4面板数据的转换*quitabyear,gen(yr)dropyr1xtreglogylogkloglyrFEvsRE?y_it=u_i+x_it*b+e_it*-Hausman检验-*基本思

19、想：如果Corr(u_i,x_it)=0,Fe,fe*随机效应模型中的时间效应xtreglogylogkloglyr*,fe*模型的筛选* 固定效应模型还是PooledOLS?xtreglogylogkloglyr*,fe/*Wald检验*/quitabid,gen(dum)/*LR检验*/reglogylogklogl/*POLS*/eststorem_olsreglogylogklogldum*,noconseststorem_felrtestm_olsm_feesttablem_*,b(%6.3f)star(0.10.050.01)* REvsPooledOLS?* H0:Var(u)=

20、0和Re都是一致的，但Re£有效*方法一：B-P检验xtreglogylogklogl,rexttest0* 如果Corr(u_i,x_it)!=0,Fe仍然有效，但Rei有偏的* 基本步骤*情形1:huasma讷正数xtreglogylogklogl,feeststorem_fextreglogylogklogl,reeststorem_rehausmanm_fem_re*情形2:quixtreglogyhinvgovopen,feeststorefequixtreglogyhinvgovopen,reeststorerehausmanfere* Hausman检验值为负怎么办？*

21、 通常是因为RE真型的基本假设Corr(x,u_i)=0无法得到满足*检验过程中两个模型的方差-协方差矩库都采用Fe真型的hausmanfere,sigmaless*两个模型的方差-协方差矩阵都采用Ret真型的hausmanfere,sigmamore* =为何有些变量会被drop掉？usenlswork.dta,cleartssetidcodeyearxtregln_wagehourstenurettl_exp,fe/*正常执行*/* 产生诉族虚拟变量tabrace,gen(dum_race)xtregln_wagehourstenurettl_expdum_race2dum_race3,f

22、e* 为何dum_race2和dum_race3会被dropped?* 固定效应模型的设定：y_it=u_i+x_it*b+e_it(1)* 由于个体效应u_i不随时间改变，* 因此若x_it包含了任何不随时间改变的变量，* 都会与u_i构成多重共线性，Stata会自动删除之。*异方差、序列相关和截面相关问题*y_it=x_it*b+u_i+e_it* 由于面板数据同时兼顾了截面数据和时间序列的特征，* 所以异方差和序列相关必然会存在于面板数据中；* 同时，由于面板数据中每个截面（公司、个人、国家、地区）之间还可能存在内在的联系，* 所以，截面相关性也是一个需要考虑的问题。*同方差假设序列无关

23、假设截面不相关假设* 此前的分析依赖三个假设条件:* (1)Vare_it=sigmaA2* (2)Corre_it,e_it-s=0* (3)Corre_it,e_jt=0*当这三个假设无法得到满足时，便分别出现异方差、序列相关和截面相关问题;我们一方面要采用各种方法来检验这些假设是否得到了满足；另一方面，也要在这些假设无法满足时寻求合理的估计方法。假设检验*=组间异方差检验（截面数据的特征）* Var（e_i）=sigma_iA2* Fe模型xtreglogylogklogl,fexttest3* Re模型* Re本身已经较大程度的考虑了异方差问题，主要体现在sigma_uA2上*=序列相

24、关检验*Fe模型* xtserialWooldridge(2002),若无序列相关，则一阶差分后残差相关系数应为-0.5xtseriallogylogkloglxtseriallogylogklogl,output* Re模型xtreglogylogklogl,rexttestl/估计方法*=异方差稳健型估计xtreglogyhinvgovopen,ferobusteststorefe_rb提供多个统计检验量*/*=截面相关检验*xttest2命令H0:所有截面残差的相关系数都相等xtreglogylogklogl,fexttest2*由于检验过程中执行了SU的计，所以要求T>Nxtre

25、glogylogkloglifid<6,fexttest2*xtcsd命令(提供了三种检验方法)xtreglogylogklogl,fextcsd,pesaran/*Pesaran(2004)*/xtcsd,friedman/*Friedman(1937)*/xtreglogylogklogl,rextcsd,pesaranxtreglogyhinvgovopen,ferobusteststorefe*结果对比esttabfe_rbfe,b(%6.3f)se(%6.3f)mtitle(fe_rbfe)*=序列相关估计* 一阶自相关xtregar,fe/re* 模型：y_it=u_i+x_

26、it结果对比xtreglogyhinvgovopen,feeststorefelocalmodels"fefe_ar1re_ar1fe_ar1_two"esttab'models',b(%6.3f)se(%6.3f)mtitle('models')r2sca(r2_wcorr)b+v_it(1)* v_it=rho*v_it-1+z_it(2)xtregarlogyhinvgovopen,feeststorefe_ar1xtregarlogyhinvgovopen,feIbi/*Baltagi-WuLBItest*/* 说明：* (1)这里的

27、Durbin-Watson=1.280677具有较为复杂的分布，* 不同于时间序列中的D-W充计量。* (2)其临界值见Bhargavaetal.(1982,TheReviewofEconomicStudies49:553-549)* (3)Baltagi-WuLBI=1.4739834基本上没有太大的参考价值，* 因为他们并未提供临界值表，而该统计量的分布又相当复杂xtregarlogyhinvgovopen,reeststorere_ar1* 两阶段估计xtregarlogyhinvgovopen,fetwostepeststorefe_ar1_two* 高阶自相关* newey2命令ne

28、wey2logyhinvgovopen,lag(2)*=组间相关(截面相关)*cluster选项usextcs.dta,clearxtreglogyhinvgovopen,fecluster(id)eststorefe_clusterxtreglogyhinvgovopen,recluster(id)eststorere_cluster*看过去很晕，采用一种综合处理:* xtgls命令* xtpcse命令* -简介(Greene,2000,chp15)*Iy_n|X_n|e_n|*Ee_i*e_i'=s_iA2*|s1A20.0|0s2A2.0|.|V=|.|.|00.snA2|*-截

29、面相关* Ee_i*e_i'=s_ijA2* |s_11s_12.s_1n|* |s_21s_22.s_2n|* |.|* V=|.|*sigmaA2* |.|* |s_n1s_n2.s_nn|*-序列相关* Ee_i*e_i'=s_iA2*M_i*|s1A2*M_10.0|0s2A2*M_2.0|.|V=|.|.|00.snA2*M_n|*GLS估计*b=X'VA-1检验截面相关xtglsmarketinveststock,panel(het)xttest2XA-1X'VA-1y*Varb=X'VA-1*XA-1*-估计和检验*=xtgls命令usei

30、nvest2.dta,clearxtglsmarketinveststock,panels(iid)/*iid,eststoreg_0regmarketinveststockeststoreg_olsxtglsmarketinveststock,panel(het)/*eststoreg_phetxtglsmarketinveststock,corr(ar1)/*数*/eststoreg_par1xtglsmarketinveststock,corr(psar1)/*/eststoreg_psar1xtglsmarketinveststock,panel(corr)/*eststoreg_pc

31、orrxtglsmarketinveststock,p(c)corr(ar1)eststoreg_all* 检验异方差xtglsmarketinveststock,panel(het)/*xttest3* 检验序列相关等同于PooledOLS*/截面异方差*/所有截面具有相同的自相关系每个截面有自己的自相关系数截面间相关且异方差*/截面异方差*/xtserialmarketinveststock结果对比xtregmarketinveststock,feeststorefelocalmodels"feg_0g_ols"esttab'models',b(%6.3

32、f)se(%6.3f)mtitle('models')r2sca(r2_w)localmodels"feg_phetg_par1g_psar1g_pcorrg_all"esttab'models',b(%6.3f)se(%6.3f)mtitle('models')r2sca(r2_w)compress* 说明：* 为何xtgls不汇报R2?* 因为此时的R殊必介于0和1之间，不具有传统线性回归模型中R2勺含义*=xtpcse命令* 默认假设：存在截面异方差和截面相关* 估计方法：OLS或Prais-Winsten回归* 有别于

33、xtgls(采用FGLSf计)* 更适于方块面板N不大(10-20),T不大(10-40)* 与xtgls的区别：估计方法不同* xtgls采用GLSS行估计，而xtpsce采用OLSuseinvest2.dta,clearxtpcseinvestmarketstockeststorepcse_full/结果对比xtreginvestmarketstock,feeststorefeOLS估计，调整异方差和截面相关后的标准误*/xtglsinvestmarketstock,panels(correlated)eststorem_xtgls/*FGLS估计，异方差和截面相关*/xtpcseinve

34、stmarketstock,correlation(arl)eststorepcse_ar1/*Prais_Winsten估计，共同的自相关系数*/xtpcseinvestmarketstock,correlation(arl)hetonlyeststorepcse_ar1/*不考虑截面相关*/localmodels"fepcse_fullm_xtglspcse_ar1pcse_ar1"esttab'models',b(%6.3f)se(%6.3f)mtitle('models')r2sca(r2_w)*xtpcse的结果与xtgls非常相似

35、，但前者可以汇报R2* 当Ntt大时，采用该方法会非常费时，* 因为方差协方差矩阵是采用OLSfc计的残差计算的usextcs.dta,clearxtdesxtpcsetlsizendtstangtobinnpr/=简介=*大约5-8分钟*/eststorextpcsextregtlsizendtstangtobinnpr,feeststorefe* 结果对比localmodels"fextpcse"esttab'models',b(%6.3f)se(%6.3f)mtitle('models')r2sca(r2_w)* 系数估计值有较大差别，

36、但符号和显著性是一致的。一*动态面板模型PartI* 动态面板模型* 简介* 一阶差分IV估计量(AndersonandHisao,1982)* 一阶差分GMM计量(ArellanoandBond,1991)* 系统GMMr计量(AB,1995;BB,1998)*模型：yit=a0*yit-1+a1*xit+a2*wit+u_i+eit*特征：解释变量中包含了被解释变量的一阶滞后项*可以是非平行面板，但要保证时间连续xit严格外生变量Ex_it,e_is=0foralltands即，所有干扰项与/b4相关wit先决变量Ew_it,e_is!=0fors<t,butEx_it,v_is=0

37、foralls>=t即，前期干扰项与当期x相关，但当期和未来期干扰项与x不相关。yit-1内生变量Ex_it,e_is!=0fors<=t即，前期和当期“当期干扰项与x相关u_i随机效应，在截面间是iid的。u_i与eit独立。*内生性问题:*(1)*(2)*若假设u_i为随机效应，则Corr(yi,t-1,u_i)!=0若假设u_i为个体效应，需要想办法去除之，因为数据为“大N、T”一阶差分：D.yi,t-1=yi,t-1-yi,t-2D.ei,t=ei,t-ei,t-1显然：Corr(D.yi,t-1,D.ei,t)!=0,差分方程存在内生问题;组内去心：ymi,t-1=yi,

38、t-1-1/(T-1)*(yi,t-1+yi,T)emi,t=ei,t-1/T*(ei,t+ei,t-1+ei,T)显然：Corr(ymi,t-1,emi,t)!=0,仍然存在内生性问题*处理办法：IV估计或GMM计，选择合适的工具变量*矩条件：Ee_it,z_it=0*=*=一阶差分IV估计量=*=AndersonandHisao(1982)=*=*基本思想：采用一阶差分去除个体效应u_i,y的滞后二阶作为D.yit-1的工具变量同时，D.yit-2也可以作为D.yit-1的工具变量useabdata.dta,cleardes/*变量的定义*/tssetidyear模型：n_it=b1*n_

39、it-1+b2*n_it-2+b3*w_it+b4*w_it-1+b5*k_it+b6*k_it-1+b7*k_it-2+b8*ys_it+b9*ys_it-1+b10*ys_it-2xtivregnL2.nwL1.wkL1.kL2.kysL1.ysL2.ysyr1981-yr1984/(L.n=L3.n),fdxtivregnL2.nL(0/1).wL(0/2).(kys)yr1981-yr1984(L.n=L3.n),fd*=* =一阶差分gmM计量=*=AB91=* L.Hansen(1982)提出GMM* ArellanoandBond(1991)* 模型：* yit=a0*yit-1

40、+a1*xit+a2*wit+u_i+vit* 假设条件：* 干扰项vit不存在序列相关；* 适用范围：* 大N,小T* 随后，我们会介绍“小N庆T型动态面板的估计方法* 基本思想：在AndersonandHisao(1982)基础上增加了更多可用的工具变量* 在t=3处，y_i1可以作为所有滞后项的工具变量* 在t=4处,y_i1,y_i2可以作为所有滞后项的工具变量* D.yit=a1*D.yit-1+a2*D.Xit+D.vitX_it=x_it,w_it* 因此，所有工具变量构成的矩阵如下：*|y_i100000.0.0D.x_i3|0y_i1y_i2000.0.0D.x_i4|000

41、y_i1y_i2y_i3.0.0D.x_i5|*Z_i=|.|000000.y_i1.yiT-2D.x_iT|Z_i的行数为T-2Z_i的列数为sum_(m=1)A(T-2)m+K,K为X的列数以T=7,K=3为例，WJZ_i的列数为(1+2+3+4+5)+3=18* 设定工具变量的基本原则：* 对-内生-变量的处理：与上述方法类似，* 即滞后两阶以上的水平变量均可作为差分方程的工具变量(GMMtype)* 对-先决-变量的处理：滞后一阶以上的水平变量均可作为工具变量(GMMtype)* 对-外生-变量的处理：自己作为自己的工具变量(StandardIV)* =例1:一阶差分估计量的基本设定* 解释变量仅包含y_it的一阶滞后项，默认设定* 干扰项同方差，一阶段估计us