If your looking at state space models for forecasting, the Hyndman's
forecasting has quite a bit.
On 8/14/07, Bernardo Ribeiro <[EMAIL PROTECTED]> wrote:
>
> Here is one of the codes.....
>
> Thanks a lot...
>
>
> module(finmetrics)
>
> genvasicek.ssf = function(param, tau=NULL, freq=1/12){
>
> ## 1. Check for valid inputs
>
> if(length(param) < 5)
> stop("Parameters must have length greater than 4")
> N = length(param) - 13
> if (length(tau) != N)
> stop("Length of Tau is inconsistent with Parameters")
>
> ## 2.Extract Parameters and Impose Constraints
>
> k1 = (param[1])
> k2 = (param[2])
> k3 = (param[3])
>
> lambda1 = param[4]
> lambda2 = param[5]
> lambda3 = param[6]
>
> Delta = (param[7])
>
> s1 = (param[8])
> s2 = (param[9])
> s3 = (param[10])
>
> r21 = param[11]
> r31 = param[12]
> r32 = param[13]
>
> st1 = (param[14])
> st2 = (param[15])
> st3 = (param[16])
> st4 = (param[17])
>
> ## 3. Compute Matrixes
>
> A1 = function(x,tau){
> (1-exp(-x*tau))/x
> }
>
> u1 = -A1(k1,tau) ## ok
> u2 = -A1(k2,tau) ## ok
> u3 = -A1(k3,tau) ## ok
> H = cbind (-u1/tau,-u2/tau,-u3/tau) ## ok
>
> v1 = (lambda1/k1)*(tau - A1(k1,tau)) + (lambda2/k2)*(tau - A1(k2,tau))
> +
> (lambda3/k3)*(tau - A1(k3,tau))
> v2 = +3 * Delta * tau # Teste
>
> v311 = ((s1*s1)/(k1*k1)) * (tau - A1(k1,tau) - A1(k1,tau) +
> A1(2*k1,tau))
> v312 = ((s1*s2*r21)/(k1*k2)) * (tau - A1(k1,tau) - A1(k2,tau) +
> A1(k1+k2,tau))
> v313 = ((s1*s3*r31)/(k1*k3)) * (tau - A1(k1,tau) - A1(k3,tau) +
> A1(k1+k3,tau))
>
> v322 = ((s2*s2)/(k2*k2)) * (tau - A1(k2,tau) - A1(k2,tau) +
> A1(2*k2,tau))
> v321 = ((s1*s2*r21)/(k1*k2)) * (tau - A1(k1,tau) - A1(k2,tau) +
> A1(k1+k2,tau))
> v323 = ((s2*s3*r32)/(k2*k3)) * (tau - A1(k2,tau) - A1(k3,tau) +
> A1(k2+k3,tau))
>
> v333 = ((s3*s3)/(k3*k3)) * (tau - A1(k3,tau) - A1(k3,tau) +
> A1(2*k3,tau))
> v331 = ((s1*s3*r31)/(k1*k3)) * (tau - A1(k3,tau) - A1(k1,tau) +
> A1(k3+k1,tau))
> v323 = ((s2*s3*r32)/(k2*k3)) * (tau - A1(k2,tau) - A1(k3,tau) +
> A1(k2+k3,tau))
>
> v = v1+v2+0.5*(v311+v312+v313+v322+v321+v323+v333+v331+v323)
>
> d = -v
>
> A = matrix(c(0,0,0,0,0,0,0,0,0),ncol=3)
> A[1,1] = 1-k1*freq
> A[1,2] = 0
> A[1,3] = 0
> A[2,1] = 0
> A[2,2] = 1-k2*freq
> A[2,3] = 0
> A[3,1] = 0
> A[3,2] = 0
> A[3,3] = 1-k3*freq
>
> ct = c(-lambda1*freq,-lambda2*freq,-lambda3*freq)
>
> Q = matrix(c(0,0,0,0,0,0,0,0,0),ncol=3)
> Q[1,1] = s1^2
> Q[1,2] = s1*s2*r21
> Q[1,3] = s1*s3*r31
> Q[2,1] = s2*s1*r21
> Q[2,2] = s2^2
> Q[2,3] = s2*s3*r32
> Q[3,1] = s1*s3*r31
> Q[3,2] = s2*s3*r32
> Q[3,3] = s3^2
>
> ## 4. Compute the State Space Form
>
> mDelta = matrix(c(ct,d),ncol=1)
>
> mPhi = rbind(A,H)
>
> mOmega1 = matrix(c(Q,0,0,0,0,0,0,0,0,0,0,0,0),ncol=7)
> mOmega2 =
> matrix(c(0,0,0,0,0,0,0,0,0,0,0,0,diag(c(st1,st2,st3,st4))),ncol=7)
>
> mOmega = rbind(mOmega1,mOmega2)
>
> mSigma = matrix(c(0,0,0,0,0,0,0,0,0,0,0,0),ncol=3)
>
> mSigma[1,1] = ((s1^2)/(2*k1))*(1-exp(-2*k1*freq))
> mSigma[2,2] = ((s2^2)/(2*k2))*(1-exp(-2*k2*freq))
> mSigma[3,3] = ((s3^2)/(2*k3))*(1-exp(-2*k3*freq))
>
> ## 7. Return State Space Form
>
> ssf.mod = list(mDelta=mDelta, mPhi=mPhi, mOmega=mOmega)
> CheckSsf(ssf.mod)
>
> }
>
> start = c(0.1,0.1,0.1,0.3,0.3,0.3,0.06,0.005,0.005,0.005,0.001,0.001,0.001
> ,
> 0.003,0.1,0.003,0.01)
>
> names(start)=c("logkappa1","logkappa2","logkappa3","lambda1","lambda2","lambda3","logDelta","logSigma1","logSigma2","logSigma3","Corr12","Corr13","Corr23","logSigmat1","logSigmat2","logSigmat3","logSigmat4")
>
> tau = c(0.25,0.5,1,5)
>
> l = c(0,0,0,0,0,0,0,0,0,-0,-1,-1,-1,0,0,0,0)
> u = c(10,10,10,10,10,10,0.5,10,10,10,1,1,1,10,10,10,10)
>
>
> ans.vasicek = SsfFit(start, fama.bliss, genvasicek.ssf, tau=tau,
> freq=1/12,
> trace=T, lower=l, upper=u, control=nlminb.control(abs.tol=1e-6,,
> rel.tol=1e-6, x.tol=1e-6, eval.max=1000,iter.max=135))
>
>
> On 8/13/07, Bernardo Ribeiro <[EMAIL PROTECTED]> wrote:
> >
> > Hey all,
> >
> > I am trying to work under a State Space form, but I didn't get the help
> > exactly.
> > Have anyone eles used this functions?
> >
> > I was used to work with S-PLUS, but I have some codes I need to adpt.
> >
> > Thanks alot,
> >
> > Bernardo
> >
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help@stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
[[alternative HTML version deleted]]
______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
