[Yuima-commits] r353 - in pkg/yuima: . R man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Nov 17 00:30:52 CET 2014
Author: lorenzo
Date: 2014-11-17 00:30:52 +0100 (Mon, 17 Nov 2014)
New Revision: 353
Modified:
pkg/yuima/DESCRIPTION
pkg/yuima/R/setCarma.R
pkg/yuima/R/setCogarch.R
pkg/yuima/R/simulate.R
pkg/yuima/R/yuima.R
pkg/yuima/man/setCarma.Rd
pkg/yuima/man/setCogarch.Rd
Log:
Show method for Cogarch Model
Modified: pkg/yuima/DESCRIPTION
===================================================================
--- pkg/yuima/DESCRIPTION 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/DESCRIPTION 2014-11-16 23:30:52 UTC (rev 353)
@@ -1,8 +1,8 @@
Package: yuima
Type: Package
Title: The YUIMA Project package for SDEs
-Version: 1.0.49
-Date: 2014-11-11
+Version: 1.0.50
+Date: 2014-11-17
Depends: methods, zoo, stats4, utils, expm, cubature, mvtnorm
Author: YUIMA Project Team
Maintainer: Stefano M. Iacus <stefano.iacus at unimi.it>
Modified: pkg/yuima/R/setCarma.R
===================================================================
--- pkg/yuima/R/setCarma.R 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/R/setCarma.R 2014-11-16 23:30:52 UTC (rev 353)
@@ -73,6 +73,7 @@
Carma.var="v",
Latent.var="x",
XinExpr=FALSE,
+ Cogarch=FALSE,
...){
# We use the same parametrization as in Brockwell (2000)
@@ -142,7 +143,11 @@
}else{
if(XinExpr==TRUE){
Int.Var<-paste(Latent.var,"0",sep="")
- mydots$xinit<-paste(Int.Var,c(0:(p-1)),sep="")
+ if(Cogarch==FALSE){
+ mydots$xinit<-paste(Int.Var,c(0:(p-1)),sep="")
+ }else{
+ mydots$xinit<-paste(Int.Var,c(1:p),sep="")
+ }
}
}
} else{
@@ -156,10 +161,17 @@
beta_coeff0<-paste("-",ar.par,sep="")
beta_coeff<-paste(beta_coeff0,p:1,sep="")
- coeff_alpha<-c(paste(ma.par1,0:q,sep=""),as.character(matrix(0,1,p-q-1)))
+ if(Cogarch==FALSE){
+ coeff_alpha<-c(paste(ma.par1,0:q,sep=""),as.character(matrix(0,1,p-q-1)))
+ }else{
+ coeff_alpha<-c(paste(ma.par1,1:(q+1),sep=""),as.character(matrix(0,1,p-q-1)))
+ }
fin_alp<-length(coeff_alpha)
-
- Y_coeff<-paste(Latent.var,0:(p-1),sep="")
+ if(Cogarch==FALSE){
+ Y_coeff<-paste(Latent.var,0:(p-1),sep="")
+ }else{
+ Y_coeff<-paste(Latent.var,1:p,sep="")
+ }
fin_Y<-length(Y_coeff)
V1<-paste(coeff_alpha,Y_coeff,sep="*")
V2<-paste(V1,collapse="+")
Modified: pkg/yuima/R/setCogarch.R
===================================================================
--- pkg/yuima/R/setCogarch.R 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/R/setCogarch.R 2014-11-16 23:30:52 UTC (rev 353)
@@ -70,9 +70,9 @@
setCogarch<-function(p,
q,
- ar.par="a",
- ma.par="b",
- loc.par="mu",
+ ar.par="b",
+ ma.par="a",
+ loc.par="a0",
Cogarch.var="g",
V.var="v",
Latent.var="y",
@@ -109,7 +109,11 @@
ar.par=ar.par,
ma.par=ma.par,
loc.par=loc.par,
- lin.par=ma.par)
+ lin.par=ma.par,
+ Cogarch=TRUE)
+ # In order to have a representation of a Cogarch(p,q) model coherent with the
+ # Chaadra Brockwell and Davis we need to modify the slot xinit and drift[1]
+
}else{
if(!is.null(mydots$xinit)){
aux.Carma<-setCarma(p=q,
@@ -118,10 +122,13 @@
ar.par=ar.par,
ma.par=ma.par,
lin.par=ma.par,
- xinit=mydots$xinit)
+ xinit=mydots$xinit,
+ Cogarch=TRUE)
+ # In order to have a representation of a Cogarch(p,q) model coherent with the
+ # Chaadra Brockwell and Davis we need to modify the slot xinit and drift[1]
}
}
-
+
newdrift<-c(0,as.character(aux.Carma at drift))
newdiffusion<-c(0,as.character(eval(aux.Carma at diffusion)))
line1<-c(paste0("sqrt(",V.var,")"),as.character(matrix(0,nrow=(q+1),ncol=1)))
Modified: pkg/yuima/R/simulate.R
===================================================================
--- pkg/yuima/R/simulate.R 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/R/simulate.R 2014-11-16 23:30:52 UTC (rev 353)
@@ -52,6 +52,13 @@
return(out)
})
+
+# We rewrite the code of simulate method. We build a new internal function aux.simulate containing
+# the old code. This function starts directly if the model is an object of yuima.model-class
+# or yuima.carma-class. If the model is an object of class yuima.cogarch, simulate method runs
+# the internal function aux.simulateCogarch that generates first a path of the underlying levy and then
+# uses this path to construct the trajectories of the Cogarch model by calling the aux.simulate function.
+
setMethod("simulate", "yuima",
function(object, nsim=1, seed=NULL, xinit, true.parameter,
space.discretized=FALSE, increment.W=NULL, increment.L=NULL,
@@ -61,188 +68,459 @@
grid = as.numeric(NULL), random = FALSE, sdelta=as.numeric(NULL),
sgrid=as.numeric(NULL), interpolation="none"){
- ##:: errors checks
-
- ##:1: error on yuima model
- yuima <- object
-
- if(missing(yuima)){
- yuima.warn("yuima object is missing.")
- return(NULL)
+ if(is(object at model,"yuima.cogarch")){
+ res<-aux.simulateCogarch(object, nsim, seed, xinit, true.parameter,
+ space.discretized, increment.W, increment.L,
+ hurst,methodfGn,
+ sampling, subsampling,
+ Initial, Terminal, n, delta,
+ grid, random, sdelta,
+ sgrid, interpolation)
+ }else{
+ res<-aux.simulate(object, nsim, seed, xinit, true.parameter,
+ space.discretized, increment.W, increment.L,
+ hurst,methodfGn,
+ sampling, subsampling,
+ Initial, Terminal, n, delta,
+ grid, random, sdelta,
+ sgrid, interpolation)
}
+ return(res)
- tmphurst<-yuima at model@hurst
-
- if(!missing(hurst)){
- yuima at model@hurst=hurst
- }
-
- if (is.na(yuima at model@hurst)){
- yuima.warn("Specify the hurst parameter.")
- return(NULL)
- }
-
- tmpsamp <- NULL
- if(is.null(yuima at sampling)){
- if(missing(sampling)){
- tmpsamp <- setSampling(Initial = Initial,
- Terminal = Terminal, n = n, delta = delta,
- grid = grid, random = random, sdelta=sdelta,
- sgrid=sgrid, interpolation=interpolation)
- } else {
- tmpsamp <- sampling
- }
- } else {
- tmpsamp <- yuima at sampling
- }
-
- yuima at sampling <- tmpsamp
-
- sdeModel <- yuima at model
- Terminal <- yuima at sampling@Terminal[1]
- n <- yuima at sampling@n[1]
- r.size <- sdeModel at noise.number
- d.size <- sdeModel at equation.number
-
- ##:2: error on xinit
- if(missing(xinit)){
- xinit <- sdeModel at xinit
- } else {
- if(length(xinit) != d.size){
- if(length(xinit)==1){
- xinit <- rep(xinit, d.size)
- } else {
- yuima.warn("Dimension of xinit variables missmatch.")
- return(NULL)
- }
- }
- }
-
- xinit <- as.expression(xinit) # force xinit to be an expression
-
+# ##:: errors checks
+#
+# ##:1: error on yuima model
+# yuima <- object
+#
+# if(missing(yuima)){
+# yuima.warn("yuima object is missing.")
+# return(NULL)
+# }
+#
+# tmphurst<-yuima at model@hurst
+#
+# if(!missing(hurst)){
+# yuima at model@hurst=hurst
+# }
+#
+# if (is.na(yuima at model@hurst)){
+# yuima.warn("Specify the hurst parameter.")
+# return(NULL)
+# }
+#
+# tmpsamp <- NULL
+# if(is.null(yuima at sampling)){
+# if(missing(sampling)){
+# tmpsamp <- setSampling(Initial = Initial,
+# Terminal = Terminal, n = n, delta = delta,
+# grid = grid, random = random, sdelta=sdelta,
+# sgrid=sgrid, interpolation=interpolation)
+# } else {
+# tmpsamp <- sampling
+# }
+# } else {
+# tmpsamp <- yuima at sampling
+# }
+#
+# yuima at sampling <- tmpsamp
+#
+# sdeModel <- yuima at model
+# Terminal <- yuima at sampling@Terminal[1]
+# n <- yuima at sampling@n[1]
+# r.size <- sdeModel at noise.number
+# d.size <- sdeModel at equation.number
+#
+# ##:2: error on xinit
+# if(missing(xinit)){
+# xinit <- sdeModel at xinit
+# } else {
+# if(length(xinit) != d.size){
+# if(length(xinit)==1){
+# xinit <- rep(xinit, d.size)
+# } else {
+# yuima.warn("Dimension of xinit variables missmatch.")
+# return(NULL)
+# }
+# }
+# }
+#
+# xinit <- as.expression(xinit) # force xinit to be an expression
+#
+#
+# par.len <- length(sdeModel at parameter@all)
+#
+# if(missing(true.parameter) & par.len>0){
+# true.parameter <- vector(par.len, mode="list")
+# for(i in 1:par.len)
+# true.parameter[[i]] <- 0
+# names(true.parameter) <- sdeModel at parameter@all
+# }
+#
+# yuimaEnv <- new.env()
+#
+# if(par.len>0){
+# for(i in 1:par.len){
+# pars <- sdeModel at parameter@all[i]
+#
+# for(j in 1:length(true.parameter)){
+# if( is.na(match(pars, names(true.parameter)[j]))!=TRUE){
+# assign(sdeModel at parameter@all[i], true.parameter[[j]],yuimaEnv)
+# }
+# }
+# #assign(sdeModel at parameter@all[i], true.parameter[[i]], yuimaEnv)
+# }
+# }
+#
+#
+# if(space.discretized){
+# if(r.size>1){
+# warning("Space-discretized EM cannot be used for multi-dimentional models. Use standard method.")
+# space.discretized <- FALSE
+# }
+# if(length(sdeModel at jump.coeff)){
+# warning("Space-discretized EM is for only Wiener Proc. Use standard method.")
+# space.discretized <- FALSE
+# }
+# }
+#
+# ##:: Error check for increment specified version.
+# if(!missing(increment.W) & !is.null(increment.W)){
+# if(space.discretized == TRUE){
+# yuima.warn("Parameter increment must be invalid if space.discretized=TRUE.")
+# return(NULL)
+# }else if(dim(increment.W)[1] != r.size){
+# yuima.warn("Length of increment's row must be same as yuima at model@noise.number.")
+# return(NULL)
+# }else if(dim(increment.W)[2] != n){
+# yuima.warn("Length of increment's column must be same as sampling at n[1].")
+# return(NULL)
+# }
+# }
+#
+# ##:: Error check for increment specified version.
+# if(!missing(increment.L) & !is.null(increment.L)){
+# if(space.discretized == TRUE){
+# yuima.warn("Parameter increment must be invalid if space.discretized=TRUE.")
+# return(NULL)
+# }else if(dim(increment.L)[1] != length(yuima at model@jump.coeff[[1]]) ){ #r.size){
+# yuima.warn("Length of increment's row must be same as yuima at model@noise.number.")
+# return(NULL)
+# }else if(dim(increment.L)[2] != n){
+# yuima.warn("Length of increment's column must be same as sampling at n[1].")
+# return(NULL)
+# }
+# }
+#
+#
+# yuimaEnv$dL <- increment.L
+#
+#
+# if(space.discretized){
+# ##:: using Space-discretized Euler-Maruyama method
+# yuima at data <- space.discretized(xinit, yuima, yuimaEnv)
+#
+# yuima at model@hurst<-tmphurst
+# return(yuima)
+# }
+#
+#
+# ##:: using Euler-Maruyama method
+# delta <- Terminal/n
+#
+# if(missing(increment.W) | is.null(increment.W)){
+#
+# if( sdeModel at hurst!=0.5 ){
+#
+# grid<-sampling2grid(yuima at sampling)
+# isregular<-yuima at sampling@regular
+#
+# if((!isregular) || (methodfGn=="Cholesky")){
+# dW<-CholeskyfGn(grid, sdeModel at hurst,r.size)
+# yuima.warn("Cholesky method for simulating fGn has been used.")
+# } else {
+# dW<-WoodChanfGn(grid, sdeModel at hurst,r.size)
+# }
+#
+# } else {
+#
+# delta<-Terminal/n
+# if(!is.Poisson(sdeModel)){ # if pure CP no need to setup dW
+# dW <- rnorm(n * r.size, 0, sqrt(delta))
+# dW <- matrix(dW, ncol=n, nrow=r.size,byrow=TRUE)
+# } else {
+# dW <- matrix(0,ncol=n,nrow=1) # maybe to be fixed
+# }
+# }
+#
+# } else {
+# dW <- increment.W
+# }
+#
+# if(is.Poisson(sdeModel)){
+# yuima at data <- simCP(xinit, yuima, yuimaEnv)
+# } else {
+# yuima at data <- euler(xinit, yuima, dW, yuimaEnv)
+# }
+#
+# for(i in 1:length(yuima at data@zoo.data))
+# index(yuima at data@zoo.data[[i]]) <- yuima at sampling@grid[[1]] ## to be fixed
+# yuima at model@xinit <- xinit
+# yuima at model@hurst <-tmphurst
+#
+# if(missing(subsampling))
+# return(yuima)
+# subsampling(yuima, subsampling)
+#
+ }
+ )
- par.len <- length(sdeModel at parameter@all)
-
- if(missing(true.parameter) & par.len>0){
- true.parameter <- vector(par.len, mode="list")
- for(i in 1:par.len)
- true.parameter[[i]] <- 0
- names(true.parameter) <- sdeModel at parameter@all
- }
-
- yuimaEnv <- new.env()
-
- if(par.len>0){
- for(i in 1:par.len){
- pars <- sdeModel at parameter@all[i]
-
- for(j in 1:length(true.parameter)){
- if( is.na(match(pars, names(true.parameter)[j]))!=TRUE){
- assign(sdeModel at parameter@all[i], true.parameter[[j]],yuimaEnv)
- }
- }
- #assign(sdeModel at parameter@all[i], true.parameter[[i]], yuimaEnv)
- }
- }
-
-
- if(space.discretized){
- if(r.size>1){
- warning("Space-discretized EM cannot be used for multi-dimentional models. Use standard method.")
- space.discretized <- FALSE
- }
- if(length(sdeModel at jump.coeff)){
- warning("Space-discretized EM is for only Wiener Proc. Use standard method.")
- space.discretized <- FALSE
- }
- }
-
- ##:: Error check for increment specified version.
- if(!missing(increment.W) & !is.null(increment.W)){
- if(space.discretized == TRUE){
- yuima.warn("Parameter increment must be invalid if space.discretized=TRUE.")
- return(NULL)
- }else if(dim(increment.W)[1] != r.size){
- yuima.warn("Length of increment's row must be same as yuima at model@noise.number.")
- return(NULL)
- }else if(dim(increment.W)[2] != n){
- yuima.warn("Length of increment's column must be same as sampling at n[1].")
- return(NULL)
- }
- }
-
- ##:: Error check for increment specified version.
- if(!missing(increment.L) & !is.null(increment.L)){
- if(space.discretized == TRUE){
- yuima.warn("Parameter increment must be invalid if space.discretized=TRUE.")
- return(NULL)
- }else if(dim(increment.L)[1] != length(yuima at model@jump.coeff[[1]]) ){ #r.size){
- yuima.warn("Length of increment's row must be same as yuima at model@noise.number.")
- return(NULL)
- }else if(dim(increment.L)[2] != n){
- yuima.warn("Length of increment's column must be same as sampling at n[1].")
- return(NULL)
- }
- }
-
-
- yuimaEnv$dL <- increment.L
-
-
- if(space.discretized){
- ##:: using Space-discretized Euler-Maruyama method
- yuima at data <- space.discretized(xinit, yuima, yuimaEnv)
-
- yuima at model@hurst<-tmphurst
- return(yuima)
- }
-
-
- ##:: using Euler-Maruyama method
- delta <- Terminal/n
-
- if(missing(increment.W) | is.null(increment.W)){
-
- if( sdeModel at hurst!=0.5 ){
-
- grid<-sampling2grid(yuima at sampling)
- isregular<-yuima at sampling@regular
-
- if((!isregular) || (methodfGn=="Cholesky")){
- dW<-CholeskyfGn(grid, sdeModel at hurst,r.size)
- yuima.warn("Cholesky method for simulating fGn has been used.")
- } else {
- dW<-WoodChanfGn(grid, sdeModel at hurst,r.size)
- }
-
- } else {
-
- delta<-Terminal/n
- if(!is.Poisson(sdeModel)){ # if pure CP no need to setup dW
- dW <- rnorm(n * r.size, 0, sqrt(delta))
- dW <- matrix(dW, ncol=n, nrow=r.size,byrow=TRUE)
- } else {
- dW <- matrix(0,ncol=n,nrow=1) # maybe to be fixed
- }
- }
-
- } else {
- dW <- increment.W
- }
-
- if(is.Poisson(sdeModel)){
- yuima at data <- simCP(xinit, yuima, yuimaEnv)
- } else {
- yuima at data <- euler(xinit, yuima, dW, yuimaEnv)
- }
-
- for(i in 1:length(yuima at data@zoo.data))
- index(yuima at data@zoo.data[[i]]) <- yuima at sampling@grid[[1]] ## to be fixed
- yuima at model@xinit <- xinit
- yuima at model@hurst <-tmphurst
-
- if(missing(subsampling))
- return(yuima)
- subsampling(yuima, subsampling)
-
- })
+aux.simulate<-function(object, nsim, seed, xinit, true.parameter,
+ space.discretized, increment.W, increment.L,
+ hurst,methodfGn,
+ sampling, subsampling,
+ Initial, Terminal, n, delta,
+ grid, random, sdelta,
+ sgrid, interpolation){
+
+ ##:: errors checks
+
+ ##:1: error on yuima model
+ yuima <- object
+
+ if(missing(yuima)){
+ yuima.warn("yuima object is missing.")
+ return(NULL)
+ }
+
+ tmphurst<-yuima at model@hurst
+
+ if(!missing(hurst)){
+ yuima at model@hurst=hurst
+ }
+
+ if (is.na(yuima at model@hurst)){
+ yuima.warn("Specify the hurst parameter.")
+ return(NULL)
+ }
+
+ tmpsamp <- NULL
+ if(is.null(yuima at sampling)){
+ if(missing(sampling)){
+ tmpsamp <- setSampling(Initial = Initial,
+ Terminal = Terminal, n = n, delta = delta,
+ grid = grid, random = random, sdelta=sdelta,
+ sgrid=sgrid, interpolation=interpolation)
+ } else {
+ tmpsamp <- sampling
+ }
+ } else {
+ tmpsamp <- yuima at sampling
+ }
+
+ yuima at sampling <- tmpsamp
+
+ sdeModel <- yuima at model
+ Terminal <- yuima at sampling@Terminal[1]
+ n <- yuima at sampling@n[1]
+ r.size <- sdeModel at noise.number
+ d.size <- sdeModel at equation.number
+
+ ##:2: error on xinit
+ if(missing(xinit)){
+ xinit <- sdeModel at xinit
+ } else {
+ if(length(xinit) != d.size){
+ if(length(xinit)==1){
+ xinit <- rep(xinit, d.size)
+ } else {
+ yuima.warn("Dimension of xinit variables missmatch.")
+ return(NULL)
+ }
+ }
+ }
+
+ xinit <- as.expression(xinit) # force xinit to be an expression
+
+
+ par.len <- length(sdeModel at parameter@all)
+
+ if(missing(true.parameter) & par.len>0){
+ true.parameter <- vector(par.len, mode="list")
+ for(i in 1:par.len)
+ true.parameter[[i]] <- 0
+ names(true.parameter) <- sdeModel at parameter@all
+ }
+
+ yuimaEnv <- new.env()
+
+ if(par.len>0){
+ for(i in 1:par.len){
+ pars <- sdeModel at parameter@all[i]
+
+ for(j in 1:length(true.parameter)){
+ if( is.na(match(pars, names(true.parameter)[j]))!=TRUE){
+ assign(sdeModel at parameter@all[i], true.parameter[[j]],yuimaEnv)
+ }
+ }
+ #assign(sdeModel at parameter@all[i], true.parameter[[i]], yuimaEnv)
+ }
+ }
+
+
+ if(space.discretized){
+ if(r.size>1){
+ warning("Space-discretized EM cannot be used for multi-dimentional models. Use standard method.")
+ space.discretized <- FALSE
+ }
+ if(length(sdeModel at jump.coeff)){
+ warning("Space-discretized EM is for only Wiener Proc. Use standard method.")
+ space.discretized <- FALSE
+ }
+ }
+
+ ##:: Error check for increment specified version.
+ if(!missing(increment.W) & !is.null(increment.W)){
+ if(space.discretized == TRUE){
+ yuima.warn("Parameter increment must be invalid if space.discretized=TRUE.")
+ return(NULL)
+ }else if(dim(increment.W)[1] != r.size){
+ yuima.warn("Length of increment's row must be same as yuima at model@noise.number.")
+ return(NULL)
+ }else if(dim(increment.W)[2] != n){
+ yuima.warn("Length of increment's column must be same as sampling at n[1].")
+ return(NULL)
+ }
+ }
+
+ ##:: Error check for increment specified version.
+ if(!missing(increment.L) & !is.null(increment.L)){
+ if(space.discretized == TRUE){
+ yuima.warn("Parameter increment must be invalid if space.discretized=TRUE.")
+ return(NULL)
+ }else if(dim(increment.L)[1] != length(yuima at model@jump.coeff[[1]]) ){ #r.size){
+ yuima.warn("Length of increment's row must be same as yuima at model@noise.number.")
+ return(NULL)
+ }else if(dim(increment.L)[2] != n){
+ yuima.warn("Length of increment's column must be same as sampling at n[1].")
+ return(NULL)
+ }
+ }
+
+
+ yuimaEnv$dL <- increment.L
+
+
+ if(space.discretized){
+ ##:: using Space-discretized Euler-Maruyama method
+ yuima at data <- space.discretized(xinit, yuima, yuimaEnv)
+
+ yuima at model@hurst<-tmphurst
+ return(yuima)
+ }
+
+
+ ##:: using Euler-Maruyama method
+ delta <- Terminal/n
+
+ if(missing(increment.W) | is.null(increment.W)){
+
+ if( sdeModel at hurst!=0.5 ){
+
+ grid<-sampling2grid(yuima at sampling)
+ isregular<-yuima at sampling@regular
+
+ if((!isregular) || (methodfGn=="Cholesky")){
+ dW<-CholeskyfGn(grid, sdeModel at hurst,r.size)
+ yuima.warn("Cholesky method for simulating fGn has been used.")
+ } else {
+ dW<-WoodChanfGn(grid, sdeModel at hurst,r.size)
+ }
+
+ } else {
+
+ delta<-Terminal/n
+ if(!is.Poisson(sdeModel)){ # if pure CP no need to setup dW
+ dW <- rnorm(n * r.size, 0, sqrt(delta))
+ dW <- matrix(dW, ncol=n, nrow=r.size,byrow=TRUE)
+ } else {
+ dW <- matrix(0,ncol=n,nrow=1) # maybe to be fixed
+ }
+ }
+
+ } else {
+ dW <- increment.W
+ }
+
+ if(is.Poisson(sdeModel)){
+ yuima at data <- simCP(xinit, yuima, yuimaEnv)
+ } else {
+ yuima at data <- euler(xinit, yuima, dW, yuimaEnv)
+ }
+
+ for(i in 1:length(yuima at data@zoo.data))
+ index(yuima at data@zoo.data[[i]]) <- yuima at sampling@grid[[1]] ## to be fixed
+ yuima at model@xinit <- xinit
+ yuima at model@hurst <-tmphurst
+
+ if(missing(subsampling))
+ return(yuima)
+ subsampling(yuima, subsampling)
+
+}
+
+aux.simulateCogarch<-function(object, nsim, seed, xinit, true.parameter,
+ space.discretized, increment.W, increment.L,
+ hurst,methodfGn,
+ sampling, subsampling,
+ Initial, Terminal, n, delta,
+ grid, random, sdelta,
+ sgrid, interpolation){
+ yuimaCogarch<-object
+ model<-yuimaCogarch at model
+ info<-model at info
+ samp <- yuimaCogarch at sampling
+ aux.Noise<-setModel(drift="1",
+ diffusion="0",
+ jump.coeff="1",
+ measure=info at measure,
+ measure.type=info at measure.type)
+
+# aux.samp<-setSampling(Initial = samp at Initial, Terminal = samp at Terminal[1], n = samp at n[1], delta = samp at delta,
+# grid=samp at grid, random = samp at random, sdelta=samp at sdelta,
+# sgrid=samp at sgrid, interpolation=samp at interpolation )
+
+ aux.samp<-setSampling(Initial = samp at Initial, Terminal = samp at Terminal[1], n = samp at n[1])
+ auxModel<-setYuima(model=aux.Noise, sampling= aux.samp)
+
+ if(length(model at parameter@measure)==0){
+ aux.incr2<-aux.simulate(object=auxModel, nsim=nsim, seed=seed,
+ space.discretized=space.discretized, increment.W=increment.W,
+ increment.L=increment.L,
+ hurst=0.5,methodfGn=methodfGn)
+ }else{
+ aux.incr2<-aux.simulate(object=auxModel, nsim=nsim, seed=seed,
+ true.parameter = true.parameter[model at parameter@measure],
+ space.discretized=space.discretized, increment.W=increment.W,
+ increment.L=increment.L,
+ hurst=0.5,methodfGn=methodfGn)
+ }
+ increment<-diff(as.numeric(get.zoo.data(aux.incr2)[[1]]))
+ # Using the simulated increment for generating the quadratic variation
+ # As first step we compute it in a crude way. A more fine approach is based on
+ # the mpv function.
+ quadratVariation <- increment^2
+ incr.L <- t(matrix(c(increment,quadratVariation),ncol=2))
+ incr.W <- matrix(0, nrow=1,ncol=length(increment))
+ # Simulate trajectories Cogarch
+ result <- aux.simulate(object=yuimaCogarch, nsim=nsim, seed=seed, xinit=xinit,
+ true.parameter = true.parameter,
+ space.discretized = space.discretized,increment.W =incr.W,
+ increment.L=incr.L,
+ hurst=hurst,methodfGn=methodfGn,
+ sampling=sampling, subsampling=subsampling,
+ Initial=Initial, Terminal=Terminal, n=n, delta=delta,
+ grid=grid, random=random, sdelta=sdelta,
+ sgrid=sgrid, interpolation=interpolation)
+ return(result)
+}
+
Modified: pkg/yuima/R/yuima.R
===================================================================
--- pkg/yuima/R/yuima.R 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/R/yuima.R 2014-11-16 23:30:52 UTC (rev 353)
@@ -293,6 +293,7 @@
is.fracdiff <- FALSE
is.jumpdiff <- FALSE
is.carma <- FALSE
+ is.cogarch <- FALSE
is.poisson <- is.Poisson(mod)
if(length(mod at drift)>0) has.drift <- TRUE
@@ -315,9 +316,14 @@
}
if( class(mod) == "yuima.carma")
is.carma <- TRUE
-
+ if( class(mod) == "yuima.cogarch"){
+ is.cogarch <- TRUE
+ is.wienerdiff <- FALSE
+ is.fracdiff <- FALSE
+ }
+
if( is.wienerdiff | is.fracdiff | is.jumpdiff ){
- if( is.wienerdiff & ! is.carma & !is.poisson){
+ if( is.wienerdiff & ! is.carma & !is.poisson & !is.cogarch){
cat("\nDiffusion process")
if( is.fracdiff ){
if(!is.na(mod at hurst)){
@@ -331,16 +337,20 @@
}
if(is.carma)
cat(sprintf("\nCarma process p=%d, q=%d", mod at info@p, mod at info@q))
+ if(is.cogarch)
+ cat(sprintf("\nCogarch process p=%d, q=%d", mod at info@p, mod at info@q))
if(is.poisson)
cat("\nCompound Poisson process")
- if( is.jumpdiff ){
+ if( (is.jumpdiff & ! is.cogarch) ){
if( (is.wienerdiff | is.carma) & !is.poisson ){
cat(" with Levy jumps")
} else {
if(!is.poisson)
cat("Levy jump process")
}
+ }else{
+ cat(" with Levy jumps")
}
cat(sprintf("\nNumber of equations: %d", mod at equation.number))
@@ -348,6 +358,8 @@
cat(sprintf("\nNumber of Wiener noises: %d", mod at noise.number))
if(is.jumpdiff)
cat(sprintf("\nNumber of Levy noises: %d", 1))
+ if(is.cogarch)
+ cat(sprintf("\nNumber of quadratic variation: %d", 1))
if(length(mod at parameter@all)>0){
cat(sprintf("\nParametric model with %d parameters",length(mod at parameter@all)))
}
Modified: pkg/yuima/man/setCarma.Rd
===================================================================
--- pkg/yuima/man/setCarma.Rd 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/man/setCarma.Rd 2014-11-16 23:30:52 UTC (rev 353)
@@ -33,7 +33,7 @@
\usage{
setCarma(p,q,loc.par=NULL,scale.par=NULL,ar.par="a",ma.par="b",
-lin.par=NULL,Carma.var="v",Latent.var="x",XinExpr=FALSE, ...)
+lin.par=NULL,Carma.var="v",Latent.var="x",XinExpr=FALSE, Cogarch=FALSE, ...)
}
\arguments{
@@ -47,6 +47,7 @@
\item{Latent.var}{a character-string that is the label of the unobserved process. Defaults to \code{"x"}.}
\item{lin.par}{a character-string that is the label of the linear coefficients. If \code{lin.par=NULL}, the default, the 'setCarma' builds the CARMA(p, q) model defined as in Brockwell (2000). }
\item{XinExpr}{a logical variable. The default value \code{XinExpr=FALSE} implies that the starting condition for \code{Latent.var} is zero. If \code{XinExpr=TRUE}, each component of \code{Latent.var} has a parameter as a initial value.}
+ \item{Cogarch}{a logical variable. The default value \code{Cogarch=FALSE} implies that the parameters are specified according to Brockwell (2000). }
\item{...}{Arguments to be passed to 'setCarma', such as the slots of \code{\link{yuima.model-class}}
\describe{
\item{\code{measure}}{Levy measure of jump variables.}
Modified: pkg/yuima/man/setCogarch.Rd
===================================================================
--- pkg/yuima/man/setCogarch.Rd 2014-11-11 02:18:07 UTC (rev 352)
+++ pkg/yuima/man/setCogarch.Rd 2014-11-16 23:30:52 UTC (rev 353)
@@ -14,25 +14,25 @@
\code{dGt = sqrt(Vt)dZt}
-\code{Vt = mu0 + (b0 Yt(0) + ... + b(q) Yt(q))}
+\code{Vt = a0 + (a1 Yt(1) + ... + b(p) Yt(p))}
-\code{dYt(0) = Yt(1) dt}
+\code{dYt(1) = Yt(2) dt}
\code{...}
- \code{dYt(p-2) = Yt(p-1) dt}
+ \code{dYt(q-1) = Yt(q) dt}
- \code{dYt(p-1) = (-a(p) Yt(0) - ... - a(1) Yt(p-1))dt + (mu0 + (b0 Yt(0) + ... + b(q) Yt(q))d[ZtZt]_{q}}
+ \code{dYt(q) = (-b(q) Yt(1) - ... - b(1) Yt(q))dt + (a0 + (a1 Yt(1) + ... + a(p) Yt(p))d[ZtZt]_{q}}
}
\usage{
-setCogarch(p, q, ar.par = "a", ma.par = "b", loc.par = "mu", Cogarch.var = "g", V.var = "v", Latent.var = "y", jump.variable = "z", time.variable = "t", measure = NULL, measure.type = NULL, XinExpr = FALSE, startCogarch = 0, work = TRUE, ...)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/yuima -r 353
More information about the Yuima-commits
mailing list