******************************************************************************** * THESE ARE THE PROGRAMS FOR MICHAEL DUEKER'S MAY/JUNE 1993 REVIEW ARTICLE. * * * ******************************************************************************** *tvpmrkf.prg* @>>>>>>>>>>>>>>>>>>>>Time varying parameter model >>>>>>>>>>>@ @>>>>>>>>>>>>>>>>>>>>with Markov Varance >>>>>>>>>>>>>>>>>>>@ @>>>>>>>>>>>>>>>>>>>>Using Optmum Proc >>>>>>>>>>>>>>>>>>>>>>@ @===========> GAUSS 2.1 + NEW APPLICATION MODULE@ new; #lineson; library optmum; #include c:\gauss\src\optmum.ext; optset; output file=c:\letter\tvpmrkf.t05 reset; load data[134,12]=c:\dat\veldata.dat; load data2[134,10]=c:\dat\newduke2.dat; m2=data2[.,4]; clear data2; /* column: 1:quarter index,2:ind. prod. 3:cpi 4:ngdp, 5:ngnp 6:gdpdef 7:gnpdef 8:m1 9:stlb 10:bogb 11:tb3 12: caputil */ t=rows(data); @ vstl= ngdp/stlb @ vstl= data[.,4]./m2; qtr_idx= data[3:t,1]; y= 100*(ln(vstl[3:t])-ln(vstl[2:t-1])); dr3= (data[2:t-1,11]-data[1:t-2,11])./data[1:t-2,11]; gmb= 100*(ln(m2[2:t-1])-ln(m2[1:t-2])); t=t-2; /* load data[114,6]=d:\gmm\kimtvp.dat; column 1: Quarter Index 2: m1===growth rate of quarterly average M1 3: dint=change in the lagged interest rate (3-month T-bill) 4: inf==lagged inflation 5: surpl==lagged full employment budget surplus 6: m1lag==lag of m1 1959.3--1987.4, ........................ qtr_idx=data[.,1]; data=data[.,2:6]; T=rows(data); */ AA=ones(t,1)~dr3~gmb; yy=y; dddd=AA; beta=invpd(AA'AA)*AA'yy; res= yy-AA*beta; sig=res'res/(t-1); beta; sig; my=yy-sumc(yy)/rows(yy); r2= 1 - res'res/(my'my); r2;; r2; @ m1, intcpt, dint, inf, surp,m1lag@ format /m1 /rd 12,6; y; @================= Initialize Global Variables============@ START=11; @1961.1--1987.4 @ PRMTR_NO=7; @ >>Number of parameters to be estimated. >>@ TVP_NO=3; @ >>Number of Time-varying Coefficients >>@ /* prmtr_in=0.323289|0.090582|0.032839|0.040333|0.063470|0.01|0.970468| 0.928421|0.958763; output on; "initial values of prmtr is(vre,intcpt,dint,inf,surp,mlag,w,p,q"; */ output on; prmtr_in=.4*sig|.2|2E-2|2E-3|1.5*sig|.90|.95; prmtr_in; PRIOR_CF=zeros(TVP_NO,1); @Initial values of reg coefs@ PRIOR_CF=inv(dddd[1:10,.]'dddd[1:10,.])*dddd[1:10,.]'yy[1:10]; PRIOR_VR=EYE(TVP_NO)*START; @Initial var-cov for reg coef @ @ Maximum Likelihood Estimation @ @==================================================@ {xout,fout,gout,hout}=optmum(&lik_fcn,tran_inv(PRMTR_IN)); @ prmtr estimates, -log lik value, Grandient, code@ PRM_FNL=TRANS(xout); F_SS =SAVE_OUT(PRM_FNL); OUTPUT FILE=c:\letter\TVPMRKF.D05 RESET; F_SS; "Calculating Hessian and its invers.....";; "THANK YOU FOR YOUR PATIENCE..............!!!"; hout=inv(hessp(&lik_fcn,xout)); grdn_fnl=gradfd(&TRANS,xout); Hsn_fnl=grdn_fnl*hout*grdn_fnl'; SD_fnl =sqrt(diag(Hsn_fnl)); @Standard errors of the estimated coefficients@ output on; "==FINAL OUTPUT========================================================"; "likelihood value is ";; fout; "Estimated parameters are:"; prm_fnl'; "var-cov matrix is:"; Hsn_fnl; "Standard errors of parameters are:"; sd_fnl; "==============================================================="; END; @========================================================================@ @========================================================================@ PROC LIK_FCN(PRMTR1); local prmtr, ppr,qpr,aaa,prob__0,prob__1,p_cf_0,p_cf_1, p_vr_0, p_vr_1, qq, lik, j_iter, h, pr__0_l,pr__1_l, pst_cf0, pst_cf1, pst_vr0, pst_vr1, f_cast0, f_cast1, ss00, ss01, ss10, ss11, pr_vl00, pr_vl01, pr_vl10, pr_vl11, pr_val,k_gn00, k_gn01, k_gn10, k_gn11, p_cf00, p_cf01, p_cf10, p_cf11, p_vr00, p_vr01, p_vr10, p_vr11, pro_00, pro_01, pro_10, pro_11, likv; EXTERNAL PROC TRANS, V_PROB; PRMTR=TRANS(PRMTR1); PPR=PRMTR[6,1]; @Pr[St=1/St-1=1]@ QPR=PRMTR[7,1]; @Pr[St=0/St-1=0]@ AAA=EYE(TVP_NO); @>>>>>>>>>>>>>>>>>>>>>>>>> INITIAL PROB. Pr[S0/Y0] @ PROB__1=(1-QPR)/(2-PPR-QPR); @Pr[St-1=1/Yt-1], STEADY STATE PROB.@ PROB__0=1-PROB__1 ; @Pr[ST-1=0/Yt-1], STEADY STATE PROB @ P_CF_0=PRIOR_CF; P_CF_1=PRIOR_CF; P_VR_0=PRIOR_VR; P_VR_1=PRIOR_VR; @ Initial values@ QQ=(PRMTR[2,1]^2~0~0)| (0~PRMTR[3,1]^2~0)| (0~0~PRMTR[4,1]^2); LIKV=0.0; J_ITER = 1; DO UNTIL J_ITER>T; @H=AA[J_ITER,.]~zz[J_ITER,.];@ H=AA[J_ITER,.]; PR__0_L= QPR*PROB__0 + (1-PPR)*PROB__1; @Pr[St=0/Yt-1]@ PR__1_L= (1-QPR)*PROB__0 + PPR*PROB__1; @Pr[St=1/Yt-1]@ @===================PREDICTION=============================@ PST_CF0 = AAA * P_CF_0; @ Prediction of coef when St-1=0 @ PST_CF1=AAA*P_CF_1; PST_VR0 = AAA * P_VR_0 * AAA' +QQ; PST_VR1 = AAA * P_VR_1 * AAA' +QQ; F_CAST0 = yy[J_ITER,.] - H * PST_CF0; @ Forecast error when St-1=0 @ F_CAST1 = yy[J_ITER,.] - H * PST_CF1; @ Forecast error when St-1=0 @ SS00= H * PST_VR0 * H' + PRMTR[1,1]^2; @ St-1=0,St=0 @ SS01= H* PST_VR0 * H' + PRMTR[5,1]^2; @ St-1=0,St=1 @ SS10= H * PST_VR1 * H' + PRMTR[1,1]^2; @ St-1=1,St=0 @ SS11= H * PST_VR1 * H' + PRMTR[5,1]^2; @ St-1=1,St=1 @ @<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<@ PR_VL00=V_PROB(F_CAST0,SS00)* QPR*PROB__0; @Pr[St,Yt/Yt-1]@ PR_VL01=V_PROB(F_CAST0,SS01)*(1-QPR)*PROB__0; PR_VL10=V_PROB(F_CAST1,SS10)*(1-PPR)*PROB__1; PR_VL11=V_PROB(F_CAST1,SS11)* PPR*PROB__1; PR_VAL=PR_VL00+PR_VL01+PR_VL10+PR_VL11; @Pr[Yt/Yt-1]@ @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ LIK=-1*LN(PR_VAL); @=======================UPDATING===============================@ K_GN00 =PST_VR0 * H' / SS00; K_GN01 =PST_VR0 * H' / SS01; K_GN10 =PST_VR1 * H' / SS10; K_GN11 =PST_VR1 * H' / SS11; P_CF00 = PST_CF0 + K_GN00 * F_CAST0; P_CF01 = PST_CF0 + K_GN01 * F_CAST0; P_CF10 = PST_CF1 + K_GN10 * F_CAST1; P_CF11 = PST_CF1 + K_GN11 * F_CAST1; P_VR00 = (EYE(TVP_NO) - K_GN00 * H ) * PST_VR0; P_VR01 = (EYE(TVP_NO) - K_GN01 * H ) * PST_VR0; P_VR10 = (EYE(TVP_NO) - K_GN10 * H ) * PST_VR1; P_VR11 = (EYE(TVP_NO) - K_GN11 * H ) * PST_VR1; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ PRO_00=PR_VL00/PR_VAL; @Pr[St,St-1/Yt]@ PRO_01=PR_VL01/PR_VAL; PRO_10=PR_VL10/PR_VAL; PRO_11=PR_VL11/PR_VAL; PROB__0=PRO_00+PRO_10; @Pr[St=0/Yt]@ PROB__1=PRO_01+PRO_11; @Pr[St=1/Yt]@ @<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<@ P_CF_0=(PRO_00*P_CF00 + PRO_10*P_CF10)/PROB__0; P_CF_1=(PRO_01*P_CF01 + PRO_11*P_CF11)/PROB__1; P_VR_0=(PRO_00*(P_VR00+(P_CF_0-P_CF00)*(P_CF_0-P_CF00)')+ PRO_10*(P_VR10+(P_CF_0-P_CF10)*(P_CF_0-P_CF10)'))/PROB__0; P_VR_1=(PRO_01*(P_VR01+(P_CF_1-P_CF01)*(P_CF_1-P_CF01)')+ PRO_11*(P_VR11+(P_CF_1-P_CF11)*(P_CF_1-P_CF11)'))/PROB__1; IF J_ITER < START; goto skip; endif; LIKV = LIKV+LIK; skip: J_ITER = J_ITER+1; ENDO; RETP(LIKV); ENDP; PROC (1) = SAVE_OUT(PRMTR1); local prmtr, ppr,qpr,aaa,prob__0,prob__1,p_cf_0,p_cf_1, p_vr_0, p_vr_1, qq, lik, j_iter, h, pr__0_l,pr__1_l, pst_cf0, pst_cf1, pst_vr0, pst_vr1, f_cast0, f_cast1, ss00, ss01, ss10, ss11, pr_vl00, pr_vl01, pr_vl10, pr_vl11, pr_val,k_gn00, k_gn01, k_gn10, k_gn11, p_cf00, p_cf01, p_cf10, p_cf11, p_vr00, p_vr01, p_vr10, p_vr11, pro_00, pro_01, pro_10, pro_11, likv,cf,vr, DTA_MAT, SS,S1,S2,f_cast,coeffs; EXTERNAL PROC TRANS, V_PROB; DTA_MAT=ZEROS(T,4); PRMTR=PRMTR1; PPR=PRMTR[6,1]; @Pr[St=1/St-1=1]@ QPR=PRMTR[7,1]; @Pr[St=0/St-1=0]@ AAA=EYE(TVP_NO); @>>>>>>>>>>>>>>>>>>>>>>>>> INITIAL PROB. Pr[S0/Y0] @ PROB__1=(1-QPR)/(2-PPR-QPR); @Pr[St-1=0/Yt-1], STEADY STATE PROB.@ PROB__0=1-PROB__1 ; @Pr[ST-1=1/Yt-1], STEADY STATE PROB @ P_CF_0=PRIOR_CF; P_CF_1=PRIOR_CF; P_VR_0=PRIOR_VR; P_VR_1=PRIOR_VR; @ Initial values@ LIKV=0.0; J_ITER = 1; DO UNTIL J_ITER>T; @H=AA[J_ITER,.]~zz[J_ITER,.];@ H=AA[J_ITER,.]; PR__0_L= QPR*PROB__0 + (1-PPR)*PROB__1; @Pr[St=0/Yt-1]@ PR__1_L= (1-QPR)*PROB__0 + PPR*PROB__1; @Pr[St=1/Yt-1]@ @===================PREDICTION=============================@ PST_CF0 = AAA * P_CF_0; @ Prediction of coef when St-1=0 @ PST_CF1=AAA*P_CF_1; @Prediction of coef when St-1=1 @ cf=prob__0*pst_cf0 +prob__1*pst_cf1; @Estimate of coefs based on past info@ cf; @ cf=pr__0_l*pst_cf0 +pr__1_l*pst_cf1;@ PST_VR0 = AAA * P_VR_0 * AAA' +QQ; PST_VR1 = AAA * P_VR_1 * AAA' +QQ; vr= prob__0*(pst_vr0+ (cf-pst_cf0)*(cf-pst_cf0)') + prob__1*(pst_vr1+ (cf-pst_cf1)*(cf-pst_cf1)'); @ vr= pr__0_l*(pst_vr0+ (cf-pst_cf0)*(cf-pst_cf0)') + pr__1_l*(pst_vr1+ (cf-pst_cf1)*(cf-pst_cf1)');@ @vr=> variance covariance matrix of coef(cf) calculated based on past information only. @ F_CAST0 = yy[J_ITER,.] - H * PST_CF0; @ Forecast error when St-1=0 @ F_CAST1 = yy[J_ITER,.] - H * PST_CF1; @ Forecast error when St-1=0 @ f_cast= pr__0_l*F_CAST0+ pr__1_l*F_CAST1; SS = H * vr * H' + (prmtr[1,1]^2 +(prmtr[5,1]^2-prmtr[1,1]^2)* pr__1_L); DTA_MAT[j_iter,.]=f_cast~ ss ~ (h*vr*h')~ (prmtr[1,1]^2+(prmtr[5,1]^2-prmtr[1,1]^2)*pr__1_l); SS00= H * PST_VR0 * H' + PRMTR[1,1]^2; @ St-1=0,St=0 @ SS01= H* PST_VR0 * H' + PRMTR[5,1]^2; @ St-1=0,St=1 @ SS10= H * PST_VR1 * H' + PRMTR[1,1]^2; @ St-1=1,St=0 @ SS11= H * PST_VR1 * H' + PRMTR[5,1]^2; @ St-1=1,St=1 @ @<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<@ PR_VL00=V_PROB(F_CAST0,SS00)* QPR*PROB__0; @Pr[St,Yt/Yt-1]@ PR_VL01=V_PROB(F_CAST0,SS01)*(1-QPR)*PROB__0; PR_VL10=V_PROB(F_CAST1,SS10)*(1-PPR)*PROB__1; PR_VL11=V_PROB(F_CAST1,SS11)* PPR*PROB__1; PR_VAL=PR_VL00+PR_VL01+PR_VL10+PR_VL11; @Pr[Yt/Yt-1]@ @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ @=======================UPDATING===============================@ K_GN00 =PST_VR0 * H' / SS00; K_GN01 =PST_VR0 * H' / SS01; K_GN10 =PST_VR1 * H' / SS10; K_GN11 =PST_VR1 * H' / SS11; P_CF00 = PST_CF0 + K_GN00 * F_CAST0; P_CF01 = PST_CF0 + K_GN01 * F_CAST0; P_CF10 = PST_CF1 + K_GN10 * F_CAST1; P_CF11 = PST_CF1 + K_GN11 * F_CAST1; P_VR00 = (EYE(TVP_NO) - K_GN00 * H ) * PST_VR0; P_VR01 = (EYE(TVP_NO) - K_GN01 * H ) * PST_VR0; P_VR10 = (EYE(TVP_NO) - K_GN10 * H ) * PST_VR1; P_VR11 = (EYE(TVP_NO) - K_GN11 * H ) * PST_VR1; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ PRO_00=PR_VL00/PR_VAL; @Pr[St,St-1/Yt]@ PRO_01=PR_VL01/PR_VAL; PRO_10=PR_VL10/PR_VAL; PRO_11=PR_VL11/PR_VAL; PROB__0=PRO_00+PRO_10; @Pr[St=0/Yt]@ PROB__1=PRO_01+PRO_11; @Pr[St=1/Yt]@ @<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<@ P_CF_0=(PRO_00*P_CF00 + PRO_10*P_CF10)/PROB__0; P_CF_1=(PRO_01*P_CF01 + PRO_11*P_CF11)/PROB__1; P_VR_0=(PRO_00*(P_VR00+(P_CF_0-P_CF00)*(P_CF_0-P_CF00)')+ PRO_10*(P_VR10+(P_CF_0-P_CF10)*(P_CF_0-P_CF10)'))/PROB__0; P_VR_1=(PRO_01*(P_VR01+(P_CF_1-P_CF01)*(P_CF_1-P_CF01)')+ PRO_11*(P_VR11+(P_CF_1-P_CF11)*(P_CF_1-P_CF11)'))/PROB__1; skip: J_ITER = J_ITER+1; ENDO; RETP(DTA_MAT~qtr_idx); ENDP; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ PROC TRANS(coef0); @ constraining values of reg. coeff.@ local coef1,m,u; u=0.9999999999; m=0.00000000001; coef1=coef0^2+0.00001; coef1[6:7,.]=m+(u-m)./ (1+exp(-1*coef0[6:7,.])); retp(coef1); endp; PROC TRAN_INV(coef0); local coef1,m,u; u=0.9999999999; m=0.00000000001; coef1=sqrt(coef0-0.00001); coef1[6:7,.]=-1*ln((u-coef0[6:7,.]) ./ (coef0[6:7,.]-m)); retp(coef1); endp; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ PROC V_PROB(EV, HE); @ CALCULATES Pr[Yt/St,Yt-1]@ LOCAL VAL; VAL=(1/SQRT(2*PI*HE))*EXP(-0.5*EV*EV/HE); RETP(VAL); ENDP; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ *gdpsim.prg* new; output file = c:\letter\gdpsim.o18 reset; screen off; @ define global variables @ alph=zeros(3,1); alph= .00985|-1|.1; bzero=zeros(2,1); bzero= .00414|-.45; tbzero= .004|-.454; @ for monte carlo draw tbzero randomly each iteration @ sigv= (1E-05)~0~0~(.05); @sigtb=.02; @ sige= (1.08)^2; tgdp0= ln(6000); gdp0= ln(6000); dm0= 100*(.00985); @ alph; bzero; tbzero; @ t=400; nc=400; outp1= zeros(t-3,nc); q=1; do while q <= nc; init=rndn(2,1); init=(2E-03|.05).*init; tbzero=tbzero+init; @ generate random disturbances @ dist=rndn(t,3); e=sqrt(sige)*dist[.,1]; v=(sigv[1,1]*dist[.,2])~(sigv[1,4]*dist[.,3]); @v1= sigtb*dist[.,4]; @ @initialize variables @ yf= zeros(t,1); yfo=zeros(t,1); y=zeros(t,1); gdp=zeros(t,1); tgdp=zeros(t,1); dm=zeros(t,1); gdp[1]=gdp0; tgdp[1]=tgdp0; dm[1]=dm0; yf[1]= (100~dm0)*bzero; yf[1]= 0; truebet=zeros(t,2); truebet[1,1]=1*tbzero[1]+v[1,1]; tbarg= zeros(t,1); tbarg[1]= tbzero[2]; truebet[1,2]= tbarg[1]; y[1]= (100~dm0)*(truebet[1,.]')+ e[1]; dm[1]= 100*(1~yf[1]/100~(-gdp0+tgdp0))*alph; tgdp[1]= tgdp0 + ln(1+alph[1]); gdp[1]= gdp0+ dm[1]/100 + y[1]/100; tgdp[1]=gdp[1]; @ yf[1]~y[1]~tgdp[1]~gdp[1]~dm[1]; end; @ X= 100*ones(t,1)~zeros(t,1); X[1,2]=dm0; beta=zeros(t,2); beta[1,.]=bzero'; eta=zeros(t,1); eta[1]= y[1]-yf[1]; R=zeros(t,4); R[1,.]=sigv; H= zeros(t,1); H[1,.]= X[1,.]*((R[1,1]~R[1,2])|(R[1,3]~R[1,4]))*X[1,.]'+sige; K=zeros(t,2); K[1,.]= (((R[1,1]~R[1,2])|(R[1,3]~R[1,4]))*(X[1,.]')/H[1])'; newR=zeros(t,4); temp= (eye(2)-K[1,.]'*X[1,.])*((R[1,1]~R[1,2])|(R[1,3]~R[1,4])); newR[1,.]= temp[1,1]~temp[1,2]~temp[2,1]~temp[2,2]; nbeta=zeros(t,2); nbeta[1,.]= beta[1,.]+ (K[1,.]'eta[1])'; X[1,2]=dm[1]; j=2; do while j <= 400; @ j~truebet[j-1,2]~nbeta[j-1,2]~dm[j-1]~tgdp[j-1,.]~gdp[j-1,.]; @ tbarg[j]=.95*tbarg[j-1] + v[j,2]; truebet[j,2]= tbarg[j]; truebet[j,1]=1*truebet[j-1,1]+v[j,1]; beta[j,1]= 1*nbeta[j-1,1]; beta[j,2]= .95*nbeta[j-1,2]; yf[j]= X[j-1,.]*(beta[j,.]'); yf[2]= 0; if j >= 3; yf[j]= (1/(j-2))*(y[j-1]-y[1]); endif; if j >= 18; yf[j]= (1/16)*(y[j-1]-y[j-17]); endif; y[j]= X[j-1,.]*(truebet[j,.]')+ e[j]; dm[j]= 100*(1~yf[j]/100~(-gdp[j-1]+tgdp[j-1]))*alph; tgdp[j]= tgdp[j-1] + ln(1+alph[1]); gdp[j]= gdp[j-1]+ dm[j]/100 + y[j]/100; X[j,2]=dm[j]; eta[j]= y[j]-yf[j]; @ if j <= 12; dm[j]=dm0; beta[j,2]=bzero[2]; eta[j]=0; endif; @ temp1= ((1~0)|(0~.95^2))*((newR[j-1,1]~newR[j-1,2])|(newR[j-1,3]~newR[j-1,4])); R[j,.]= (temp1[1,1]+1E-05)~temp1[1,2]~temp1[2,1]~(temp1[2,2]+.05); H[j,.]= X[j-1,.]*((R[j,1]~R[j,2])|(R[j,3]~R[j,4]))*X[j-1,.]'+sige; K[j,.]= (((R[j,1]~R[j,2])|(R[j,3]~R[j,4]))*X[j-1,.]'./H[j,.])'; temp= (eye(2)-K[j,.]'*X[j-1,.])*((R[j,1]~R[j,2])|(R[j,3]~R[j,4])); newR[j,.]= temp[1,1]~temp[1,2]~temp[2,1]~temp[2,2]; nbeta[j,.]= beta[j,.]+ (K[j,.]'eta[j])'; j=j+1; endo; outp1[1:t-4,q:q+1]= gdp[1:t-4]~(dm[4:t-1]+dm[3:t-2]+dm[2:t-3]+dm[1:t-4]); q=q+2; endo; outp1; end; *GDPSTAT.PRG* new; output file = c:\letter\gdpstat.o18 reset; t=396; tgdp=zeros(t,1); tgdp0=ln(6000); alph= .00985; tgdp[1]= tgdp0+alph; i=2; do while i <= t; tgdp[i]=tgdp[i-1]+alph; i=i+1; endo; load xx[t,400]=c:\letter\gdpsim.o18; avegap= zeros(200,1); msegap=zeros(200,1); vargap=zeros(200,1); corgap=zeros(200,1); maxval=zeros(200,1); minval=zeros(200,1); i=1; do while i <= t; r=.5*(i+1); avegap[r]= (1/t)*sumc(xx[1:t,i]-tgdp[1:t]); msegap[.5*(i+1)]= (1/t)*sumc((xx[1:t,i]-tgdp[1:t]).^2); vargap[.5*(i+1)]= (1/t)*sumc((xx[1:t,i]-tgdp[1:t]-avegap[.5*(i+1)]).^2); v1= xx[1:t,i]-tgdp[1:t]; k=1; qn=20; n=rows(v1); boxq=zeros(1,1); do while k <= qn; j=1+k; m=n-k; boxq[1]=boxq[1]+((v1[1:m]'v1[j:n]/(n-k))/(sumc(v1^2)/n))^2/(n-k); k=k+1; endo; boxq=n*(n+2)*boxq; corgap[.5*(i+1)]= boxq; maxval[.5*(i+1)]=maxc(xx[1:t,i]-tgdp[1:t]); minval[.5*(i+1)]=-maxc(tgdp[1:t]-xx[1:t,i]); i=i+2; endo; ? "Nominal GDP statistics"; ? "average gap=" (1/200)*sumc(avegap); ? "mean sq error=" (1/200)*sumc(msegap); ? "variance of target error=" (1/200)*sumc(vargap); ? "serial correlation in gap=" (1/200)*sumc(corgap); ? "extreme values=" (1/200)*sumc(maxval)~(1/200)*sumc(minval); i=2; do while i <= t; r=.5*(i); avegap[r]= (1/t)*sumc(xx[1:t,i]); msegap[.5*(i)]= (1/t)*sumc((xx[1:t,i]).^2); vargap[.5*(i)]= (1/t)*sumc((xx[1:t,i]-avegap[.5*(i+1)]).^2); v1= xx[1:t,i]-alph; k=1; qn=20; n=rows(v1); boxq=zeros(1,1); do while k <= qn; j=1+k; m=n-k; boxq[1]=boxq[1]+((v1[1:m]'v1[j:n]/(n-k))/(sumc(v1^2)/n))^2/(n-k); k=k+1; endo; boxq=n*(n+2)*boxq; corgap[.5*(i)]= boxq; maxval[.5*(i)]=maxc(xx[1:t,i]); minval[.5*(i)]=-maxc(-xx[1:t,i]); i=i+2; endo; ? "MONEY GROWTH STATISTICS"; ? "average gap=" (1/200)*sumc(avegap); ? "mean sq error=" (1/200)*sumc(msegap); ? "variance of target error=" (1/200)*sumc(vargap); ? "serial correlation in gap=" (1/200)*sumc(corgap); ? "extreme values=" (1/200)*sumc(maxval)~(1/200)*sumc(minval); end;