/* THIS IS FREESHARE GAUSS CODE: THE USUAL CAVEATS APPLY                                        **
** THIS IS A TEMPLATE PROGRAM FILE THAT CALLS THE DPD PROCEDURE FILE                            **
** RESEARCHERS USING THIS CODE SHOULD QUOTE THE ORIGINAL PAPER:                                 **
** HARRIS, M.N. and MATYAS, L. (2004),                                                          **
** "A Comparative Analysis of Different IV and GMM Estimators of Dynamic Panel Data Models",    **    
** International Statistical Review, vol. 72, 3, pp. 397-408.                                   **
*/
 
                       /* GROWTH DPD PROGRAM */

#include c:\mark\gauss\procs\dpd.prc; 
#include c:\mark\gauss\procs\pd_thres.prc;

library co;
coset;
closeall;

cdo = cdir(0);
load path = ^cdo;
save path = ^cdo;

cdo = cdir(0);
load path = ^cdo;
save path = ^cdo;


infile  	= "c:\\mark\\data\\gauss\\fin_bal";
open f1 	= ^infile varindxi;
vars    	= getname(infile);
xx      	= readr(f1,rowsf(f1));
f1      	= close(f1);

xx          = delif(xx,xx[.,iNOT_IN]);

timev   = xx[.,iYEAR];
t       = maxc(timev) - minc(timev) + 1;
tvec	= timev;


    /* DEFINE DEPENDENT VARIABLE AND OTHER VARS */

y       = xx[.,iGROWTH];

"Number of initial observations";;rows(y);


lnlly	    = ln(xx[.,iLLY]);
lnbtot      = ln(xx[.,iBTOT]);
lnpriv      = ln(xx[.,iPRIVO]);
lngdp_i     = xx[.,iLINITIAL];
lngov	    = ln(xx[.,igov]);
lnbmp	    = ln(xx[.,iBMP]+1);
lntrade     = ln(xx[.,iTRADE]);
lnplly      = lnlly .* xx[.,iLPI];
lnpbtot     = lnbtot .* xx[.,iLPI];
lnppriv     = lnpriv .* xx[.,iLPI];
lnpgdp      = lngdp_i .* xx[.,iLPI];
lnpgov      = lngov .* xx[.,iLPI];
lnpbmp      = lnbmp .* xx[.,iLPI];
lnptrade    = lntrade .* xx[.,iLPI];
lnpschl     = xx[.,iSCHOOL] .* xx[.,iLPI];


/*
DEFINE THE FULL SET OF EXPLANATORY VARIABLES
INCLUDING ADDITIONAL INSTRUMENTS (note 1st column = id, next year, then constant)
NOTE ID VECTOR NEED NOT BE CONTIGUOUS
*/

tvec    = xx[.,iYEAR] - minc(xx[.,iYEAR]) + 1;
idvec   = xx[.,iID] - minc(xx[.,iID]) + 1;
T       = maxc(tvec) - minc(tvec) + 1;
N       = maxc(idvec) - minc(idvec) + 1;

x		= ones(rows(xx),1)~xx[.,iLPI iI]~LNGDP_I~lnppriv~lnpriv;
x		= ones(rows(xx),1)~xx[.,iLPI iI]~LNGDP_I~lnplly~lnlly~(xx[.,iYEAR] .eq 1903)~(xx[.,iYEAR] .eq 1904)~
			(xx[.,iYEAR] .eq 1905)~(xx[.,iYEAR] .eq 1906)~(xx[.,iYEAR] .eq 1907);
x		= ones(rows(xx),1)~xx[.,iLPI iI]~LNGDP_I~lnplly~lnlly;

k		= cols(x);
z       = xx[.,iMONEY iM_1] ./ 100;
kz      = cols(z)+1;
xz 		= idvec~tvec~x~z;
xz 		= idvec~tvec~x;
yxztmp	= y~xz;
/*
yxztmp  = delif(yxztmp,yxztmp[.,3] .eq minc(yxztmp[.,3]));      @ DELETE 1ST OBS if USE LAGGED MONEY @ 
*/
yxztmp	= packr(yxztmp);

idvec   = yxztmp[.,2] - minc(yxztmp[.,2]) + 1;
tvec    = yxztmp[.,3] - minc(yxztmp[.,3]) + 1;
T       = maxc(tvec) - minc(tvec) + 1;
N       = maxc(idvec) - minc(idvec) + 1;


/* check for balanced panel */

yxz = zeros(1,cols(yxztmp));

for i (1,N,1);

   	sel_i     	= yxztmp[.,2] .eq i;

    if sumc(sel_i) eq T;
        yxz = yxz|selif(yxztmp,sel_i);
    endif;

endfor;

yxz     = yxz[2:rows(yxz),.];
y       = yxz[.,1];
idvec   = yxz[.,2];
tvec    = yxz[.,3] - minc(yxz[.,3]) + 1;
xz      = yxz[.,4:cols(yxz)];
NT      = rows(yxz);


/* CREATE CONTIGUOUS ID VECTOR */

idvec1   = 1;

for irep (2,NT,1);

    if idvec[irep] .eq idvec[irep-1];
        idvec1 = idvec1|idvec1[irep-1];
    else;
        idvec1 = idvec1|(idvec1[irep-1]+1);
    endif;               

endfor;

idvec   = idvec1;
yxz     = y~idvec~xz;
NT      = rows(yxz);
T       = maxc(tvec) - minc(tvec) + 1;
N       = maxc(idvec) - minc(idvec) + 1;


/* SPLIT SPLINE VARIABLE AS PER OLS1.OUT IE TWO BREAKPOINTS */

gama    = 0.1484~0.2087;
infl    = xz[.,2];
infl_s  = _splinex(infl,gama,0);
xz      = xz[.,1]~xz[.,3:cols(xz)]~infl_s;


"Number of balanced observations";;NT;

            /* NAME X'S AND Z'S & SAVE */


let vnames = GROWTH ID YEAR CONSTANT
    	       LPI I LNGDP_I LNPLLY LNLLY;


let vnames = GROWTH ID TIME CONSTANT 
    	       LPI I LNGDP_I LNPLLY LNLLY M M1;


let vnames = GROWTH ID TIME CONSTANT
    	       I LNGDP_I LNPLLY LNLLY LPI_LO LPI_MED LPI_HI ;


        if not saved(y~idvec~tvec~xz,"fin_dat",vnames);
            errorlog "Disk full";stop;
            end;
        endif;



/* DESCRIPTIVE STATS */

__miss = 2; 
{ vnam,mean,var,std,min,max,valid,mis } = dstat("fin_dat",vnames); 



clear f1, xx;


    /* LINEAR (linear), WANSBEEK-BEKKER (wb) & GMM (gmm) ESTIMATORS? */

let ivtype     = linear;


    /* TIME-VARYING OR TIME INVARIANT?                      */
    /* Supply row vector which corresponds to 1 if          */
    /* XZ variable is time variant, 0 otherwise.            */
    /* Note 1st column of "XZ" is taken to be the constant  */

timev = 0~ones(1,5);

timev = 0~ones(1,k-1)~ones(1,kz-1);

timev = 0~ones(1,rows(vnames)-4);



                    @ 1st element is constant @


    /* EXPLANATORY VARIABLES?                                       */
    /* Supply a row vector corresponding to "XZ" variables above    */
    /* indicating they are explanatory variables (=1), the rest     */
    /* by default will be used as instruments - in levels in the    */
    /* levels model and in differences in the differences model     */


xvars = 1~ones(1,k-1)~zeros(1,kz-1);         @ 1st element is constant @
xvars = 1~ones(1,k-1);         @ 1st element is constant @
xvars = 1~ones(1,rows(vnames)-4);


    /* ENDOGENOUS VARIABLES?                                        */
    /* Supply a row vector corresponding to "XZ" variables above    */
    /* indicating they are explanatory variables AND endogenous     */
    /* 1 indicates explanatory variable AND endogenous              */

endvars = 0~1~zeros(1,k-2)~zeros(1,kz-1);       @ 1st element is constant @
endvars = 0~zeros(1,k-1);                       @ 1st element is constant @
endvars = 0~zeros(1,rows(vnames)-4);


    /* GMM CONDITION 6 (SEE BELOW)                                  */
    /* Supply a row vector corresponding to "ZX" variables above    */
    /* indicating they are to be used in GMM condition 6, a 1       */
    /* indicates they are to be used                                */
    /* Note: only used in GMM estimation and 1st element is the     */
    /* constant and always set = 0                                  */

gmm6 = zeros(1,cols(xz));


    /* GMM CONDITION 7 (SEE BELOW)                              */
    /* Supply a row vector corresponding to "ZX" variables above*/
    /* indicating they are to be used in GMM condition 7, a 1   */
    /* indicates they are to be used                            */
    /* Note: only used in GMM estimation and 1st element is the */
    /* constant and always set = 0                              */

gmm7 = zeros(1,cols(xz));


    /* GMM CONDITIONS                                           */
    /* Supply a row vector corresponding to which GMM           */
    /* conditionsyou want to use, in order:                     */
    /* 1: E(v_it = a_i + u_it)      = 0                         */
    /* 2: V(v_it)                   = sig_a^2 + sig_u^2         */
    /* 3: E(v_it)(v_is)             = sig_a^2                   */
    /* 4: E[(v_it - v_is)(y_i,t-l)] = 0                         */
    /* 5: E(v_it)(x_it)             = 0                         */
    /* note 5 is 1 - endvars                                    */
    /* 6: E(v_it)(x_i0)             = 0                         */
    /* note 6 works in conjunction with gmm6                    */
    /* 7: E(v_it)(x_is)             = 0                         */
    /* note 7 works in conjunction with gmm7                    */

gmm_con = 1~1~1~0~0~0~0;


    /* OPTIMAL GMM STANDARD ERRORS                              */
    /* Supply a scalar as                                       */
    /* 1 = numerical standard errors (recommended)              */
    /* 0 = analytical stadard errors                            */

gmm_se = 1;

    {ols,R2_ols,NTols,
    within,R2_wn,
    fgls,  R2_fgls,
    dols,  R2_dols,NTdols,
    aha,   R2_aha, ahaid, AR2aha, NTah,
    ahb,   R2_ahb, ahbid, AR2ahb, 
    ahc,   R2_ahc, ahcid, AR2ahc, 
    ara,   R2_ara, araid, AR2ara,
    arb,   R2_arb, arbid, AR2arb,
    arc,   R2_arc, arcid, AR2arc, 
    aba,   R2_aba, abaid, AR2aba, 
    abc,   R2_abc, abcid, AR2abc, 
    bb,    R2_bb,  bbid,  AR2bb, 

    wba,   R2_wba, wbaid,
    wbb,   R2_wbb, wbbid,

    gmmi,  R2_gmmi,pwaldi,
    gmm,   R2_gmm, pwald}
        = dpd("fin_dat",ivtype,timev,xvars,endvars,gmm_con,gmm6,gmm7,gmm_se);


/* WRITE RESULTS AS A GAUSS DATA SET */

output file = c:\mark\gauss\output\macro\fin_bal110.out reset;
output on;

selx_lv     = zeros(1,3)~xvars;
selx_df     = zeros(1,3)~(xvars .* timev);
vname_l     = selif(vnames,selx_lv');
if rows(vname_l) gt 2;
    vname_lv    = vname_l[1]|"ldv"|vname_l[2:rows(vname_l)];
else;
    vname_lv    = vname_l[1]|"ldv";
endif;

vname_d    = selif(vnames,selx_df');
if rows(vname_d) ge 1;
    vname_df    = "ldv"|vname_d;
else;
    vname_df    = "ldv";
endif;


format 10,4;
print;
print;
"OLS";
print;
print;
i=1;
do while i le rows(ols);
    $vname_lv[i];;ols[i,.];
i=i+1;
endo;
print;
"Adjusted R^2";;R2_ols;
print;
"NT";;NTols;
print;
print;


print;
print;
"DELTA OLS";
print;
print;
i=1;
do while i le rows(dols);
    $vname_df[i];;dols[i,.];
i=i+1;
endo;
print;
"Adjusted R^2";;R2_dols;
print;
"NT";;NTdols;
print;
print;

print;
print;
"WITHIN";
print;
print;
i=1;
do while i le rows(within);
    $vname_df[i];;within[i,.];
i=i+1;
endo;
print;
"Adjusted R^2";;R2_wn;
print;
"NT";;NTols;
print;
print;

print;
print;
"FGLS";
print;
print;
i=1;
do while i le rows(fgls);
    $vname_lv[i];;ols[i,.];
i=i+1;
endo;
print;
"Adjusted R^2";;R2_fgls;
print;
"NT";;NTols;
print;
print;


let mtype = linear wb gmm;
mtype = indcv(ivtype,mtype);


if mtype eq 1;

    print;
    print;
    "ANDERSON-HSIAO A";
    print;
    print;
    i=1;
    do while i le rows(aha);
        $vname_df[i];;aha[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_aha;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    ahaid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2aha;
    print;
    "NT";;NTah;
    print;
    print;

    print;
    print;
    "ANDERSON-HSIAO B";
    print;
    print;
    i=1;
    do while i le rows(ahb);
        $vname_df[i];;ahb[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_ahb;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    ahbid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2ahb;
    print;
    "NT";;NTah;
    print;
    print;
    
    print;
    print;
    "ANDERSON-HSIAO C";
    print;
    print;
    i=1;
    do while i le rows(ahc);
        $vname_df[i];;ahc[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_ahc;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    ahcid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2ahc;
    print;
   "NT";;NTah;
    print;
    print;

    print;
    print;
    "ARELLANO A";
    print;
    print;
    i=1;
    do while i le rows(ara);
        $vname_df[i];;ara[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_ara;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    araid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2ara;
    print;
    "NT";;NTdols;
    print;
    print;

    print;
    print;
    "ARELLANO B";
    print;
    print;
    i=1;
    do while i le rows(arb);
        $vname_df[i];;arb[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_arb;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    arbid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2arb;
    print;
    "NT";;NTdols;
    print;
    print;

    print;
    print;
    "ARELLANO C";
    print;
    print;
    i=1;
    do while i le rows(arc);
        $vname_df[i];;arc[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_arc;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    arcid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2arc;
    print;
    "NT";;NTdols;
    print;
    print;
    
    print;
    print;
    "ARELLANO-BOND A";
    print;
    print;
    i=1;
    do while i le rows(aba);
        $vname_df[i];;aba[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_aba;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    abaid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2aba;
    print;
    "NT";;NTdols;
    print;
    print;

    print;
    print;
    "ARELLANO-BOND C";
    print;
    print;
    i=1;
    do while i le rows(abc);
        $vname_df[i];;abc[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_abc;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    abcid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2abc;
    print;
    "NT";;NTdols;
    print;
    print;

    print;
    print;
    "BLUNDELL-BOND";
    print;
    print;
    i=1;
    do while i le rows(bb);
        $vname_lv[i];;bb[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_bb;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    bbid;
    print;
    "M2 Test for AR(2) Residuals: Standard Normal";
    print;
    AR2bb;
    print;
    "NT";;NTdols*2;
    print;
    print;

elseif mtype eq 2;

    print;
    print;
    "WANBEEK-BEKKER A";
    print;
    print;
    i=1;
    do while i le rows(wba);
        $vname_df[i];;wba[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_wba;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    wbaid;
    print;
    "NT";;NTols;
    print;
    print;

    print;
    print;
    "WANBEEK-BEKKER B";
    print;
    print;
    i=1;
    do while i le rows(wbb);
        $vname_df[i];;wbb[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_wbb;
    print;
    "Sargan Test for Over-Identifying restrictions (Chi^2)";
    print;
    "     test        df's   P-value";
    print;
    wbbid;
    print;
    "NT";;NTols;
    print;
    print;

else;

    print;
    print;
    "GMM INDENTITY";
    print;
    print;
    i=1;
    do while i le rows(gmmi);
        $vname_lv[i];;gmmi[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_gmmi;
    print;
    "NT";;NTols;
    print;
    print;

    print;
    print;
    "OPTIMAL GMM";
    print;
    print;
    i=1;
    do while i le rows(gmm);
        $vname_lv[i];;gmm[i,.];
    i=i+1;
    endo;
    print;
    "Adjusted R^2";;R2_gmm;
    print;
    "P-value of Hausman Test for Over-Identifying Restrictions";pwald;
    print;
    "Conditions used:";gmm_con';
    print;
    "Endogenous X vector used:";endvars';
    print;
    "NT";;NTols;
    print;
    print;

endif;


output off;
stop;


