diff options
Diffstat (limited to 'R/estUtil.R')
| -rw-r--r-- | R/estUtil.R | 27 |
1 files changed, 16 insertions, 11 deletions
diff --git a/R/estUtil.R b/R/estUtil.R index 72b5143..21e132f 100644 --- a/R/estUtil.R +++ b/R/estUtil.R @@ -17,7 +17,7 @@ levbmqdt <- function(...,obj,theta0,N,opt){ update <- 1 # variable to count the number of times objective function is called countObj <- 0 - sumSqRatio <- 1 + sumSqDiff <- 1 repeat{ i=i+1 @@ -27,20 +27,20 @@ levbmqdt <- function(...,obj,theta0,N,opt){ repeat{ # Update Parameters H <- t(l$grad)%*%l$grad + d*diag(dim(theta0)[1]) - Hinv <- solve(H) - theta <- theta0 + Hinv%*%t(l$grad)%*%e - + Hinv <- solve(H); g <- t(l$grad)%*%e + theta <- theta0 + Hinv%*%g + # Update residuals e <- l$Y-l$X%*%theta sumsq <- sum(e^2) - sumSqRatio <- (sumsq0-sumsq)/sumsq0 + sumSqDiff <- sumsq0-sumsq countObj <- countObj + 1 - if(abs(sumSqRatio) < tol) break + if(abs(sumSqDiff)/100 < tol) break # If sum square error with the updated parameters is less than the # previous one, the updated parameters become the current parameters # and the damping coefficient is reduced by a factor of mu - if(sumSqRatio > 0){ + if(sumSqDiff > 0){ d <- d/mu theta0 <- theta sumsq0 <- sumsq @@ -50,10 +50,10 @@ levbmqdt <- function(...,obj,theta0,N,opt){ } } - if((abs(sumSqRatio) < tol)||(i == maxIter)) break + if((abs(sumSqDiff)/100 < tol)||(i == maxIter)) break } - if(abs(sumSqRatio) < tol){ + if(abs(sumSqDiff)/100 < tol){ WhyStop <- "Tolerance" } else{ WhyStop <- "Maximum Iteration Limit" @@ -61,7 +61,7 @@ levbmqdt <- function(...,obj,theta0,N,opt){ e <- e[1:N,] sigma2 <- sum(e^2)/df - vcov <- Hinv + vcov <- Hinv*sigma2 list(params=theta,residuals=e,vcov=vcov,sigma = sqrt(sigma2), termination=list(WhyStop=WhyStop,iter=i,FcnCount = countObj)) @@ -81,7 +81,7 @@ levbmqdt <- function(...,obj,theta0,N,opt){ #' @param LMstep Size of the Levenberg-Marquardt step #' #' @export -optimOptions <- function(tol=1e-3,maxIter=20,LMinit=0.1,LMstep=2){ +optimOptions <- function(tol=1e-4,maxIter=20,LMinit=0.1,LMstep=2){ return(list(tol=tol,maxIter= maxIter, adv= list(LMinit=LMinit, LMstep=LMstep))) } @@ -196,4 +196,9 @@ bjGrad <- function(theta,e,dots){ } return(l) +} + +inv <- function(H){ + chdecomp <- chol(H) + chol2inv(H) }
\ No newline at end of file |