* * TVARYING.RPF * Example of time-varying parameters VAR. * * RATS User's Guide, Example from Section 7.14.2,7.14.3 * cal(q) 1959 open data haversample.rat data(format=rats) 1959:1 2006:4 gdph fm1 log gdph log fm1 compute kstart=1960:1,kend=2006:4 equation kalmeq gdph # constant gdph{1 to 4} fm1{1 to 4} declare vector b(9) bmean(9) ewise bmean(i)=(i==2) associate(perm) kalmeq b * linreg gdph # constant gdph{1 to 4} compute gseesq=.9*%seesq linreg fm1 # constant fm1{1 to 4} compute mseesq=.9*%seesq * * Set up the system. Use GSEESQ as the measurement equation variance * system kalmeq kfset(constant,noscale,likelihood=likely) xxx # gseesq tvarying tvx end(system) * * The hyperparameters (PIx) below are the following: * PI5 = Overall tightness. This corresponds to TIGHTNESS^2 on * SPECIFY. * PI2 = Relative tightness on other variables. This corresponds * to w^2 on the symmetric prior. * PI3 = Relative tightness on the constant. (No correspondence). * PI7 = Relative tightness on time variation. * PI8 = Shrinkage factor toward mean (between 0 and 1 , 1 = no shrinkage) * * This uses a HARMONIC decay with DECAY=.5 * compute pi5 = .20^2 compute pi2 = .5^2 compute pi3 = 10000. compute pi7 = 0.00000001 compute pi8 = 1.00 * dim xxx(9,9) tvx(9,9) * compute likely=%const(0.0) compute b=bmean compute xxx=%const(0.0) compute xxx(1,1)=pi5*pi3*gseesq do i=1,4 compute decayfac=i^(-2*.5) compute xxx(i+1,i+1)=pi5*decayfac compute xxx(i+5,i+5)=pi5*pi2*decayfac*gseesq/mseesq end do i compute tvx=pi7*xxx * do time=kstart,kend if time==kstart kalman(start=time) else { compute b=pi8*b+(1-pi8)*bmean compute xxx=(pi8^2)*xxx kalman(print=(time==kend)) } end do time * * Compute the concentrated log likelihood * compute nobs=kend-kstart+1 compute loglikely=-.5*(likely(1,1)+nobs*log(likely(2,1)/nobs)) compute pi5 = .20^2 , pi2 = .5^2 , pi3 = 10000. compute pi7 = 0.00000001 , pi8 = 1.00 nonlin pi8 pi7 dec real loglikely dim xxx(9,9) tvx(9,9) find(method=genetic,iterations=200) maximum loglikely { if pi7<0.or.pi8>1.00 { compute loglikely=%na next } compute likely=%const(0.0) compute b=bmean compute xxx=%const(0.0) compute xxx(1,1)=pi5*pi3*gseesq do i=1,4 compute decayfac=i^(-2*.5) compute xxx(i+1,i+1)=pi5*decayfac compute xxx(i+5,i+5)=pi5*pi2*decayfac*gseesq/mseesq end do i compute tvx=pi7*xxx * do time=kstart,kend if time==kstart kalman(start=time) else { compute b=pi8*b+(1-pi8)*bmean compute xxx=(pi8^2)*xxx kalman } end do time compute nobs=kend-kstart+1 compute loglikely=-.5*(likely(1,1)+nobs*log(likely(2,1)/nobs)) } end find