I was trying to use mle()
to minimize the -log likelihood function, however it always pops out an error message like this:
Error in minuslogl() (from kf.r@7524iUV#2) : argument "psi" is missing, with
no default
My functions is:
KFlogL2 <- function(psi){
.........
return (res)
}
> psi.test
[1] 2.00 0.20 0.20 0.02 0.20 0.20 0.50
> mle(KFlogL2,start = list(psi = psi.test),method = 'L-BFGS-B',lower =
> lb,upper = ub)
> mle(KFlogL2,start = list(psi = c(2,.2,.2,.02,.2,.2,.5)),method =
> 'L-BFGS-B',lower = lb,upper = ub)
I tried both ways above but neither works...
Anyone knows why?
For those who would like to introduce optim()
, I've actually tried that,
however my Hessian matrix, after taking inverse and extracting diagonals -
turns to be negative! The convergence is fine. I have no idea what that
means. So I wanted to try with mle()
to see if things can change.
My Optimization function:
KFopt <- optim(par = psi.test, fn = KFlogL, method = 'L-BFGS-B',
lower = lb, upper = ub, hessian = TRUE)
My KF function:
KFlogL2 <- function(psi){
k <- psi[1]
sigmax <- psi[2]
lambdax <- psi[3]
mu <- psi[4]
sigmae <- psi[5]
rnmu <- psi[6]
pxe <- psi[7]
m <- length(init.state)
N <- ncol(y.use)
nobs <- nrow(y.use)
s <- rep(.01,N)
dt <- 7/360
cc <- c(0,mu * dt)
T <- diag(m)
T[1,1] <- exp(-k * dt)
xx <- (1-exp(-2 * k * dt)) * sigmax * sigmax / 2 / k
xy <- (1-exp(-k * dt)) * pxe * sigmax * sigmae / k
yx <- (1-exp(-k * dt)) * pxe * sigmax * sigmae / k
yy <- sigmae * sigmae * dt
Q <- matrix(c(xx,xy,yx,yy),m,m)
R <- diag(m)
# measurement equation
d <- rep(0,N)
Z <- matrix(0,N,m)
for (i in 1:N){
matur.i <- matur.use[i]
p1 <- (1-exp(-2 * k * matur.i)) * sigmax * sigmax / 2 / k
p2 <- sigmae * sigmae * matur.i
p3 <- 2 * (1-exp(k * matur.i)) * pxe * sigmax * sigmae / k
d[i] <- rnmu * matur.i - (1-exp(-k * matur.i)) * lambdax / k + 1/2 * (p1 + p2 + p3)
Z[i,] <- exp(-k * matur.i)
}
H <- diag(s)
# Kalman Filter
save.ytt1 <- save.vtt <- save.vt <- matrix(0,nobs,N)
save.att1 <- save.att <- matrix(0,nobs,m)
save.Ptt1 <- save.Ptt <- matrix(0,nobs,m*m)
save.Ftt1 <- matrix(0,nobs,N*N)
save.detF <- save.vFv <- rep(0,nobs)
Ptt <- init.dist
att <- init.state
for (iter in 1:nobs){
Ptt1 <- T %*% Ptt %*% t(T) + R %*% Q %*% t(R)
Ftt1 <- Z %*% Ptt1 %*% t(Z) + H
att1 <- T %*% att + cc
yt <- y.use[iter,]
ytt1 <- Z %*% att1 + d
vt <- yt - ytt1
temp <- Ptt1 %*% crossprod(Z,solve(Ftt1))
att <- att1 + temp %*% vt
Ptt <- Ptt1 - temp %*% Z %*% Ptt1
ytt <- Z %*% att + d
vtt <- yt - ytt
save.vtt[iter,] <- vtt
save.vt[iter,] <- vt
save.att[iter,] <- att
save.Ptt1[iter,] <- as.numeric(Ptt1)
save.Ptt[iter,] <- as.numeric(Ptt)
save.detF[iter] <- det(Ftt1)
save.vFv[iter] <- crossprod(vt,solve(Ftt1,vt))
}
logL <- - (N + 1) * nobs/2 * log(2 * pi) - 1/2 * sum(log(save.detF)) - 1/2 * sum(save.vFv)
res <- -logL
return (res)
}
I don't know how to attach a file here, so here is just a glimpse:
3.130263167 3.130700134 3.058707073 3.012589391 2.999724295 2.991724252
3.102342009 3.09421922 2.999724295 2.952824773 2.940747965 2.932259851
3.118392286 3.125882958 3.006177531 2.949164638 2.926917958 2.913979772
3.106378794 3.072693315 2.991724252 2.950211758 2.937573359 2.926917958
3.104138147 3.111735949 3.031099417 2.991724252 2.979602892 2.969388298
3.107273648 3.113959655 3.033509638 2.992226134 2.979602892 2.969388298
3.085572978 3.086486637 3.059176446 3.033509638 3.026746327 3.020912572
3.100092289 3.099641737 3.064325065 3.033028058 3.025776395 3.019936962
3.073618812 3.07176696 3.061520014 3.054001182 3.052112607 3.050220459
3.058707073 2.961658293 3.054001182 3.039270576 3.03543364 3.032546247
3.003204288 3.006672214 3.016024977 3.007660844 3.003700443 3.001217204
2.96217549 2.959068289 3.013080912 3.017493765 3.017493765 3.017493765
3.013572192 3.0194488 3.022860941 3.014554028 3.006672214 3.004692015
3.007166651 3.010620886 3.042138646 3.038791763 3.033028058 3.0301337
2.865623588 2.872434057 2.972463647 2.994731773 2.998727783 3.002707887
2.85474458 2.858766418 2.983659692 2.989714201 2.987700102 2.985681938
2.882003508 2.944438979 2.999226163 2.999226163 2.99723115 2.995232149
2.931726944 2.938103161 3.012097628 3.006672214 3.009635179 3.00914196
2.904712875 2.905807566 2.983659692 2.984165637 2.983659692 2.983659692
2.971439581 2.975019232 3.007660844 3.001217204 2.996731774 2.993730271
2.862772146 2.872434057 2.993229143 2.998727783 2.994731773 2.991724252
2.890371758 2.889260029 2.992727765 3.005187432 3.003700443 3.003700443
2.797890905 2.814210397 2.930126516 2.959586827 2.97603964 2.986691529
2.855895328 2.862772146 2.947591898 2.965788397 2.973997781 2.980110893
2.744060639 2.750470917 2.913437031 2.949164638 2.964241606 2.975019232
2.840247371 2.841414913 2.936512914 2.962692419 2.97603964 2.986186861
2.817203515 2.821378886 2.915606229 2.935982269 2.947067102 2.957511061
2.836150204 2.829677689 2.938103161 2.961140829 2.973997781 2.986186861
2.903068589 2.925846146 3.032546247 3.039749159 3.043569603 3.045474365
2.987700102 3.023347441 3.060583246 3.060114532 3.060114532 3.060114532
3.023833704 3.029650492 3.085572978 3.079153882 3.078693794 3.078693794
3.387774361 3.343215099 3.273742726 3.23317313 3.201119103 3.168003494
3.284663565 3.274121299 3.221272949 3.185939325 3.16463081 3.147594623
3.348499593 3.346741196 3.298426104 3.244153632 3.212455257 3.195811885
3.32251486 3.327909586 3.24998682 3.191710157 3.162940193 3.144152279
3.377587516 3.371425223 3.291010423 3.224857897 3.191710157 3.173459961
3.410817625 3.426215145 3.308716529 3.229222117 3.193763124 3.176803048
3.510948246 3.508855256 3.390473418 3.300640127 3.259249719 3.240245851
3.646232879 3.625140613 3.476614021 3.354455119 3.292126287 3.258865473
3.538928277 3.524888854 3.454738149 3.328626689 3.249211025 3.206803244
3.706964922 3.698829785 3.541828511 3.400863993 3.316365446 3.272606147
3.672241813 3.660737148 3.507357577 3.375537635 3.296207168 3.25617161
3.399529325 3.379973745 3.237501289 3.141994781 3.101442728 3.075928816
3.553346059 3.542118073 3.380994674 3.270329106 3.213662258 3.176803048
3.478467017 3.485232111 3.325395668 3.223266173 3.165475048 3.127199036
3.505557397 3.50013733 3.348850901 3.247268899 3.195402469 3.162093811
3.384390263 3.364187556 3.256556892 3.185112195 3.143289838 3.118834471
3.498021566 3.492256113 3.348148161 3.259633817 3.210843653 3.181796817
3.369018483 3.422958873 3.295466427 3.211246798 3.171364842 3.148024084
3.277144733 3.273742726 3.169265324 3.113515309 3.097837496 3.08876714