| Library: | gam |
| See also: | gamfit gintest gintestpl intest interact intestpl pcad |
| Quantlet: | TimeVarAddModel2 | |
| Description: | estimates a dynamic factor model from the form: yt = m0(u) + bt1*m1(u) + bt2*m2(u)... btL*mL(u), where m0 to mL are 2-dimensional invariant basis functions on the grid u and bt0=1. bt1 to btL are scalar weights depending on time T. After estimation, the functions m are orthogonalized under the empirical norm of the observed data points and ordered according to the maximal variance of b. The model only uses fits in the local neighbourhood of the design points. Thus, it is suited for data with degenerated design, like implied volatility data. |
| Usage: | {beta,mhat} = TimeVarAddModel2(x,L,h,u,rmax,L2tol) | |
| Input: | ||
| x | nx4 matrix, columns 1 and 2 contain the grid points, column 3 the observed surface value of the function, column 4 contains a numeric vector with different T dates. T must be at least three times L (T >= 3L) | |
| L | integer, number of basis functions to be estimated (1 < L <= 10) | |
| h | 2x1 vector, contains the bandwidths. First entry refers to column 1 in x, second entry to column 2 in x. | |
| u | mx2 matrix, grid points on which the functions are to be estimated, (the stepwidth in each direction must be constant) | |
| rmax | optional scalar, maximum number of iterations in the estimation, default=301. | |
| L2tol | optional scalar, minimal L2 tolerance of the estimation criterion, default=10^-5. | |
| Output: | ||
| beta | Tx(L+1) matrix, column 1 contains date vector, column 2 weights bt1 corresponding to function m1, column 3 weights bt2 corresponding to m2, etc. | |
| mhat | Mx(L+2) matrix, columns 1 and 2 contain the grid u, column 3 contains m0, column 4 m1, column 5 m2, etc. | |
library("gam")
; generates the equally spaced grid on [0,1]x[0,1] with stepwidth 0.1
usize = 11
u = grid(0|0,0.1|0.1,usize|usize)
Ru = matrix(rows(u))
; simulates the basis functions for L=2
m0 = 5*Ru
m1 = 10*(sin(u[,1]*pi) + cos(u[,2])*pi)
m2 = 3*u[,1] + 2*u[,2]
; generates time series of coefficients
T = 200
; AR(1) with starting at 3 and weight 0.8
b1 = genar(normal(T), 3, 0.8 )
; AR(1) with starting at 1 and weight 0.2
b2 = genar(normal(T), 1, 0.2 )
; blows up the variables to generate the artificial data set
b1b = kron(b1, Ru)
b2b = kron(b2, Ru)
m0b = kron(m0, matrix(T))
m1b = kron(matrix(T), m1)
m2b = kron(matrix(T), m2)
; generates a date vector from 1 to 200 and blows it up
datum = kron(1:T, Ru)
y = kron(matrix(T), u)~(m0b + b1b.*m1b + b2b.*m2b)~datum
; adds noise on the function values
y[,3] = y[,3] + 0.05*sqrt(var(y[,3]))*normal(rows(y))
; estimates the model
{beta, mhat} = TimeVarAddModel2(y, 3, 0.3|0.3, u, 30, 0.001)
setsize(800,300)
d = createdisplay(1,3)
gs1 = grsurface(mhat[,1:2]~mhat[,3])
gs2 = grsurface(mhat[,1:2]~mhat[,4])
gs3 = grsurface(mhat[,1:2]~mhat[,5])
show(d,1,1, gs1)
show(d,1,2, gs2)
show(d,1,3, gs3)
rot =(#(0.70711, -0.13795, -0.69351)~#(-0.70710, -0.13793, -0.69353)~#(0.00001, 0.98079, -0.19508))'
setgopt(d,1,1,"title", "maht0", "rotcos", rot)
setgopt(d,1,2,"title", "mhat1", "rotcos", rot)
setgopt(d,1,3,"title", "mhat2", "rotcos", rot)
; inspects day 3 for example
f = mhat[,1:2]~sum(beta[3,2:cols(beta)].*mhat[,3:cols(mhat)],2)
f = grsurface(f)
g = paf(y[,1:3], y[,4] == 3)
g = setmask(g,"fillcircle","red","size","small")
setsize(400,400)
plot(f,g)
setgopt(plotdisplay,1,1,"title","Estimated functions vs. observations", "rotcos", rot)
In this example an artificial model spanning 200 days is created and estimated. The estimated basis functions are plotted. An additional plot shows the estimated function together with the original observations at time 3.