Time Varying Parameter Example

[TomDoan]

This is an example of a time-varying coefficients regression, taken from Lutkepohl's "New Introduction to Multiple Time Series Analysis", pp 637-640. It estimates the model several different ways, including Gibbs sampling and maximum likelihood.

Code: Select all

*
* Example 18.5 from pp 637-640
*
open data e1.dat
calendar(q) 1960
data(format=prn,org=columns,skips=6) 1960:1 1982:4 invest income cons
*
set dlogy = log(income/income{1})
set dlogc = log(cons/cons{1})
*
linreg(define=ceqn,robust) dlogc
# constant dlogy dlogy{1} dlogc{1} dlogy{2} dlogc{2}
*
* These will be used for graphs
*
dec vect[series] b(6) lower(6) upper(6)
*
* The magnitudes of many of the variances are quite small. This causes
* some numerical problems for variational methods (such as BFGS). To
* counter this, we estimate everything in log form.
*
dec vect lsigsqw(6)
*
* These are Lutkepohl's values. It's not clear what was used as the
* prior mean and variance (gamma0-bar and sigma0 in the text), so these
* are similar but not exactly the same. The estimates themselves are
* quite different from ML estimates that we get below, so these are
* probably heavily dependent upon those unknown priors.
*
compute lsigsqw=%log(||2.04e-5,.14e-2,.46e-2,.45e-2,.51e-2,.62e-2||)
compute lsigsqv=log(3.91e-5)
dlm(c=%eqnxvector(ceqn,t),sw=%diag(%exp(lsigsqw)),sv=exp(lsigsqv),exact,y=dlogc,$
  method=solve,type=smooth) 1960:4 * xstates vstates
do i=1,6
   set b(i) 1960:4 * = xstates(t)(i)
   set lower(i) 1960:4 * = b(i)-2.0*sqrt(vstates(t)(i,i))
   set upper(i) 1960:4 * = b(i)+2.0*sqrt(vstates(t)(i,i))
end do i
*
procedure DoGraph
option string   footer
option equation equation
*
spgraph(vfields=3,hfields=2)
do i=1,%eqnsize(equation)
   graph(header=%eqnreglabels(equation)(i)) 3
   # b(i)
   # lower(i) / 2
   # upper(i) / 2
end do i
spgraph(done)
end DoGraph
*
@DoGraph(footer="At Lutkepohl's Values",equation=ceqn)
*
* These are the ML estimates with exact treatment for initial conditions.
*
nonlin lsigsqv lsigsqw
dlm(c=%eqnxvector(ceqn,t),sw=%diag(%exp(lsigsqw)),sv=exp(lsigsqv),exact,y=dlogc,$
  method=bfgs,iters=500,type=smooth) 1960:4 * xstates vstates
*
do i=1,6
   set b(i) 1960:4 * = xstates(t)(i)
   set lower(i) 1960:4 * = b(i)-2.0*sqrt(vstates(t)(i,i))
   set upper(i) 1960:4 * = b(i)+2.0*sqrt(vstates(t)(i,i))
end do i
@DoGraph(footer="At ML Values",equation=ceqn)
*
* ML estimates of discount model. This tends to have very wide bands at
* the front end of the data, becoming much narrower as more data are
* obtained.
*
compute delta=.95
nonlin lsigsqv delta
dlm(c=%eqnxvector(ceqn,t),discount=delta,sv=exp(lsigsqv),exact,y=dlogc,$
  method=bfgs,type=smooth) 1960:4 * xstates vstates
*
do i=1,6
   set b(i) 1960:4 * = xstates(t)(i)
   set lower(i) 1960:4 * = b(i)-2.0*sqrt(vstates(t)(i,i))
   set upper(i) 1960:4 * = b(i)+2.0*sqrt(vstates(t)(i,i))
end do i
@DoGraph(footer="At ML Estimates of Discount Model",equation=ceqn)
*
* Gibbs sampling
*
do i=1,6
   set b(i)     1960:4 * = 0.0
   set lower(i) 1960:4 * = 0.0
   set upper(i) 1960:4 * = 0.0
end do i
*
* Prior for variances on the coefficient drift. These are inverse gamma
* with a "mean" of sumsq0/df0. Because of the scale difference between
* it and the other coefficients, this is probably a bit too loose on the
* CONSTANT.
*
dec vect sigsqw(6)
compute sigsqw=%exp(lsigsqw)
compute sigsqv=exp(lsigsqv)
compute sumsq0=.0005
compute df0   =5
*
compute nburn=5000,ndraws=5000
infobox(action=define,progress,lower=-nburn,upper=ndraws) $
   "Gibbs Sampling"
do draw=-nburn,ndraws
   *
   * Draw w's and v's given the variances.
   *
   dlm(c=%eqnxvector(ceqn,t),sw=%diag(sigsqw),sv=sigsqv,exact,y=dlogc,$
     type=csimulate,what=what,vhat=vhat) 1960:4 * xstates
   *
   * Draw sigsqv given v's. (This has a non-informative prior).
   *
   sstats 1960:4 * vhat(t)(1)^2>>sumsqv
   compute sigsqv=sumsqv/%ranchisqr(%nobs)
   *
   * Draw sigsqw given w's.
   *
   do i=1,6
      sstats 1960:4 * what(t)(i)^2>>sumsqw
      compute sigsqw(i)=(sumsqw+sumsq0)/%ranchisqr(%nobs+df0)
   end do i
   *
   * Update information if we're out of the burn-in period
   *
   infobox(current=draw)
   if draw>0 {
      do i=1,6
         set b(i) 1960:4 *     = b(i)+xstates(t)(i)
         set lower(i) 1960:4 * = lower(i)+xstates(t)(i)^2
      end do i
   }
end do draw
infobox(action=remove)
*
do i=1,6
   set b(i) 1960:4 * = b(i)/ndraws
   set lower(i) 1960:4 * = b(i)-2.0*sqrt(lower(i)/ndraws-b(i)^2)
   set upper(i) 1960:4 * = 2*b(i)-lower(i)
end do i
@DoGraph(footer="Gibbs Sampling",equation=ceqn)

Data file:

e1.dat

(2.64 KiB) Downloaded 1058 times

[tclark]

Tom --

If I may ask a clarifying question about the Gibbs sampling estimation in the example you provided a few weeks ago.... In the recent literature on VARs with time-varying parameters (e.g., Cogley and Sargent (2005 Review of Economic Dynamics) and Primiceri (2005 RESTUD)), a number of papers have used the Carter and Kohn (1994 Biometrika) approach of filtering backward and drawing coefficients at each t-1 conditioned on period t information. Do I correctly understand that, in your Gibbs sampling code, the CSIMULATE option on DLM uses the Durbin-Koopman 2002 approach of filtering backward and drawing coefficients? Thanks very much.

[TomDoan]

That's right. Carter-Kohn and Durbin-Koopman are different algorithms for implementing the same type of conditional simulation.

[Henrique Andrade]

Dear Tom, my model is a Taylor Rule for the Brazilian economy:

interest(t) = b0 + b1*interest(t-1) + b2*inflation(t) + b3*gap(t) + b4(t)*gap_usa(t) + e(t)

I'm trying to check if the impact of external product gap, b4(t), is increasing overtime. I adapted your example using my data (4 variables with monthly observations) and everything worked just fine. Unfortunately, as I am a beginner in the RATS world (and in the state-space models too), a question had remained. The values of gamma0-bar and sigma0. In the example you use:

Code: Select all

compute lsigsqw=%log(||2.04e-5,.14e-2,.46e-2,.45e-2,.51e-2,.62e-2||)

How can I make the shocks of "interest", "inflation", and "gap" have zero variance in the state equation? I think I need to change this line to something like this:

Code: Select all

compute lsigsqw=%log(||0,0,0,0,X||)

Am I right?

Best regards,
Henrique

[TomDoan]

That is allowing no drift in any of the coefficients except the gap. Usually, I also recommend that you allow drift in the constant (if the drifting coefficient has a non-zero mean) but the gap probably hangs close enough to zero that you don't necessarily need that.

[Henrique Andrade]

That is allowing no drift in any of the coefficients except the gap. Usually, I also recommend that you allow drift in the constant (if the drifting coefficient has a non-zero mean) but the gap probably hangs close enough to zero that you don't necessarily need that.

Dear Tom, thanks a lot for your help. Now I have just one more doubt: In your example, how did you define the variance values?

Code: Select all

compute lsigsqw=%log(||2.04e-5,.14e-2,.46e-2,.45e-2,.51e-2,.62e-2||)

Best,

[TomDoan]

That is allowing no drift in any of the coefficients except the gap. Usually, I also recommend that you allow drift in the constant (if the drifting coefficient has a non-zero mean) but the gap probably hangs close enough to zero that you don't necessarily need that.

Dear Tom, thanks a lot for your help. Now I have just one more doubt: In your example, how did you define the variance values?

Code: Select all

compute lsigsqw=%log(||2.04e-5,.14e-2,.46e-2,.45e-2,.51e-2,.62e-2||)

Best,

I didn't. Luetkepohl did. Those are given in the book. And he didn't have any record of how the model was actually estimated. (Standard maximum likelihood gives quite different variances).