/* @>>>>>>>>>>>>>>>>>>>>Time varying parameter model >>>>>>>>>>>@ @>>>>>>>>>>>>>>>>>>>>with Markov Varance >>>>>>>>>>>>>>>>>>>@ @>>>>>>>>>>>>>>>>>>>>Using Optmum Proc >>>>>>>>>>>>>>>>>>>>>>@ @ May 9, 1990 version @ cjkim@korea.ac.kr */ new; library optmum,pgraph; load data[120,6]=tvpmrkf.prn; /* column 1: m1===growth rate of quarterly average M1 2: dint=change in the lagged interest rate (3-month T-bill) 3: inf==lagged inflation 4: surp==full employment budget surplus (detrended) 5: m1lag==lag of m1 6: quarter index 1959.3--1987.4, ........................ */ qtr_idx=data[.,6]; T=rows(data); yy=data[.,1]; xt=ones(t,1)~data[.,2:5]; output file=tvpmrkf.out reset; output off; @================= Initialize Global Variables============@ START=11; @1962.1--1989.2 @ PRMTR_NO=9; @ >>Number of parameters to be estimated. >>@ TVP_NO=5; @ >>Number of Time-varying Coefficients >>@ pr_1=sqrt(0.2); pr_2=sqrt(0.2); pr_3=sqrt(0.2); pr_4=sqrt(0.2);pr_5=sqrt(0.2);pr_6=sqrt(0.2); pr_7=sqrt(1.050); pr_8=0.9;pr_9=0.9; @initial values of parameters@ PRMTR_IN=pr_1|pr_2|pr_3|pr_4|pr_5|pr_6|pr_7|pr_8|pr_9; prmtr_in=0.36072|0.07080|0.079|0.0791|0.070|0.072|0.89439| 0.93266|0.95990; output on; "initial values of prmtr is(vre,intcpt,dint,inf,surp,mlag,w,p,q"; prmtr_in'; output off; PRIOR_CF=zeros(TVP_NO,1); @Initial values of reg coefs@ prior_vr=eye(tvp_no)*100; @ Maximum Likelihood Estimation @ @==================================================@ {xout,fout,gout,hout}=optmum(&lik_fcn,tran_inv(PRMTR_IN)); @ prmtr estimates, -log lik value, Grandient, Hessian@ "Calculating Hessian..... Please be patient!!!!"; hout0=hessp(&lik_fcn,xout); hout=inv(hessp(&lik_fcn,xout)); PRM_FNL=TRANS(xout); @ Estimated coefficients, constrained@ 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; "==============================================================="; OUTPUT OFF; F_SS=SAVE_OUT(xout); OUTPUT FILE=TVPMRKF.DTA RESET; F_SS; OUTPUT OFF; xy(seqa(1,1,rows(f_ss)),f_ss[.,1]); xy(seqa(1,1,rows(f_ss)),f_ss[.,2:4]); 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[8,1]; @Pr[St=1/St-1=1]@ QPR=PRMTR[9,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~0)| (0~PRMTR[3,1]^2~0~0~0)| (0~0~PRMTR[4,1]^2~0~0)| (0~0~0~PRMTR[5,1]^2~0)| (0~0~0~0~PRMTR[6,1]^2); LIKV=0.0; J_ITER = 1; DO UNTIL J_ITER>T; H=xt[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[7,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[7,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 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; EXTERNAL PROC TRANS, V_PROB; PRMTR=TRANS(PRMTR1); DTA_MAT=ZEROS(T,4); PPR=PRMTR[8,1]; @Pr[St=1/St-1=1]@ QPR=PRMTR[9,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@ QQ=(PRMTR[2,1]^2~0~0~0~0)|(0~PRMTR[3,1]^2~0~0~0)|(0~0~PRMTR[4,1]^2~0~0)| (0~0~0~PRMTR[5,1]^2~0)|(0~0~0~0~PRMTR[6,1]^2); LIKV=0.0; J_ITER = 1; DO UNTIL J_ITER>T; H=xt[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=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;@ f_cast= prob__0*F_CAST0+ prob__1*F_CAST1; SS = H * vr * H' + (prmtr[1,1]^2 +(prmtr[7,1]^2-prmtr[1,1]^2)* pr__1_L); DTA_MAT[j_iter,.]=f_cast~ ss ~ (h*vr*h')~ (prmtr[1,1]^2+(prmtr[7,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[7,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[7,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[start:t,.]~qtr_idx[start:t,.]); ENDP; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@ PROC TRANS(coef0); @ constraining values of reg. coeff.@ local coef1,m,u; coef1=coef0; coef1[8:9,.]=1./ (1+exp(-1*coef0[8:9,.])); retp(coef1); endp; PROC TRAN_INV(coef0); local coef1,m,u; coef1=coef0; coef1[8:9,.]=-1*ln((1-coef0[8:9,.]) ./ (coef0[8:9,.])); 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; @>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>@