options crt memory=20 ;
? ==================================================
? [Dynamic Panel Data Estimation]
? [Using GMM(DIF) or GMM (SYS)]
? [without Time Dummies or with Time Dummies]
? ==================================================
? We can select the options to estimate
? dynamic panel data model in different modes.
? ==================================================
? Designed by Yoshitsugu Kitazawa
? First Design: 2003.12.1 Last Update: 2003.12.6
? E-mail: kitazawa@ip.kyusan-u.ac.jp
? ==================================================
? ** This script realizes almost the same results
? ** as Ox version DPD except for the serial correlation
? ** tests.
? ** This script is written on the accumulation of
? ** the past technical supports by
? ** Professor Bronwyn Hall and Dr. Clint Cummins.
? ==================================================
? * T=6:  All time periods used for the estimation
? * N=47: Number of individuals
? * t = 1 to T and i = 1 to N
?
? * Estimation Model without and with Time Dummy
? For all i (individual) and t (time), 
?    y_{it} = alp * y_{i,t-1} + bet * x_{it} + (TD_t) + u_{it}
? and
?    u_{it} = eta_i + v_{it}
? where y_{it} is the dependent variable, x_{it} is
? the explanatory variable, alp and beta are the 
? parameter of interest to be estimated, eta_i is
? the fixed effect (the nuisance parameter), and
? v_{it} is the disturbance. In this model, 
? E[eta_i x_{it}] =\ 0 for all T, and x_{it} is
? endogenously determined so that E[v_{it} x_{is}]=0
? for s=1 to t-1.
? Further, it is supposed that y_{it} and x_{it} are
? mean-stationary. In addition, (TD_t)s are the
? time dummies to be estimated, when (TD_t)s are
? specified in the model.
?
? * Specify the options at Specification of Options below.
? lt = 6 : nnumber of all time periods used in the
?          estimation in this case.
? 
? * Specify the options at Specification of Options below.
? [1] tm=0 and om=1 : without Time Dummies and using GMM(DIF)
? [2] tm=0 and om=2 : without Time Dummies and using GMM(SYS)
? [3] tm=1 and om=1 : with Time Dummies and using GMM(DIF)
? [4] tm=1 and om=2 : with Time Dummies and using GMM(SYS)
?
? * References
? Arellano (2003): "Panel Data Econometrics", Oxford
? University Press, Great Britain.
? Arellano and Bond (1991): "Some test of specification
? for panel data: Monte Carlo evidence and an application
? to employment equations", Review of Economic Studies,
? 58, 277-297.
? Arellano and Bover (1995): "Another look at the
? instrumental variables estimation of error-components
? models", Journal of Econometrics, 68, 29-51.
? Blundell and Bond (1998): "Initial conditions and
? moment restrictions in dynamic panel data models",
? Journal of Econometrics, 87, 115-143.
? Doornik, Arellano, and Bond (2002): "Panel Data
? estimation using DPD for Ox",
? http://www.nuff.ox.ac.uk/Users/Doornik/papers/dpd.pdf
?
? * See on detail
? http://www.ip.kyusan-u.ac.jp/J/keizai/kitazawa/SOFT/TSP_DPD1/index.htm
? =====================================================

input "cov1spdm.tsp" ;
input "cov1ssys.tsp" ;
input "cov2stpm.tsp" ;
input "evdgen.tsp" ;
input "lm2test2.tsp" ;
input "mattrim1.tsp" ;


? ===== Specification of Options ===================
? --------------------------------------------------
? *** Specify here [the number of time periods],
? *** [the usage of time dummies], and
? *** [the selection of estimation method].
? --------------------------------------------------
? Number of all time periods used in the estimation
set lt = 6 ;
? Time Dummy (TD) Specification
set tm = 1 ;	? 0: Without TD,	1: With TD
? Selection of estimation method
set om = 2 ;	? 1: GMM(DIF),		2: GMM(SYS)
? ==================================================

print lt tm om ;

if(.not. (tm=0 .or. tm=1)) ; then ; do ;
  title "Error: tm must be 0 or 1" ;
  stop ;
enddo ;

if(.not. (om=1 .or. om=2)) ; then ; do ;
  title "Error: om must be 1 or 2" ;
  stop ;
enddo ;

set ltm1 = lt - 1 ;
set ltm2 = lt - 2 ;


? Level Eq.
list(first=2, last=lt) tmplis ;
dot(index=i) tmplis ;
  set i=i+1 ;
  frml ux.  y. - ( alp * y_m1 + bet * x. + tdx. ) ;
  set j=i-1 ; list(first=j, last=j) jlis ; dot jlis ;
    frml y_m1 y. ; enddot ;
  frml tdx. td. ;
  eqsub(name=ut., print) ux. y_m1 tdx. ;
  frml tdx. 0 ;
  eqsub(name=un., print) ux. y_m1 tdx. ;
enddot ;


? Differenced Eq.
list(first=3, last=lt) tmplis ;
dot(index=i) tmplis ;
  set i=i+2 ;
  ? 
  frml dux.  ut. - un_m1 ;
  set j=i-1 ; list(first=j, last=j) jlis ; dot jlis ;
    frml un_m1 un. ; eqsub un_m1 un. ; enddot ;
  eqsub(name=dutn., print) dux. ut. un_m1 ;
  ? 
  frml dux.  ut. - ut_m1 ;
  set j=i-1 ; list(first=j, last=j) jlis ; dot jlis ;
    frml ut_m1 ut. ; eqsub ut_m1 ut. ; enddot ;
  eqsub(name=dut., print) dux. ut. ut_m1 ;
  ? 
  frml dux.  un. - un_m1 ;
  set j=i-1 ; list(first=j, last=j) jlis ; dot jlis ;
    frml un_m1 un. ; eqsub un_m1 un. ; enddot ;
  eqsub(name=dun., print) dux. un. un_m1 ;
enddot ;

title "Equations used for the estimation" ;

if(om=1 .and. tm=0) ; then ; do ;
  frml du2 0 ;
  list(first=3, last=lt) tmplis ;
  dot tmplis ;
    frml du. dun. ;
    eqsub(print) du. dun. ;
  enddot ;
enddo ;

if(om=1 .and. tm=1) ; then ; do ;
  frml du2 0 ;
  list(first=3, last=lt) tmplis ;
  dot tmplis ;
    frml du. dutn. ;
    eqsub(print) du. dutn. ;
  enddot ;
enddo ;

if(om=2 .and. tm=0) ; then ; do ;
  frml du2 0 ;
  list(first=3, last=lt) tmplis ;
  dot tmplis ;
    frml du. dun. ;
    eqsub(print) du. dun. ;
  enddot ;
  list(first=2, last=lt) tmplis ;
  dot tmplis ;
    frml u. un. ;
    eqsub(print) u. un. ;
  enddot ;
enddo ;

if(om=2 .and. tm=1) ; then ; do ;
  frml du2 0 ;
  list(first=3, last=lt) tmplis ;
  dot tmplis ;
    frml du. dut. ;
    eqsub(print) du. dut. ;
  enddot ;
  list(first=2, last=lt) tmplis ;
  dot tmplis ;
    frml u. ut. ;
    eqsub(print) u. ut. ;
  enddot ;
enddo ;


list(first=lt, last=lt) tmplis ;
dot tmplis ;
  list eqlisa du2-du. ;
  list eqlisb u2-u. ;
enddot ;

if(om=1) ; then ;
  do ; list eqlis  eqlisa ; enddo ;
if(om=2) ; then ;
  do ; list eqlis  eqlisa eqlisb ; enddo ;


list(first=ltm2, last=ltm2) tmplis ;
dot tmplis ; list insta  y1-y. x1-x. ; enddot ;
list(first=ltm1, last=ltm1) tmplis ;
dot tmplis ; list instb  dy2-dy. dx2-dx. ; enddot ;

if(om=1) ; then ;
  do ; list zlis  insta c ; enddo ;
if(om=2) ; then ;
  do ; list zlis  insta instb c ; enddo ;


mform(type=triang, ncol=ltm2) one22ut = 1 ;
mat one22id = ident(ltm2) ;
mform(type=sym, ncol=ltm2) zer22 = 0 ;
mform(nrow=1, ncol=ltm2) zer12 = 0 ;
mform(nrow=1, ncol=ltm2) one12 = 1 ;
mat zer12t = zer12' ;
mform(nrow=1, ncol=1) zer11 = 0 ;
mform(nrow=1, ncol=1) one11 = 1 ;

if(om=1 .and. tm=0) ; then ; do ;
  mmake(vert) zmaf zer12t  zer12t  zer11 ;  
  mmake(vert) zma  one22ut one22ut zer12 ;
  mmake zmask zmaf zma ;
enddo ;
if(om=1 .and. tm=1) ; then ; do ;
  mmake(vert) zmaf zer12t  zer12t  zer11 ;  
  mmake(vert) zma  one22ut one22ut one12 ;
  mmake zmask zmaf zma ;
enddo ;
if(om=2 .and. tm=0) ; then ; do ;
  mmake(vert) zmaf zer12t  zer12t  zer12t  zer12t  zer11 ;  
  mmake(vert) zma  one22ut one22ut zer22   zer22   zer12 ;
  mmake(vert) zmbf zer12t  zer12t  zer12t  zer12t  zer11 ;
  mmake(vert) zmb  zer22   zer22   one22id one22id zer12 ;
  mmake zmask zmaf zma zmbf zmb ;
enddo ;
if(om=2 .and. tm=1) ; then ; do ;
  mmake(vert) zmaf zer12t  zer12t  zer12t  zer12t  zer11 ;  
  mmake(vert) zma  one22ut one22ut zer22   zer22   zer12 ;
  mmake(vert) zmbf zer12t  zer12t  zer12t  zer12t  one11 ;
  mmake(vert) zmb  zer22   zer22   one22id one22id one12 ;
  mmake zmask zmaf zma zmbf zmb ;
enddo ;

print zmask ;

mat neqx = ncol(zmask) ;


if(om=1 .and. tm=1) ; then ; do ;
  list(first=3, last=lt) tmplis ;
  dot tmplis ;  frml zr_td.  td. - 0 ;  enddot ;
  list(first=lt, last=lt) tmplis ;
  dot tmplis ;  list zr_tds  zr_td3-zr_td. ;  enddot ;
enddo ;
if(om=2 .and. tm=1) ; then ; do ;
  list(first=2, last=lt) tmplis ;
  dot tmplis ;  frml zr_td.  td. - 0 ;  enddot ;
  list(first=lt, last=lt) tmplis ;
  dot tmplis ;  list zr_tds  zr_td2-zr_td. ;  enddot ;
enddo ;


freq n ;
smpl 1 47 ;
read(file="ydata.txt") code y1-y6 ;
read(file="xdata.txt") code x1-x6 ;

? print y1-y6 ;
? print x1-x6 ;

genr y_lag1 = y1 ;
list(first=2, last=lt) tmplis ;
dot tmplis ;
  genr dy. = y. - y_lag1 ;
  genr y_lag1 = y. ;
enddot ;

genr x_lag1 = x1 ;
list(first=2, last=lt) tmplis ;
dot tmplis ;
  genr dx. = x. - x_lag1 ;
  genr x_lag1 = x. ;
enddot ;


param alp bet ;
if(om=1 .and. tm=1) ; then ; do ;
  list(first=lt, last=lt) tmplis ;
  dot tmplis ; param td3-td. ; enddot ;
enddo ;
if(om=2 .and. tm=1) ; then ; do ;
  list(first=lt, last=lt) tmplis ;
  dot tmplis ; param td2-td. ; enddot ;
enddo ;

title "<One step estimation>" ;
if(om=1) ; then ; do ;
  cov1spdm neqx zlis zmask ani1 ;
enddo ;
if(om=2) ; then ; do ;
  cov1ssys neqx zlis zmask ani1 ;
enddo ;
mat ani1 = ani1 / 2 ;
gmm(covoc=ani1, noiteroc, nooptcov, hetero, nma=0,
    maxit=100, silent, inst=zlis, mask=zmask) eqlis ;
print alp bet ;
if(tm=1) ; then ; do ;
  if(om=1) ; then ; do ;
    list(first=lt, last=lt) tmplis ;
    dot tmplis ;  print td3-td. ;  enddot ;
  enddo ;
  if(om=2) ; then ; do ;
    list(first=lt, last=lt) tmplis ;
    dot tmplis ;  print td2-td. ;  enddot ;
  enddo ;
enddo ;
tstats(name=@rnms) @coef @vcov ;
if(tm=1) ; then ; do ;
  title "Wald Test: H_0: Time Dummies are all zero" ;
  analyz zr_tds ;
enddo ;
title "Sargan Test for 1 step estimation" ;
copy @res restmp ;
mat lc = ncol(@res) ;
if(om=1) ; then ; do ;
  set fc = 2 ;
enddo ;
if(om=2) ; then ; do ;
  set fc = lc / 2 + 1 ;
enddo ;
mattrim1 @res fc lc @res ;
evdgen sr1 ;
copy restmp @res ;
set sargan = @gmmovid / sr1 ;
print sargan ; cdf(chi, df=@novid) sargan ;
title "Serial Correlation Tests, AR(1) and AR(2)" ;
copy @res restmp ;
set fc = 2 ;
mat nc = ncol(@res) ;
if(om=1) ; then ; do ;
  set lc = nc ;
enddo ;
if(om=2) ; then ; do ;
  set lc = nc / 2 ;
enddo ;
mattrim1 @res fc lc @res ;
regopt(noprint) all ;  lm2test2 ;  regopt ;
copy restmp @res ;

title "<Two step estimation>" ;
regopt(noprint) all ;
cov2stpm eqlis zlis zmask ani2 ;
regopt ;
gmm(covoc=ani2, noiteroc, optcov, hetero, nma=0,
    maxit=100, silent, inst=zlis, mask=zmask) eqlis ;
print alp bet ;
if(tm=1) ; then ; do ;
  if(om=1) ; then ; do ;
    list(first=lt, last=lt) tmplis ;
    dot tmplis ;  print td3-td. ;  enddot ;
  enddo ;
  if(om=2) ; then ; do ;
    list(first=lt, last=lt) tmplis ;
    dot tmplis ;  print td2-td. ;  enddot ;
  enddo ;
enddo ;
tstats(name=@rnms) @coef @vcov ;
if(tm=1) ; then ; do ;
  title "Wald Test: H_0: Time Dummies are all zero" ;
  analyz zr_tds ;
enddo ;
title "Sargan Test for 2 step estimation" ;
set sargan = @gmmovid ;
print sargan ; cdf(chi, df=@novid) sargan ;
title "Serial Correlation Tests, AR(1) and AR(2)" ;
copy @res restmp ;
set fc = 2 ;
mat nc = ncol(@res) ;
if(om=1) ; then ; do ;
  set lc = nc ;
enddo ;
if(om=2) ; then ; do ;
  set lc = nc / 2 ;
enddo ;
mattrim1 @res fc lc @res ;
regopt(noprint) all ;  lm2test2 ;  regopt ;
copy restmp @res ;

stop ;
end ;