Problem solved.
It had something to do with calling expm(-array(c(K_1, 0, 0, K_2),
c(2,2))*h).
2011/9/22 Kristian Lind <kristian.langgaard.l...@gmail.com>
> Dear R users,
>
> When running the program below I receive the following error message:
> fit <- optim(parm, objective, yt = tyield, hessian = TRUE)
> Error in as.vector(data) :
> no method for coercing this S4 class to a vector
>
> I can't figure out what the problem is exactly. I imagine that it has
> something to do with "tyield" being a matrix. Any help on explaining what's
> going on and how to solve this is much appreciated.
>
> Thank you,
>
> Kristian
>
> library(FKF) #loading Fast Kalman Filter package
> library(Matrix) # matrix exponential package
>
> K_1 = 0.1156
> K_2 = 0.17
> sigma_1 = 0.1896
> sigma_2 = 0.2156
> lambda_1 = 0
> lambda_2 = -0.5316
> theta_1 = 0.1513
> theta_2 = 0.2055
>
> #test data
> tyield <- matrix(data = rnorm(200), nrow =2, ncol =100)
>
> # defining dimensions
> m <- 2 # m is the number of state variables
> n <- 100 # is the length of the observed sample
> d <- 2 # is the number of observed variables.
> theta <- c(theta_1, theta_2)
> h <- t <- 1/52 # time between observations
>
> ## creating state space representation of 2-factor CIR model follwing
> Driessen and Geyer et al.
> CIR2ss <- function(K_1, K_2, sigma_1, sigma_2, lambda_1, lambda_2, theta_1,
> theta_2) {
> ## defining auxilary parameters
> phi_11 <- sqrt((K_1+lambda_1)^2+2*sigma_1^2)
> phi_21 <- sqrt((K_1+lambda_2)^2+2*sigma_2^2)
> phi_12 <- K_1+lambda_1+phi_11
> phi_22 <- K_2+lambda_2+phi_12
> phi_13 <- -2*K_1*theta_1/sigma_1^2
> phi_23 <- -2*K_2*theta_2/sigma_2^2
> phi_14 <- 2*phi_11+phi_21*(exp(phi_11*t)-1)
> phi_24 <- 2*phi_12+phi_22*(exp(phi_12*t)-1)
> phi <- array(c(phi_11, phi_21, phi_12, phi_22, phi_13, phi_23, phi_14,
> phi_24), c(4,2))
> a <- array(0, c(d,n))
> for(t in n:1){
> a[,n-(t+1)] <-
> -phi_13/(n-(t+1))*log(2*phi_11*exp(phi_12*(n-(t+1))/2)/phi_14)-phi_23/(n-(t+1))*log(2*phi_21*exp(phi_22*(n-(t+1))/2)/phi_24)
> }
> b <- array(c(1,0,0,1,0), c(d,m,n))
> j <- -array(c(K_1, 0, 0, K_2), c(2,2))*h
> explh <- expm(j)
> Tt <- array(explh, c(m,m,n)) #array giving the factor of the transition
> equation
> Zt <- b #array giving the factor of the measurement equation
> ct <- a #matrix giving the intercept of the measurement equation
> dt <- (diag(m)-expm(-array(c(K_1, 0, 0, K_2), c(2,2))*h))*theta #matrix
> giving the intercept of the transition equation
> GGt <- array(c(1,0,0,1), c(d,d,n)) #array giving the variance of the
> disturbances of the measurement equation
> HHt <- array(c(1,0,0,1), c(m,m,n)) #array giving the variance of the
> innovations of the transition equation
> a0 <- c(0, 0) #vector giving the initial value/estimation of the state
> variable
> P0 <- matrix(1e6, nrow = 2, ncol = 2) # matrix giving the variance of
> a0
> return(list(a0 = a0, P0 = P0, ct = ct, dt = dt, Zt = Zt, Tt = Tt, GGt =
> GGt,
> HHt = HHt))
> }
>
> ## Objective function passed to optim
> objective <- function(parm, yt) {
> sp <- CIR2ss(parm["K_1"], parm["K_2"], parm["sigma_1"], parm["sigma_2"],
> parm["lambda_1"], parm["lambda_2"],
> parm["theta_1"], parm["theta_2"])
> ans <- fkf(a0 = sp$a0, P0 = sp$P0, dt = sp$dt, ct = sp$ct, Tt = sp$Tt,
> Zt = sp$Zt, HHt = sp$HHt, GGt = sp$GGt, yt = yt)
> return(-ans$loglik)
> }
>
> parm <- c(K_1 = 0.1156, K_2 = 0.17, sigma_1 = 0.1896, sigma_2 = 0.2156,
> lambda_1 = 0, lambda_2 = -0.5316, theta_1 = 0.1513, theta_2 =
> 0.2055) # initial parameters
>
> ##optimizing objective function
> fit <- optim(parm, objective, yt = tyield, hessian = TRUE)
> print(fit)
>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org 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.
