* * Example of bivariate HP filter. This can be adapted to more than two * variables, by changing the n=2 and redefining the YF vector. * open data oecdsample.rat calendar(q) 1981 data(format=rats) 1981:1 2006:4 canexpgdpchs usaexpgdpch * set canrgdps = log(canexpgdpchs) set usargdps = log(usaexpgdpch) * compute n=2 * dec frml[vect] yf frml yf = ||canrgdps,usargdps|| * compute lambda=1600.0 * * Standard local trend model for the common growth component * compute [rect] a=||1.0,1.0|0.0,1.0|| compute [symm] sw=%diag(||0.0,1.0/lambda||) * compute [symm] sv=%identity(n) * * Each series is modeled as y(i)=a(i)+gamma(i)*growth+v(i) * * Because there are n+1 "intercepts" in the complete model (1 for each * data series + 1 in the growth component), the series intercepts are * constrained to sum to zero. (This is accomplished by using the %perp of * a vector of 1's premultiplying an n-1 vector of free coefficients). * The "sv" matrix in the DLM gives the weighting scheme applied to the * errors in the observables. If this is the identity matrix, the series * get equal weights. * dec vect gamma(n) dec vect theta(n-1) * nonlin gamma theta compute gamma=%const(1.0),theta=%const(0.0) * compute oneperp=%perp(%fill(1,n,1.0)) * dec frml[rect] cf frml cf = tr(gamma)~~%zeros(1,n) * * This has some convergence problems. The likelihood is extremely flat * in the parameters, so almost identical trends can be generated many * ways. * dlm(type=smooth,y=yf,mu=oneperp*theta,a=a,sw=sw,sv=sv,c=cf,\$ presample=diffuse,var=concentrated,pmethod=simplex,piters=20,\$ method=bfgs,iters=500) / xstates * set canf = %dot(%xrow(oneperp,1),theta)+%dot(%xcol(cf(t),1),xstates(t)) set usaf = %dot(%xrow(oneperp,2),theta)+%dot(%xcol(cf(t),2),xstates(t)) spgraph(vfields=2,footer="Data with HP Trend") graph(header="Canada") 2 # canrgdps # canf graph(header="US") 2 # usargdps # usaf spgraph(done) * set growth = xstates(t)(1) graph(footer="Common Growth Component") # growth