Time Varying Parameter model

[ivory4]

Single equation Time varying parameter model, estimated using kalman filter.
y_t=x_tb_t+u_t,
b_t=b_t-1+v_t,
Var(u_t)=R,
Var(v_t)=Q,
User Guide provides an example for VAR, is there one for a single equation?

[TomDoan]

From the Tsay textbook examples - tsayp510.prg and tsayp533.prg. From the West and Harrison textbooks examples, westp081.prg.

[ivory4]

I have some questions about the code.
1. What is the function for line: nonlin leps leta?
I have the RM at hand and it says that nonlin set a list of free parameter before using DLM, and what does this mean?
2. sv=1.0 in this code, are they pegging the sv or using 1.0 as initial value?
When I change the sv, the estimation is changed, but if it is just initial value, the estimation should not change much.


Thanks for reading and replying.

Code: Select all

*
* Example 11.2, pp 533-535
*
open data m-fac9003.txt
calendar(m) 1990
data(format=free,org=columns) 1990:1 2003:12 aa age cat f fdx gm hpq kmb mel nyt pg trb txn sp
set gm = gm*.01
set sp = sp*.01
*
* Fixed coefficients model done with linreg
*
linreg gm
# constant sp
*
* Done with DLM. sv=1.0 with VAR=CONCENTRATE concentrates out the variance of the
* measurement equation. FREE=2 has no effect on the estimates, but corrects the
* likelihood function so it matches LINREG's by computing the likelihood
* conditional on the (freely estimated) coefficients.
*
dlm(y=gm,c=||1.0,sp||,sv=1.0,var=concentrate,exact,free=2,type=smoothed) / xstates vstates
*
* With the smoothed estimates, the state means and variances are the same for all
* time periods. Note that the %VARIANCE (the concentrated estimate of the
* measurement equation variance) isn't corrected for degrees of freedom, while the
* calculation used in the text is.
*
disp xstates(2003:12)
disp tr(%sqrt(%xdiag(vstates(2003:12))*%variance))
*
* Now estimate the component variances. Note that we're pegging sv=1.0 and using
* var=concentrate to  concentrate out the unknown measurement variance. The other
* variance estimates will have to be multiplied by %VARIANCE to get their true
* estimated values. In component models like this, it can be very difficult to
* estimate all the variances directly because they are so highly correlated. In
* this case, the estimates given in the book aren't really the correct likelihood
* maximizers. The likelihood is a bit higher with a somewhat higher variance on
* beta.
*
nonlin leps leta
compute leps=leta=0.0
dlm(y=gm,c=||1.0,sp||,sv=1.0,var=concentrate,sw=%diag(||exp(leta),exp(leps)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates
disp "Component Std Dev"
disp "Measurement" @20 sqrt(%variance)
disp "Alpha" @20 sqrt(%variance*exp(leta))
disp "Beta" @20 sqrt(%variance*exp(leps))
*
set alpha = xstates(t)(1)
set beta  = xstates(t)(2)
set expret = %dot(||1.0,sp||,xstates)
*
spgraph(hfields=2,vfields=2,footer="Figure 11.5 Time plots from a time-varying CAPM")
graph(hlabel="Excess return")
# gm
graph(hlabel="Expected return")
# expret
graph(hlabel="Alpha")
# alpha
graph(hlabel="Beta")
# beta
spgraph(done)

[TomDoan]

The first DLM has the coefficients "fixed" - they follow a random walk with zero increment variance (there's no SW option). The only reason the estimates change over the sample is that you get more information. When you Kalman smooth (as opposed to filter), the coefficients will be the same across the sample. In the second case, they follow a RW with a non-zero variance; the "target" can change over the course of the sample. As you can see from the formula for SW, leta is the log of the increment variance for the constant, and leps is the log of the increment variance for the beta coefficient. The second DLM estimates those by maximum likelihood.

Code: Select all

nonlin leps leta
compute leps=leta=0.0
dlm(y=gm,c=||1.0,sp||,sv=1.0,var=concentrate,sw=%diag(||exp(leta),exp(leps)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates

There are actually three variances in this model - the two in the SW matrix and variance in the measurement error. sv=1.0,var=concentrate does the estimation by concentrating out a scale factor in the variances. This is described in the User's Guide. If you change this to (say) sv=2.0,var=concentrate, you'll get (exactly) the same results for xstates; the SW variances will go up by a factor of 2, which means that leps and leta will increase by log(2). There is, however, no reason to change sv=1.0 if you use var=concentrate. The alternative way of handling this is to put the variance of the measurement error itself into the parameter set:

Code: Select all

nonlin lsigma leps leta
compute leps=leta=lsigma=0.0
dlm(y=gm,c=||1.0,sp||,sv=exp(lsigma),sw=%diag(||exp(leta),exp(leps)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates

[ivory4]

In my model, xstates should be 375x5, then I set as following
set intercept = xstates(t)(1)
set beta = xstates(t)(2)

But it tells me "SX22. Expected Type SERIES, Got REAL Instead"
What is the correct way of reading the estimated time varying coefficients out?

The part of code after inputting data is as following.

And also if I would like to graph those estimates, is that the correct way?

Code: Select all

nonlin lsigma leps leta
compute leps1=leps2=leps3=leps4=leps5=lsigma=0.0
dlm(y=ue,c=||1.0,def,ygap,ue1,ue2||,sv=exp(lsigma),sw=%diag(||exp(leps1),exp(leps2),exp(leps3),exp(leps4),exp(leps5)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates
disp "Component Std Dev"
disp "Measurement" @20 sqrt(%variance)

disp "Beta" @20 sqrt(%variance*exp(leps2))

set intercept = xstates(t)(1)
set beta  = xstates(t)(2)
set gama  = xstates(t)(3)
set row1  = xstates(t)(4)


spgraph(hfields=2,vfields=2,footer="Figure ")
graph(hlabel="intercept")
# intercept
graph(hlabel="beta")
# beta
graph(hlabel="gamma")
# gama
graph(hlabel="row1")
# row1
spgraph(done)

[TomDoan]

In my model, xstates should be 375x5, then I set as following
set intercept = xstates(t)(1)
set beta = xstates(t)(2)

But it tells me "SX22. Expected Type SERIES, Got REAL Instead"
What is the correct way of reading the estimated time varying coefficients out?

Based upon the error message, it looks as if BETA is used as a REAL earlier in the program.

And also if I would like to graph those estimates, is that the correct way?

That looks fine.

[ivory4]

Code: Select all

cal 1972 8 12
allocate 0 2003:12
open data "C:*******.xls"
data(format=xls,org=col) / ue def ygap ue1 ue2


nonlin leps1 leps2 leps3 leps4 leps5 lsigma
compute leps1=leps2=leps3=leps4=leps5=lsigma=0.0
dlm(y=ue,c=||1.0,def,ygap,ue1,ue2||,sv=exp(lsigma),sw=%diag(||exp(leps1),exp(leps2),exp(leps3),exp(leps4),exp(leps5)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates
disp "Component Std Dev"
disp "Measurement" @20 sqrt(%variance)
disp "Alpha" @20 sqrt(%variance*exp(leps1))
disp "Beta" @20 sqrt(%variance*exp(leps2))
disp "Gap"  @20 sqrt(%variance*exp(leps3))
disp "Row1" @20 sqrt(%variance*exp(leps4))
disp "Row2" @20 sqrt(%variance*exp(leps5))

set intercept = xstates(t)(1)
set beta  = xstates(t)(2)
set gama  = xstates(t)(3)
set row1  = xstates(t)(4)

spgraph(hfields=2,vfields=2,footer="Figure ")
graph(hlabel="intercept")
# intercept
graph(hlabel="beta")
# beta
graph(hlabel="gamma")
# gama
graph(hlabel="row1")
# row1
spgraph(done)

Error message
## SX22. Expected Type SERIES, Got VECTOR Instead
>>>>set intercept <<<<

This means that I am not recording the results in a correct way.

[TomDoan]

Error message
## SX22. Expected Type SERIES, Got VECTOR Instead
>>>>set intercept <<<<

This means that I am not recording the results in a correct way.

If you tried writing this with INTERCEPT as a VECTOR and then tried to re-run it with it as a SERIES without clearing the program first, that's what will happen. Before re-running a program after you make a change like that, hit the "Clear Program" toolbar button first (far right button). That clears out the symbol table so you don't have conflicts.

Code: Select all

dlm(y=ue,c=||1.0,def,ygap,ue1,ue2||,sv=exp(lsigma),sw=%diag(||exp(leps1),exp(leps2),exp(leps3),exp(leps4),exp(leps5)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates
disp "Component Std Dev"
disp "Measurement" @20 sqrt(%variance)
disp "Alpha" @20 sqrt(%variance*exp(leps1))

As you have this written, the measurement error variance will be exp(lsigma). It's %variance only if you peg SV at 1.0 and use VAR=CONCENTRATE. If you have modelled all the variances (which you've done here), take the %variance* out of all of those other calculations like the value for ALPHA. You just want sqrt(exp(leps1)).

[ivory4]

Code: Select all

nonlin leps1 leps2 leps3 leps4 leps5 lsigma
compute leps1=leps2=leps3=leps4=leps5=lsigma=log(0.01)
dlm(y=ue,c=||1.0,def,ygap,ue1,ue2||,sv=exp(lsigma),sw=%diag(||exp(leps1),exp(leps2),exp(leps3),exp(leps4),exp(leps5)||),$
   exact,method=bfgs,type=smoothed) / xstates vstates

If I set type=smoothed, the states are too close to each other, like
1998:06 0.095926629962
1998:07 0.095926629961
What would be a solution?

I tried to use estimated hyperparameter, I need to manually add those number to COMPUTE and run dlm again?

[TomDoan]

That would happen if the variance on one of the "w" terms goes (effectively) to zero, which would probably mean an estimated log variance of maybe -15 or smaller. With that many coefficients, I'm not surprised that at least one of the coefficients would behave that way.

Note, however, that I don't think you get good results freely estimating the variances in a model like this. There's a tendency for the maximum likelihood estimates to have a few of the regression coefficients moving dramatically to fit the data almost exactly. Something more tightly parameterized (for instance, having the SW matrix be a single scale parameter x a fixed matrix) tends to work better.

[ivory4]

That would happen if the variance on one of the "w" terms goes (effectively) to zero, which would probably mean an estimated log variance of maybe -15 or smaller. With that many coefficients, I'm not surprised that at least one of the coefficients would behave that way.

That is correct, the variance is very very small.

Note, however, that I don't think you get good results freely estimating the variances in a model like this. There's a tendency for the maximum likelihood estimates to have a few of the regression coefficients moving dramatically to fit the data almost exactly.

And I think the intercept is somehow large to be intuitively consistent with my model and I think it take away the variance in the coefficients.

Something more tightly parameterized (for instance, having the SW matrix be a single scale parameter x a fixed matrix) tends to work better.

In the workbook, 3.2 method 4?
Measurement variance concentrated out; ratio in log form?

[TomDoan]

Time-varying coefficients models are particularly tricky because there isn't likely to be as much prior knowledge to set the variances as there will be for the unobservable components models. Even in the case of a time-varying beta, where you at least have a good idea that the value is likely to be somewhere in the range of 0 to 2, the combination of the variances might produce results which look "wrong". This is from Tsay's textbook (tsayp505.prg, the data file is attached):

Code: Select all

open data m-fac9003.txt
calendar(m) 1990
data(format=free,org=columns) 1990:1 2003:12 aa age cat f fdx $
  gm hpq kmb mel nyt pg trb txn sp
*
compute sv=0.1**2,seta=0.02**2,seps=0.04**2
*
* The above values are quite far from the ones which maximize the
* likelihood. Note that the "C" matrix is actually a time-dependent
* formula, since it includes a reference to the series "sp"
*
dlm(y=gm,c=||1.0,sp||,sv=sv,sw=%diag(||seta,seps||),exact,type=smoothed) / xstates vstates
*
linreg gm
# constant sp
*
set tvbeta = xstates(t)(2)
graph(vgrid=%beta(2),footer="Smoothed estimates of time-varying betas")
# tvbeta

The seps (which is the variance on the change in the beta coefficient) is .04^2. The coefficient follows a random walk, so that would give a variance of .04^2 x 168=.27 over the course of the sample. That doesn't seem at all unreasonable. However, if you look at the results, the beta's are flying all over the place; making the estimation worthless. What's wrong? The measurement equation variance of .1^2. The data set has returns multiplied by 100 to make percentages. The .1^2 (and the value for seta as well) are scaled for data in fractions. With the values of the data in the data set being enormous relative to that .1^2, it becomes much "cheaper" to explain the data by moving the beta around than to assign the error to the measurement equation. Hence the wild betas. If you correct the scales on the sv and seta to:

Code: Select all

compute sv=10.0**2,seta=2.0**2,seps=0.04**2

you get much more reasonable results, though you might find that to be "too" stiff.

West and Harrison offer what may be a better way to handle TVC models, although it doesn't seem to have really caught on. Instead of adding variance to the state equation at each time period, it multiplies the last period variance by a constant somewhat larger than 1. (It's specified as a "discount" parameter; the reciprocal of that multiplier). You can try this instead, and see the sensitivity to the choice of discount.

Code: Select all

dlm(y=gm,c=||1.0,sp||,sv=sv,discount=.99,exact,type=smoothed) / xstates vstates

m-fac9003.txt

(16.41 KiB) Downloaded 1253 times

[ivory4]

Time-varying coefficients models are particularly tricky because there isn't likely to be as much prior knowledge to set the variances as there will be for the unobservable components models. Even in the case of a time-varying beta, where you at least have a good idea that the value is likely to be somewhere in the range of 0 to 2, the combination of the variances might produce results which look "wrong". This is from Tsay's textbook (tsayp505.prg, the data file is attached):

Code: Select all

open data m-fac9003.txt
calendar(m) 1990
data(format=free,org=columns) 1990:1 2003:12 aa age cat f fdx $
  gm hpq kmb mel nyt pg trb txn sp
*
compute sv=0.1**2,seta=0.02**2,seps=0.04**2
*
* The above values are quite far from the ones which maximize the
* likelihood. Note that the "C" matrix is actually a time-dependent
* formula, since it includes a reference to the series "sp"
*
dlm(y=gm,c=||1.0,sp||,sv=sv,sw=%diag(||seta,seps||),exact,type=smoothed) / xstates vstates
*
linreg gm
# constant sp
*
set tvbeta = xstates(t)(2)
graph(vgrid=%beta(2),footer="Smoothed estimates of time-varying betas")
# tvbeta

The seps (which is the variance on the change in the beta coefficient) is .04^2. The coefficient follows a random walk, so that would give a variance of .04^2 x 168=.27 over the course of the sample. That doesn't seem at all unreasonable. However, if you look at the results, the beta's are flying all over the place; making the estimation worthless. What's wrong? The measurement equation variance of .1^2. The data set has returns multiplied by 100 to make percentages. The .1^2 (and the value for seta as well) are scaled for data in fractions. With the values of the data in the data set being enormous relative to that .1^2, it becomes much "cheaper" to explain the data by moving the beta around than to assign the error to the measurement equation. Hence the wild betas. If you correct the scales on the sv and seta to:

Code: Select all

compute sv=10.0**2,seta=2.0**2,seps=0.04**2

you get much more reasonable results, though you might find that to be "too" stiff.

West and Harrison offer what may be a better way to handle TVC models, although it doesn't seem to have really caught on. Instead of adding variance to the state equation at each time period, it multiplies the last period variance by a constant somewhat larger than 1. (It's specified as a "discount" parameter; the reciprocal of that multiplier). You can try this instead, and see the sensitivity to the choice of discount.

Code: Select all

dlm(y=gm,c=||1.0,sp||,sv=sv,discount=.99,exact,type=smoothed) / xstates vstates

m-fac9003.txt

These two are very useful and generate "better reasonable" results.

Thank you for information. I may follow up more.

[ivory4]

## SR5. Using Range -32765 to -32765 of Series BETA. (Incorrect start,end or SMPL)

Why this may happen from time to time? When I use "Clear program", sometimes it could be solved sometimes not and I have to restart the application.

The program I am running is only slightly different from what we discussed here.

[TomDoan]

It's hard to tell without seeing the exact program. However, if you run part of a program, you might start at a spot which leaves some variable uninitialized or set to a value from later in the program when it needs the values from the start.