[Yuima-commits] r595 - pkg/yuima/R

[noreply at r-forge.r-project.org]

Author: eguchi
Date: 2017-03-27 07:53:52 +0200 (Mon, 27 Mar 2017)
New Revision: 595

Added:
   pkg/yuima/R/IC.R
Log:
added

Added: pkg/yuima/R/IC.R
===================================================================
--- pkg/yuima/R/IC.R	                        (rev 0)
+++ pkg/yuima/R/IC.R	2017-03-27 05:53:52 UTC (rev 595)
@@ -0,0 +1,62 @@
+## information criteria
+
+IC <- function(yuima, data, start, lower, upper, joint = FALSE, rcpp = FALSE,...){
+  if(is(yuima,"yuima") == TRUE){
+    
+    if(is.matrix(data) == FALSE){
+      sub.zoo.data <- list(zoo(x = data, order.by = yuima at sampling@grid[[1]]))
+      names(sub.zoo.data)[1] <- "Series 1"
+    }else{
+      sub.zoo.data <- list()
+      for(i in 1:ncol(data)){
+        sub.zoo.data <- c(sub.zoo.data, list(zoo(x = data[,i], order.by = yuima at sampling@grid[[1]])))
+        names(sub.zoo.data)[i] <- paste("Series", i)
+      }
+    }
+    yuima at data@zoo.data <- sub.zoo.data
+    alpha <- yuima at model@parameter at drift
+    beta <- yuima at model@parameter at diffusion
+    Terminal <- yuima at sampling@Terminal
+    n <- yuima at sampling@n
+    
+    para.num.init  <- match(yuima at model@parameter at all, names(start))
+    para.num.low  <- match(yuima at model@parameter at all, names(lower))
+    para.num.upp  <- match(yuima at model@parameter at all, names(upper))
+    para.start <- NULL
+    para.lower <- NULL
+    para.upper <- NULL
+    for(j in 1:length(yuima at model@parameter at all)){
+      para.start <- c(para.start, list(start[[para.num.init[j]]]))
+      para.lower <- c(para.lower, list(lower[[para.num.low[j]]]))
+      para.upper <- c(para.upper, list(upper[[para.num.upp[j]]]))
+    }
+    names(para.start) <- yuima at model@parameter at all
+    names(para.lower) <- yuima at model@parameter at all
+    names(para.upper) <- yuima at model@parameter at all
+
+    mle <- qmle(yuima, start = para.start, lower = para.lower, upper = para.upper, method = "L-BFGS-B", joint = joint, rcpp = rcpp)
+    hess <- list(mle at details$hessian)
+    hess.diff <- subset(hess[[1]], rownames(hess[[1]])%in%beta, select=beta)
+    hess.drif <- subset(hess[[1]], rownames(hess[[1]])%in%alpha, select=alpha)
+    
+    esti <- coef(mle)
+    names(esti) <- c(beta, alpha)
+    cic <- summary(mle)@m2logL+2*(length(alpha)+length(beta))
+    bic <- summary(mle)@m2logL+length(alpha)*log(Terminal[1])+length(beta)*log(n[1])
+    if(det(hess.diff) > 0 && det(hess.drif) > 0){
+      qbic <- summary(mle)@m2logL+log(det(hess.diff))+log(det(hess.drif))
+    }else{
+      warning(cat("det(hessian of diffusion coefficient)<=0 or det(hessian of drift coefficient)<=0 \n"))
+      qbic <- summary(mle)@m2logL+length(alpha)*log(Terminal)+length(beta)*log(n)
+    }
+    
+    final.res <- list(par = esti, BIC = bic, QBIC = qbic, CIC = cic)
+    return(final.res)
+  }else{
+    yuima.stop("This function is for yuima-class.")
+  }
+}
+
+
+
+


Property changes on: pkg/yuima/R/IC.R
___________________________________________________________________
Added: svn:executable
   + *