@ GMMMHF.EXP @ @ Another example file for GMM.SET to illustrate a more complicated application. See GMM.EXP for a simpler example.@ @ This is a minor modification of a program used by Cooley and Ogaki (1993), "A Time Series Analysis of Real Wages, Consumption, and Asset Returns: A Cointegration-Euler Equation Approach," Rocheseter Center for Economic Research Working Paper No. 285R, University of Rochester, for a asset pricing model with habit formation. See "GMM: A User's Guide" Section 5.A. @ @ To run this example file, type RUN MINQUAD.SET [F4][F2] RUN GMM.SET [F4][F2] RUN GMMHF.EXP [F4][F2] When you run this program with the default hflag value of 2 (see MINQUAD.SET for explanations of hflag), the initial idendity weighing matrix estimation will not converge. After about 30 iterations, hit 1 on your key board. Then the MINQUAD program will change the hflag value to 1 to use the DFP update method, and you will get a convergence. After that, hit 2 on your key board to set hflag to 2 again and use the outer products method, which is preferable for most cases. @ #include "gradqg.prc"; @ ======================== User Definition Area ============================ @ @ ********************************************************************** @ @ ----------- Step 1: Prepare Output File ---------------------------@ @ ********************************************************************** @ output file=gmmq.out reset; @ Specify the name of the output file. @ datestr(0); @ print date using a procedure in TIME.ARC @ timestr(0); @ print time using a procedure in TIME.ARC @ ?"GMMHF.OUT Another example for GMM.SET NDS, VWR and T-Bill, Parameter Names are 'delta a1 alpha"; @ *********************************************************************** @ @------------ Step 2: Define Global Variables ------------------------- @ @ *********************************************************************** @ @--------------------------------- Specify the following globals --------------------------------- @ nob=176; @nob=176 47:1-90:4@ izlag=1; taubeg=1+izlag; tauend=nob-2; ? "izlag,taubeg,tauend=" izlag~taubeg~tauend; tend=tauend-taubeg+1; bgm=0.99|0|1; @ bgm=b(1)|..|b(kgm); Initial values of the coefficients @ kgm=rows(bgm); @ # of parameters @ ndbt=2; @ # of disturbance terms in w(t) @ nmr=ndbt*(1+2*izlag); mas=1; @ The order of MA of the disturbance @ @---------------------------------@ @ Section B. User can leave the variables in this section as they are, and go to STEP 3. Only advanced users should modify this section.@ @---------------------------------@ gradname=&GRADQG; @ Specify the name of proc that calculates the gradient dgT(b)/db. @ const=1/sqrt(tend); @ Scaling Multiplier when W0=eye(L) @ const2=1/sqrt(tend); @ Scaing Multiplier when W0 is not eye(L) @ @ These scaling multipliers should be set so that the value of function (vof) in nonlinear search of MINQUAD.SET is near 1. If vof is too close to 0, the search will not work properly. const2=1/sqrt(tend) will generally be good. @ w0flag=0; @ scalar; If w0flag=0, W0=I is used as the initial weighing matrix W0. If w0flag=1, initial bgm is used to calculate initial W0. If w0flag=2, W0 in the memory is used as initial W0. If w0flag=3, W0 and bgm in the memory are used to give the first GMM result. @ maxitegm=5; @ Sets maximum # of iteration over weighting matrix, W0. Set w0flag=0 and maxitegm=2 to execute usual 2-Stage GMM. @ zero=1E-2; @ Iteration over W0 continues until the maximum difference of the current and the previous W0 in absolute value becomes less than 'zero', or the # of iteration exceeds maxitegm. @ calwflag=0; @ This variable is used to choose the method to calculate the distance matrics, W0. If calwflag=0, Durbin's method will be used when W0 is singular. If calwflag=1, the QS kernel estimator (nonprewhitened or prewhitened) will be used when W0 is singular. If calwflag=3, Durbin's method will be used. If calwflag=4, the QS kernel estimator (nonprewhitened or prewhitened) will be used. @ ordard=mas+3; @ Order of AR representation for Durbin's method @ @The following four variables control the QS kernel estimator.@ st=0; @ A scalar to control the bandwidth parameter for the QS kernel. if st=0, then an automatic bandwidth estimator is used. if st/=0, then st is used as the bandwidth parameter.@ wav=ones(nmr,1); @nmr by 1 vector; weights given to the a-th element of z(t)w(t) for the automatic bandwidth estimator (used only if st=0). @ maxd=0; @ scalar: if maxd=0, then nonprewhitened HAC with the QS kernel is used. if maxd>0, then the elements of DeltaLS with the absolute value greater than maxd is replaced by maxd. See Andrews and Monahan's (1990) footnote 4. @ bst=10e+5; @ scalar: When automatic bandwidth parameter is calculated to be bigger than bst, bst is used. @ @ See MINQUAD.SET for the following globals @ hflag=2; dfpflag=1; sstol=1e-25; @ See MAXMUM.DOC on MODULE9 of GAUSS for the following globals @ gradtol=1e-5; btol=1e-5; typf=1; typb=1; @ *************************************************************** @ @ ----------- Step 3: READING IN DATA --------------------------- @ @ *************************************************************** @ load nw[178,1]=qwages92.dat; load c1n[176,2]=qnrnd91.dat; load c1s[176,2]=qnrs91.dat; load mpopl[541,2]=mpop92.dat; load mret[780,6]=mret.dat; load nrfr[781,11]=rskfr91.dat; nrfr=exp(nrfr[251:781,8].*nrfr[251:781.,10]/36500); @exp(r*NDM/365), gross return@ nrfr=reshape(nrfr,177,3); nrfr=nrfr[.,3]; vwrv=mret[253:780,2]; @vwrv ex post monthly retrun, 253 is 1947 Jan.@ vwrv=1+vwrv; @gross retrun@ vwrv=reshape(vwrv,176,3); vwr=vwrv[.,1].*vwrv[.,2].*vwrv[.,3]; @gross quarterly return@ mpopl=mpopl/1000000; c1v=c1n+c1s; @ NDS@ @ 2. Transform the data if necessary @ @ real per capita consumption @ mpopl=mpopl[.,2]; @mpopl[.,2]; 16+@ popul=reshape(mpopl[1:nob*3],nob,3); popul=meanc(popul'); @ave. over each quarter@ c1=c1v[1:nob,2]./popul[1:nob]; p1=c1v[1:nob,1]./c1v[1:nob,2]; c1str=c1[2:nob,.]./c1[1:nob-1,.]; c2str=c1[3:nob,.]./c1[1:nob-2,.]; cm1str=c1[1:nob-1,.]./c1[2:nob,.]; vwr=vwr[2:nob,.].*p1[1:nob-1,.]./p1[2:nob,.]; rfr=nrfr[2:nob,.].*p1[1:nob-1,.]./p1[2:nob,.]; clear popul; @Conventional Instruments@ zpv=ones(tend,1)~c1str[taubeg-1:tauend-1,.]~vwr[taubeg-1:tauend-1,.]; i=1; do until i>(izlag-1); zpv=zpv~c1str[(taubeg-1-i):(tauend-1-i),.]~ vwr[(taubeg-1-i):(tauend-1-i),.]; i=i+1; endo; zpr=ones(tend,1)~c1str[taubeg-1:tauend-1,.]~rfr[taubeg-1:tauend-1,.]; i=1; do until i>(izlag-1); zpr=zpr~c1str[(taubeg-1-i):(tauend-1-i),.]~ rfr[(taubeg-1-i):(tauend-1-i),.]; i=i+1; endo; rfr=rfr[taubeg:tauend,.]; vwr=vwr[taubeg:tauend,.]; maxcr=maxc(-c1[taubeg:tauend+2]./c1[taubeg-1:tauend+1]); lpm=(1e+20)*ones(tend,nmr); proc hu(b); @ scale by 1+b[1]*b[2] for identification@ local wv,wr,s1str,s2str,hr,bca; bca=b[3]; if b[2]<=maxcr; hr=lpm; goto lhue; endif; s1str=((c1[taubeg+1:tauend+1]+b[2]*c1[taubeg:tauend])./(c1[taubeg:tauend]+ b[2]*c1[taubeg-1:tauend-1]))^(-bca); s2str=((c1[taubeg+2:tauend+2]+b[2]*c1[taubeg+1:tauend+1])./( c1[taubeg:tauend]+b[2]*c1[taubeg-1:tauend-1]))^(-bca); wv=(b[1]*(s1str+b[1]*b[2]*s2str).*vwr-1-b[1]*b[2]*s1str)/(1+b[1]*b[2]); wr=(b[1]*(s1str+b[1]*b[2]*s2str).*rfr-1-b[1]*b[2]*s1str)/(1+b[1]*b[2]); hr=((zpv.*wv)~(zpr.*wr)); if b[2]>1; hr=hr*((b[2]-1)^2+1); endif; if b[3]>10; hr=hr*((b[3]-5)^2+1); endif; lhue: retp(hr); endp; @=========== The User does not have to change the code below. ============== @ @--------------------- Print results of estimation---------------------@ proc (0)=prntrslt(x,s); ? " --------------- Minimization Results ---------------------- Step size: " s " Value of the objective function: " vof " Minimiser is" x' " -----------------------------------------------------------"; endp; @ ************************************************************************ @ @ THE PROGRAM STARTS @ chi=gmm(gradname,tend,nmr,mas,rows(bgm),zero,maxitegm,st,wav,maxd,bst); output off; @ ------------------------ End of Program ------------------------------ @