不可不知的4个经典案例数据：面板VAR模型

不可不知的4个案例数据实现面板VAR模型

---

案例1：使用nlswork2数据

本部分完整操作请查看文章：一文读懂面板向量自回归PVAR学习手册

简介：我们通过分析年工作时间和小时收入之间的关系来说明pvar命令集的使用，Holtz-Eakin, Newey和Rosen(1988)在他们关于面板向量自回归的开创性论文中分析了这一关系。为了将我们的新程序与Stata内置的var命令集进行比较，我们还将新的pvar命令集应用于Lutkephol(1993)的 West Germany 时间序列数据。

数据介绍：我们使用了Stata提供的1968年至1978年全国纵向调查中14-26岁女性的子样本。我们的子样本包括2039名女性，她们在至少三轮调查中报告了工资(工资)和年度工作时间(小时)，其中两轮调查是连续进行的。Holtz-Eakin等人使用了相同的调查，但不同的时间段和不同的工作人员子样本，因此结果可能不是直接可比的。使用前四个滞后时间和工资作为工具，使用pvarsoc计算一到三阶面板VAR模型选择。

模型代码集合如下：

sjlog using pvar11, replace
webuse nlswork2, clear
xtdescribe
generate ln_wks = ln(wks_work)
sjlog close, replace

sjlog using pvar12, replace
pvar ln_wks ln_wage, fd
estimates store fd_2 
pvar ln_wks ln_wage, fd instlags(2/3) 
estimates store fd_2t3
pvar ln_wks ln_wage, fod 
estimates store fod_1
pvar ln_wks ln_wage, fod instlags(1/2) 
estimates store fod_1t2
sjlog close, replace

sjlog using pvar13, replace
estimates table fd_2 fd_2t3 fod_1 fod_1t2 , b(%3.2f) se(%3.2f) stats(N J J_pval) modelwidth(8)
sjlog close, replace

sjlog using pvar14, replace
pvar ln_wks ln_wage, fd instlags(2/3) gmmstyle
estimates store fd_2t3g
pvar ln_wks ln_wage, fod instlags(1/2) gmmstyle
estimates store fod_1t2g
estimates table fd_2t3g fod_1t2g, b(%4.2f) se(%4.2f) stats(N J J_pval) modelwidth(8)
sjlog close, replace

sjlog using pvar15, replace
pvarsoc ln_wks ln_wage, maxl(3) pvaropts(instlags(1/3) fod gmmstyle)
sjlog close, replace

sjlog using pvar16, replace
pvarsoc ln_wks ln_wage, maxl(3) pvaropts(instlags(2/4) fod gmmstyle)
sjlog close, replace

案例2：psidextract数据

简介：我们使用psidextract数据，我们复制了Holtz-Eakin、Newey和Rosen(1988)中表2所示的简化形式的面板VAR，使用了Stata中收入与动态(PSID)面板研究数据中对1976-1982年528名男性的观察。

在他们最初的分析中，Holtz-Eakin, Newey和Rosen(1988)使用了1968年至1981年期间观察到的898名男性样本，使用了年工作时间和年平均时薪。所以我们的样品成分和时间段与原始的文章略有不同。在本例中，假设对数转换工资率(lwage)和工作周数(lwks)是每个变量的三个滞后函数。我们还假设工资率和工作周的系数在整个样本中是普遍的，并且系统的个体异质性被个体固定效应所捕获。变量fem是一个表示被调查者性别的二元变量。

本部分完整操作请查看文章：

Stata：面板向量自回归操作教程

完整代码为：

// Example 1: PSID
sjlog using pvar1, replace
webuse psidextract
generate lwks = ln(wks)
pvar lwks lwage if fem == 0, lags(3) 
sjlog close, replace

sjlog using pvar2, replace
pvargranger
sjlog close, replace

sjlog using pvar3, replace
pvarstable
sjlog close, replace

sjlog using pvar4, replace
pvarstable, graph
sjlog close, replace
graph export pvar1.eps, replace

sjlog using pvar5, replace
pvarirf, oirf mc(200) byoption(yrescale) porder(lwage lwks)
sjlog close, replace
graph export pvar2.eps, replace

capture erase fevd_ci.dta
sjlog using pvar6, replace
pvarfevd, mc(200) porder(lwage lwks) save("fevd_ci.dta")
sjlog close, replace

sjlog using pvar7, replace
pvarsoc lwks lwage if fem == 0, pvaropts(instlags(1/4))
sjlog close, replace

sjlog using pvar8, replace
pvar lwks lwage if fem == 0, lags(1) instlags(1/4) 
pvargranger
sjlog close, replace

sjlog using pvar9, replace
xtunitroot ht lwks if fem == 0
xtunitroot ht lwage if fem == 0
sjlog close, replace

sjlog using pvar10, replace
generate gwage = (exp(lwage)-exp(l.lwage))/exp(l.lwage)
generate gwks = (wks - l.wks)/l.wks
xtunitroot ht gwks if fem == 0
xtunitroot ht gwage if fem == 0
pvarsoc gwks gwage if fem == 0, pvaropts(instlags(1/4)) 
pvar gwks gwage if fem == 0, lags(1) instlags(1/4)
pvargranger
sjlog close, replace

案例3：lutkepohl2数据

代码为：

/ Example comparing -pvar- and -var- using time-series data
webuse lutkepohl2
gen id = 1
xtset id qtr

* Estimate VAR model
qui var dln_inv dln_inc dln_consump, lags(1)
  est store var1
  
* Estimate pVAR models  
qui pvar dln_inv dln_inc dln_consump, lags(1)
  est store pvar1_1
qui pvar dln_inv dln_inc dln_consump, lags(1) instl(1/5) gmms
  est store pvar1_5
  
* Comare VAR estimates  
est table var1 pvar1_1 pvar1_5, se stat(N tmin tmax) drop(_cons)

* Create -var- IRF file
est restore var1
irf create varirf, step(10) set(myirf1, replace)

* Compare FEVD
irf table fevd, response(dln_inv) noci
pvarfevd, response(dln_inv)

* Compare IRF - all variables
irf graph oirf, byop(yrescale)
est restore pvar1_1
  pvarirf, byop(yrescale) oirf mc(200)  

* Compare IRF - dln_inv:dln_inv  
irf graph oirf, i(dln_inv) r(dln_inv) title("VAR1") saving(var1_irf, replace)

est restore pvar1_1
pvarirf, imp(dln_inv) res(dln_inv) oirf mc(200) title("PVAR1_1") 
graph save pvar1_1_irf, replace
  
est restore pvar1_5
pvarirf, imp(dln_inv) res(dln_inv) oirf mc(200) title("PVAR1_5") 
graph save pvar1_5_irf, replace

graph combine var1_irf.gph pvar1_1_irf.gph pvar1_5_irf.gph

案例4：grunfeld数据

本部分完整操作请查看文章：Stata：面板VAR模型（pvar2命令）

*# 3.1、导入数据，然后修改变量名称，设定声明

   webuse grunfeld, clear
   rename company id
   xtset id year

**# 3.2、拟合面板VAR

   pvar2 kstock invest mvalue

**# 3.3、拟合面板VAR,滞后2期的

   pvar2 kstock invest mvalue,lag(2)


**# 3.4 面板VAR滞后阶数选择
****    估计一个PVAR(5)模型并选择最优滞后阶数

   pvar2 kstock invest mvalue, lag(5) soc


   
**# 3.5 面板VAR冲击响应函数IRF
   * 用IRFs图进行系统gmm估计(冲击6期)，蒙特卡罗模拟200次(默认值)得到误差边界:

   pvar2 kstock invest mvalue, lag(2) irf(6)
  

**# 3.6 面板VAR方差分解
   pvar2 kstock invest mvalue, irf(5) decomp(10) nograph


**# 3.7 面板VAR格兰杰检验
    pvar2 kstock invest mvalue, lag(3) granger