[Yuima-commits] r425 - in pkg/yuima: R man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Mar 29 13:35:39 CEST 2016
Author: lorenzo
Date: 2016-03-29 13:35:38 +0200 (Tue, 29 Mar 2016)
New Revision: 425
Added:
pkg/yuima/man/cogarch.est.incr.rd
pkg/yuima/man/cogarch.est.rd
Removed:
pkg/yuima/man/cogarch.gmm.incr.rd
pkg/yuima/man/cogarch.gmm.rd
Modified:
pkg/yuima/R/ClassCogarch.R
pkg/yuima/R/DiagnosticCogarch.R
pkg/yuima/R/MM.COGARCH.R
pkg/yuima/R/PseudoLogLikCOGARCH.R
pkg/yuima/R/cogarchNoise.R
pkg/yuima/R/qmle.R
pkg/yuima/R/simulate.R
pkg/yuima/man/cogarchNoise.Rd
pkg/yuima/man/gmm.rd
pkg/yuima/man/qmle.Rd
Log:
Modified Methods and Class for COGARCH model
Modified: pkg/yuima/R/ClassCogarch.R
===================================================================
--- pkg/yuima/R/ClassCogarch.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/ClassCogarch.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -1,8 +1,8 @@
# In order to deal the cogarch model in the \texttt{yuima} package
-# We define a new class called \texttt{yuima.cogarch} and its structure
+# We define a new class called \texttt{yuima.cogarch} and its structure
# is similar to those used for the carma model.
# The class \texttt{yuima.cogarch} extends the \texttt{yuima.model} and has
-# an additional slot that contains informations about the model stored in an object of
+# an additional slot that contains informations about the model stored in an object of
# class \texttt{cogarch.info}.
# The class \texttt{cogarch.info} is build internally by the function \texttt{setCogarch} and it is a
# the first slot of an object of class \texttt{yuima.cogarch}.
@@ -30,34 +30,37 @@
# Class 'gmm.cogarch'
-setClass("cogarch.gmm",representation(
- model = "yuima.cogarch",
+setClass("cogarch.est",representation(
+ model = "yuima",
objFun="character"),
contains="mle"
)
-setClass("summary.cogarch.gmm",representation( objFun = "ANY",
- objFunVal = "ANY" ),
- contains="summary.mle"
+setClass("summary.cogarch.est",
+ representation(objFun = "ANY",
+ objFunVal = "ANY"),
+ contains="summary.mle"
)
-setClass("cogarch.gmm.incr",representation(Incr.Lev = "ANY",
- logL.Incr = "ANY"
+setClass("cogarch.est.incr",
+ representation(Incr.Lev = "ANY",
+ logL.Incr = "ANY"
),
-contains="cogarch.gmm"
+contains="cogarch.est"
)
-setClass("summary.cogarch.gmm.incr",representation(logL.Incr = "ANY",
- MeanI = "ANY",
- SdI = "ANY",
- logLI = "ANY",
- TypeI = "ANY",
- NumbI = "ANY",
- StatI ="ANY"),
- contains="summary.cogarch.gmm"
+setClass("summary.cogarch.est.incr",
+ representation(logL.Incr = "ANY",
+ MeanI = "ANY",
+ SdI = "ANY",
+ logLI = "ANY",
+ TypeI = "ANY",
+ NumbI = "ANY",
+ StatI ="ANY"),
+ contains="summary.cogarch.est"
)
@@ -71,7 +74,7 @@
# ),
# contains="mle"
# )
-#
+#
# setClass("cogarch.gmm",representation(
# model = "yuima.cogarch",
# objFun="character"),
@@ -79,4 +82,4 @@
# )
# setClass("gmm.cogarch",
# contains="mle"
-# )
\ No newline at end of file
+# )
Modified: pkg/yuima/R/DiagnosticCogarch.R
===================================================================
--- pkg/yuima/R/DiagnosticCogarch.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/DiagnosticCogarch.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -1,56 +1,56 @@
# We write a Diagnostic function that evaluates the following quantity
-Diagnostic.Cogarch <- function(yuima.cogarch, param = list(),
+Diagnostic.Cogarch <- function(yuima.cogarch, param = list(),
matrixS = NULL , mu = 1,
display =TRUE){
-
-
+
+
if(missing(yuima.cogarch))
yuima.stop("yuima.cogarch or yuima object is missing.")
-
- if((length(param)==0 && !is(yuima.cogarch, "cogarch.gmm"))){
+
+ if((length(param)==0 && !is(yuima.cogarch, "cogarch.est"))){
yuima.stop("missing values parameters")
}
-
-
+
+
if(is(yuima.cogarch,"yuima")){
model<-yuima.cogarch at model
}else{
if(is(yuima.cogarch,"yuima.cogarch")){
model<-yuima.cogarch
}else{
- if(is(yuima.cogarch,"cogarch.gmm")){
- model<-yuima.cogarch at model
+ if(is(yuima.cogarch,"cogarch.est")){
+ model<-yuima.cogarch at model@model
if(length(param)==0){
param<-coef(yuima.cogarch)
}
}
}
}
-
+
if(!is.COGARCH(model)){
yuima.warn("The model does not belong to the class yuima.cogarch")
}
-
+
info <- model at info
numb.ar <- info at q
ar.name <- paste(info at ar.par,c(numb.ar:1),sep="")
numb.ma <- info at p
ma.name <- paste(info at ma.par,c(1:numb.ma),sep="")
loc.par <- info at loc.par
-
+
if(is.list(param)){
param <- unlist(param)
- }
-
+ }
+
xinit.name0 <- model at parameter@xinit
idx <- na.omit(match(c(loc.par, ma.name), xinit.name0))
xinit.name <- xinit.name0[-idx]
-
+
fullcoeff <- c(ar.name, ma.name, loc.par,xinit.name)
nm<-names(param)
oo <- match(nm, fullcoeff)
- if(!is(yuima.cogarch,"cogarch.gmm")){
+ if(!is(yuima.cogarch,"cogarch.est")){
if(length(na.omit(oo))!=length(fullcoeff))
yuima.stop("some named arguments in 'param' are not arguments to the supplied yuima.cogarch model")
}
@@ -78,7 +78,7 @@
ev.dum<-matrix(0,1,info at q)
ev.dum[1,info at q] <- 1
av.dum <-matrix(0,info at q,1)
- av.dum[c(1:info at p),1]<-acoeff
+ av.dum[c(1:info at p),1]<-acoeff
if(display==TRUE){
cat(paste0("\n COGARCH(",info at p,info at q,") model \n"))
}
@@ -113,7 +113,7 @@
if(is.complex(lambda.eig) && all(Im(lambda.eig)==0) && all(Re(lambda.eig) < 0)){
massage <- "\n the Variance is a positive process \n"
res.pos <- TRUE
- }
+ }
}
if((info at q==2 && info at p==2 )){
if(is.numeric(lambda.eig)|| (is.complex(lambda.eig) && all(Im(lambda.eig)==0))){
@@ -140,7 +140,7 @@
if(display==TRUE){
cat(massage)
}
-
+
res1<-StationaryMoments(cost,b,acoeff,mu)
res <- list(meanVarianceProc=res1$ExpVar,
meanStateVariable=res1$ExpStatVar,
@@ -157,19 +157,19 @@
StationaryMoments<-function(cost,b,acoeff,mu=1){
# We obtain stationary mean of State process
- # stationary mean of the variance
+ # stationary mean of the variance
# E\left(Y\right)=-a_{0}m_{2}\left(A+m_{2}ea'\right)^{-1}e
- q<-length(b)
+ q<-length(b)
a <- e <- matrix(0,nrow=q,ncol=1)
e[q,1] <- 1
p<-length(acoeff)
a[1:p,1] <- acoeff
B_tilde <- MatrixA(b[c(q:1)])+mu*e%*%t(a)
-
+
if(q>1){
invB<-rbind(c(-B_tilde[q,-1],1)/B_tilde[q,1],cbind(diag(q-1),matrix(0,q-1,1)))
}else{invB<-1/B_tilde}
-
+
ExpStatVar <- -cost*mu*invB%*%e
ExpVar <- cost+t(a)%*%ExpStatVar
res <- list(ExpVar=ExpVar, ExpStatVar=ExpStatVar)
Modified: pkg/yuima/R/MM.COGARCH.R
===================================================================
--- pkg/yuima/R/MM.COGARCH.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/MM.COGARCH.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -546,31 +546,151 @@
# Build an object of class mle
+# if(Est.Incr=="NoIncr"){
+# res<-new("cogarch.est", call = call, coef = coef, fullcoef = unlist(coef),
+# vcov = vcov, min = min, details = list(),
+# method = character(),
+# model = setYuima(model=model),
+# objFun = objFun
+# )
+# }
+# if(Est.Incr=="Incr"||Est.Incr=="IncrPar"){
+# L.Incr<-cogarchNoise(yuima.cogarch=model, data=observ,
+# param=as.list(coef), mu=1)
+# ttt<-observ at zoo.data[[1]]
+# tt<-index(ttt)
+# L.Incr_Fin <- zoo(L.Incr$incr.L,tt[(1+length(tt)-length(L.Incr$incr.L)):length(tt)])
+# }
+# if(Est.Incr=="Incr"){
+# # Build an object of class cogarch.gmm.incr
+# res<-new("cogarch.est.incr", call = call, coef = coef, fullcoef = unlist(coef),
+# vcov = vcov, min = min, details = list(),
+# method = character(),
+# Incr.Lev = L.Incr_Fin,
+# model = setYuima(data = L.Incr$Cogarch,
+# model=model), nobs=as.integer(length(L.Incr)+1),
+# logL.Incr = numeric(),
+# objFun= objFun
+# )
+# }
+# if(Est.Incr=="IncrPar"){
+# #estimationLevy
+#
+# fixedCon <- constdum(fixed, meas.par)
+# lowerCon <- constdum(lower, meas.par)
+# upperCon <- constdum(upper, meas.par)
+# if(aggregation==TRUE){
+# if(floor(n/index(observ at zoo.data[[1]])[n])!=env$deltaData){
+# yuima.stop("the n/Terminal in sampling information is not an integer. Aggregation=FALSE is recommended")
+# }
+# inc.levy1<-diff(cumsum(c(0,L.Incr$incr.L))[seq(from=1,
+# to=yuima at sampling@n[1],
+# by=env$deltaData
+# )])
+# }else{
+# inc.levy1 <- L.Incr$incr.L
+# }
+#
+# result.Lev <- gmm.Est.Lev(Increment.lev=c(0,inc.levy1),
+# param0=start[meas.par],
+# fixed = fixedCon[meas.par],
+# lower=lowerCon[meas.par],
+# upper=upperCon[meas.par],
+# measure=model at measure,
+# measure.type=model at measure.type,
+# aggregation=aggregation,
+# dt=1/env$deltaData
+# )
+#
+# if(is.null(result.Lev)){
+# res<-new("cogarch.est.incr", call = call, coef = coef, fullcoef = unlist(coef),
+# vcov = vcov, min = min, details = list(),
+# method = character(),
+# Incr.Lev=L.Incr_Fin,
+# model = setYuima(data = L.Incr$Cogarch,
+# model=model),
+# nobs=as.integer(length(L.Incr)+1),
+# logL.Incr = numeric(),
+# objFun= objFun
+# )
+#
+# }
+# else{
+# Inc.Parm<-result.Lev$estLevpar
+# IncVCOV<-result.Lev$covLev
+# if(length(meas.par)==length(Inc.Parm)){
+# names(Inc.Parm)<-meas.par
+# rownames(IncVCOV)<-as.character(meas.par)
+# colnames(IncVCOV)<-as.character(meas.par)
+# }
+# name.parm.cog<-names(coef)
+# coef<-c(coef,Inc.Parm)
+#
+# names.par<-names(coef)
+# cov<-NULL
+# cov<-matrix(NA,length(names.par),length(names.par))
+# rownames(cov)<-names.par
+# colnames(cov)<-names.par
+#
+# cov[unique(name.parm.cog),unique(name.parm.cog)]<-vcov
+# cov[names(Inc.Parm),names(Inc.Parm)]<-IncVCOV
+# cov<-cov
+#
+# res<-new("cogarch.est.incr", call = call, coef = coef, fullcoef = unlist(coef),
+# vcov = cov, min = min, details = list(),
+# method = character(),
+# Incr.Lev=L.Incr_Fin,
+# model = setYuima(data = L.Incr$Cogarch,
+# model=model), nobs=as.integer(length(L.Incr)+1),
+# logL.Incr = tryCatch(-result.Lev$value,error=function(theta){NULL}),
+# objFun= objFun
+# )
+#
+# }
+#
+#
+# }
+
+ res <- ExtraNoiseFromEst(Est.Incr,
+ call, coef, vcov, min, details = list(),
+ method = character(), model, objFun, observ,
+ fixed, meas.par, lower, upper, env, yuima, start, aggregation)
+ return(res)
+
+}
+
+ExtraNoiseFromEst <- function(Est.Incr,
+ call, coef, vcov, min, details = list(),
+ method = character(), model, objFun, observ,
+ fixed, meas.par, lower, upper, env, yuima, start, aggregation){
+ n <- length(index(observ at original.data))
if(Est.Incr=="NoIncr"){
- res<-new("cogarch.gmm", call = call, coef = coef, fullcoef = unlist(coef),
- vcov = vcov, min = min, details = list(),
- method = character(),
- model = model,
- objFun = objFun
- )
+ res<-new("cogarch.est", call = call, coef = coef,
+ fullcoef = unlist(coef),
+ vcov = vcov, min = min, details = details,
+ method = method,
+ model = setYuima(model=model),
+ objFun = objFun
+ )
}
if(Est.Incr=="Incr"||Est.Incr=="IncrPar"){
L.Incr<-cogarchNoise(yuima.cogarch=model, data=observ,
- param=as.list(coef), mu=1)
+ param=as.list(coef), mu=1)
ttt<-observ at zoo.data[[1]]
tt<-index(ttt)
- L.Incr_Fin <- zoo(L.Incr,tt[(1+length(tt)-length(L.Incr)):length(tt)])
+ L.Incr_Fin <- zoo(L.Incr$incr.L,tt[(1+length(tt)-length(L.Incr$incr.L)):length(tt)])
}
if(Est.Incr=="Incr"){
- # Build an object of class cogarch.gmm.incr
- res<-new("cogarch.gmm.incr", call = call, coef = coef, fullcoef = unlist(coef),
- vcov = vcov, min = min, details = list(),
- method = character(),
- Incr.Lev = L.Incr_Fin,
- model = model, nobs=as.integer(length(L.Incr)+1),
- logL.Incr = numeric(),
- objFun= objFun
- )
+ # Build an object of class cogarch.gmm.incr
+ res<-new("cogarch.est.incr", call = call, coef = coef, fullcoef = unlist(coef),
+ vcov = vcov, min = min, details = list(),
+ method = character(),
+ Incr.Lev = L.Incr_Fin,
+ model = setYuima(data = L.Incr$Cogarch,
+ model=model), nobs=as.integer(length(L.Incr)+1),
+ logL.Incr = numeric(),
+ objFun= objFun
+ )
}
if(Est.Incr=="IncrPar"){
#estimationLevy
@@ -580,36 +700,38 @@
upperCon <- constdum(upper, meas.par)
if(aggregation==TRUE){
if(floor(n/index(observ at zoo.data[[1]])[n])!=env$deltaData){
- yuima.stop("the n/Terminal in sampling information is not an integer. Aggregation=FALSE is recommended")
+ yuima.stop("the n/Terminal in sampling information is not an integer. aggregation=FALSE is recommended")
}
- inc.levy1<-diff(cumsum(c(0,L.Incr))[seq(from=1,
- to=yuima at sampling@n[1],
- by=env$deltaData
+ inc.levy1<-diff(cumsum(c(0,L.Incr$incr.L))[seq(from=1,
+ to=yuima at sampling@n[1],
+ by=env$deltaData
)])
}else{
- inc.levy1 <- L.Incr
+ inc.levy1 <- L.Incr$incr.L
}
result.Lev <- gmm.Est.Lev(Increment.lev=c(0,inc.levy1),
- param0=start[meas.par],
- fixed = fixedCon[meas.par],
- lower=lowerCon[meas.par],
- upper=upperCon[meas.par],
- measure=model at measure,
- measure.type=model at measure.type,
- aggregation=aggregation,
- dt=1/env$deltaData
- )
+ param0=start[meas.par],
+ fixed = fixedCon[meas.par],
+ lower=lowerCon[meas.par],
+ upper=upperCon[meas.par],
+ measure=model at measure,
+ measure.type=model at measure.type,
+ aggregation=aggregation,
+ dt=1/env$deltaData
+ )
if(is.null(result.Lev)){
- res<-new("cogarch.gmm.incr", call = call, coef = coef, fullcoef = unlist(coef),
- vcov = vcov, min = min, details = list(),
- method = character(),
- Incr.Lev=L.Incr_Fin,
- model = model, nobs=as.integer(length(L.Incr)+1),
- logL.Incr = numeric(),
- objFun= objFun
- )
+ res<-new("cogarch.est.incr", call = call, coef = coef, fullcoef = unlist(coef),
+ vcov = vcov, min = min, details = list(),
+ method = character(),
+ Incr.Lev=L.Incr_Fin,
+ model = setYuima(data = L.Incr$Cogarch,
+ model=model),
+ nobs=as.integer(length(L.Incr)+1),
+ logL.Incr = numeric(),
+ objFun= objFun
+ )
}
else{
@@ -633,21 +755,21 @@
cov[names(Inc.Parm),names(Inc.Parm)]<-IncVCOV
cov<-cov
- res<-new("cogarch.gmm.incr", call = call, coef = coef, fullcoef = unlist(coef),
+ res<-new("cogarch.est.incr", call = call, coef = coef, fullcoef = unlist(coef),
vcov = cov, min = min, details = list(),
method = character(),
Incr.Lev=L.Incr_Fin,
- model = model, nobs=as.integer(length(L.Incr)+1),
+ model = setYuima(data = L.Incr$Cogarch,
+ model=model), nobs=as.integer(length(L.Incr)+1),
logL.Incr = tryCatch(-result.Lev$value,error=function(theta){NULL}),
objFun= objFun
- )
+ )
}
- }
- return(res)
-
+ }
+ return(res)
}
constdum<-function(fixed, meas.par){
@@ -977,7 +1099,7 @@
}
-setMethod("plot",signature(x="cogarch.gmm.incr"),
+setMethod("plot",signature(x="cogarch.est.incr"),
function(x, type="l" ,...){
Time<-index(x at Incr.Lev)
Incr.L<-coredata(x at Incr.Lev)
@@ -993,14 +1115,14 @@
)
-setMethod("summary", "cogarch.gmm",
+setMethod("summary", "cogarch.est",
function (object, ...)
{
cmat <- cbind(Estimate = object at coef, `Std. Error` = sqrt(diag(object at vcov)))
obj<- object at min
labFun <- object at objFun
- tmp <- new("summary.cogarch.gmm", call = object at call, coef = cmat,
+ tmp <- new("summary.cogarch.est", call = object at call, coef = cmat,
m2logL = 0,
#model = object at model,
objFun = labFun,
@@ -1017,7 +1139,7 @@
#
# }
-setMethod("show", "summary.cogarch.gmm",
+setMethod("show", "summary.cogarch.est",
function (object)
{
@@ -1041,7 +1163,7 @@
)
-setMethod("summary", "cogarch.gmm.incr",
+setMethod("summary", "cogarch.est.incr",
function (object, ...)
{
cmat <- cbind(Estimate = object at coef, `Std. Error` = sqrt(diag(object at vcov)))
@@ -1051,14 +1173,14 @@
obj<- object at min
labFun <- object at objFun
- tmp <- new("summary.cogarch.gmm.incr", call = object at call, coef = cmat,
+ tmp <- new("summary.cogarch.est.incr", call = object at call, coef = cmat,
m2logL = m2logL,
objFun = labFun,
objFunVal = obj,
MeanI = mean(data),
SdI = sd(data),
logLI = object at logL.Incr,
- TypeI = object at model@measure.type,
+ TypeI = object at model@model at measure.type,
NumbI = length(data),
StatI =summary(data)
)
@@ -1066,11 +1188,16 @@
}
)
#
-setMethod("show", "summary.cogarch.gmm.incr",
+setMethod("show", "summary.cogarch.est.incr",
function (object)
{
+ if(!is.null(object at logLI)){
- cat("Two Stages GMM estimation \n\nCall:\n")
+ cat("Two Stages PSEUDO-LogLik estimation \n\nCall:\n")
+ }else{
+
+ cat("Two Stages GMM estimation \n\nCall:\n")
+ }
print(object at call)
cat("\nCoefficients:\n")
print(coef(object))
Modified: pkg/yuima/R/PseudoLogLikCOGARCH.R
===================================================================
--- pkg/yuima/R/PseudoLogLikCOGARCH.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/PseudoLogLikCOGARCH.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -3,7 +3,7 @@
#
PseudoLogLik.COGARCH <- function(yuima, start, method="BFGS", fixed = list(),
- lower, upper, Est.Incr, call, grideq, ...){
+ lower, upper, Est.Incr, call, grideq, aggregation,...){
if(is(yuima,"yuima")){
model <- yuima at model
@@ -55,49 +55,49 @@
I <- diag(info at q)
assign("I",I, envir = my.env)
-
+ param1 <- param[c(my.env$loc.par, ma.names, ar.names)]
out<-NULL
if(length(lower)==0 && length(upper)>0 && length(fixed)==0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, upper=upper, env = my.env,...)
}
if(length(lower)==0 && length(upper)==0 && length(fixed)>0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, fixed=fixed, env = my.env,...)
}
if(length(lower)>0 && length(upper)==0 && length(fixed)==0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, lower=lower, env = my.env,...)
}
if(length(lower)>0 && length(upper)>0 && length(fixed)==0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, upper = upper,
lower=lower, env = my.env,...)
}
if(length(lower)==0 && length(upper)>0 && length(fixed)>0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, upper = upper,
fixed = fixed, env = my.env,...)
}
if(length(lower)>0 && length(upper)==0 && length(fixed)>0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, lower = lower,
fixed = fixed, env = my.env,...)
}
if(length(lower)>0 && length(upper)>0 && length(fixed)>0){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, lower = lower,
fixed = fixed, upper = upper,
env = my.env,...)
@@ -105,12 +105,16 @@
if(is.null(out)){
- out <- optim(par=param, fn=minusloglik.COGARCH1,
+ out <- optim(par=param1, fn=minusloglik.COGARCH1,
method = method, env = my.env,...)
}
# control= list(maxit=100),
# Write the object mle with result
+# my.env1 <- my.env
+# my.env1$grideq <- TRUE
+ resHessian<-tryCatch(optimHess(par = out$par, fn = minusloglik.COGARCH1,
+ env = my.env),error=function(theta){NULL})
bvect<-out$par[ar.names]
@@ -120,32 +124,44 @@
a0<-out$par[info at loc.par]
- if(length(meas.par)!=0){
- idx.dumm<-match(meas.par,names(out$par))
- out$par<-out$par[- idx.dumm]
- }
+# if(length(meas.par)!=0){
+# idx.dumm<-match(meas.par,names(out$par))
+# out$par<-out$par[- idx.dumm]
+# }
+#
+# idx.dumm1<-match(start.state,names(out$par))
+ coef <- out$par
+ #coef <- out$par
+ vcov<-matrix(NA, length(coef), length(coef))
- idx.dumm1<-match(start.state,names(out$par))
- coef <- out$par[-idx.dumm1]
- vcov<-matrix(NA, length(coef), length(coef))
+
names_coef<-names(coef)
- colnames(vcov) <- names_coef
- rownames(vcov) <- names_coef
+ colnames(vcov)[1:length(names_coef)] <- names_coef
+ rownames(vcov)[1:length(names_coef)] <- names_coef
+
+ if(!is.null(resHessian)){
+ vcov[c(1:length(names_coef)),c(1:length(names_coef))]<- solve(resHessian)
+ }
+
mycoef <- start
# min <- out$value
objFun <- "PseudoLogLik"
# min <- numeric()
min <- out$value
- res<-new("cogarch.gmm", call = call, coef = coef, fullcoef = unlist(coef),
- vcov = vcov, min = min, details = list(),
- method = character(),
- model = model,
- objFun = objFun
- )
+# res<-new("cogarch.est", call = call, coef = coef, fullcoef = unlist(coef),
+# vcov = vcov, min = min, details = list(),
+# method = character(),
+# model = model,
+# objFun = objFun
+# )
+ env <- list(deltaData = yuima at sampling@delta)
+ res <- ExtraNoiseFromEst(Est.Incr,
+ call, coef, vcov, min, details = list(),
+ method = character(), model, objFun, observ=yuima at data,
+ fixed, meas.par, lower, upper, env, yuima, start, aggregation)
-
return(res)
}
@@ -166,7 +182,7 @@
stateMean <- a0/(bq-a1)*as.matrix(c(1,numeric(length=(env$q-1))))
- param[env$start.state]<-stateMean
+ #param[env$start.state]<-stateMean
state <- stateMean
# state <- param[env$start.state]
DeltaG2 <- env$Obs
Modified: pkg/yuima/R/cogarchNoise.R
===================================================================
--- pkg/yuima/R/cogarchNoise.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/cogarchNoise.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -8,14 +8,14 @@
cogarchNoise<-function(yuima.cogarch, data=NULL, param, mu=1){
if(missing(yuima.cogarch))
yuima.stop("yuima.cogarch or yuima object is missing.")
-
+
if(length(param)==0)
yuima.stop("missing values parameters")
-
+
if(!is.COGARCH(yuima.cogarch)){
yuima.warn("The model does not belong to the class yuima.cogarch")
}
-
+
if(is(yuima.cogarch,"yuima")){
model<-yuima.cogarch at model
if(is.null(data)){
@@ -30,44 +30,44 @@
observ<-data
}
}
-
+
info <- model at info
numb.ar <- info at q
ar.name <- paste(info at ar.par,c(numb.ar:1),sep="")
numb.ma <- info at p
ma.name <- paste(info at ma.par,c(1:numb.ma),sep="")
loc.par <- info at loc.par
-
+
nm <- c(names(param))
param<-as.numeric(param)
names(param)<-nm
-
+
xinit.name0 <- model at parameter@xinit
idx <- na.omit(match(c(loc.par, ma.name), xinit.name0))
xinit.name <- xinit.name0[-idx]
fullcoeff <- c(ar.name, ma.name, loc.par,xinit.name)
-
-
-
+
+
+
# oo <- match(nm, fullcoeff)
-#
+#
# if(any(is.na(oo)))
# yuima.stop("some named arguments in 'param' are not arguments to the supplied yuima.cogarch model")
-
+
acoeff <- param[ma.name]
b <- param[ar.name]
cost<- param[loc.par]
Data<-as.matrix(onezoo(observ)[,1])
freq<-round(frequency(onezoo(observ)[,1]))
- res<-auxcogarch.noise(cost,b,acoeff,mu,Data,freq)
+ res<-auxcogarch.noise(cost,b,acoeff,mu,Data,freq,model at solve.variable)
return(res)
}
-auxcogarch.noise<-function(cost,b,acoeff,mu,Data,freq){
-
+auxcogarch.noise<-function(cost,b,acoeff,mu,Data,freq,lab){
+
res<-StationaryMoments(cost,b,acoeff,mu)
ExpY0<-res$ExpStatVar
-
- q<-length(b)
+
+ q<-length(b)
a <- e <- matrix(0,nrow=q,ncol=1)
e[q,1] <- 1
p<-length(acoeff)
@@ -75,14 +75,16 @@
B <- MatrixA(b[c(q:1)])
DeltaG <- c(0,diff(Data))
squaredG <- DeltaG^2
-
+
Process_Y <- ExpY0
+# Process_Y1 <- ExpY0
# Process_Y <- as.matrix(50.33)
var_V<-cost + sum(acoeff*Process_Y)
delta <- 1/freq
- for(t in c(2:(length(Data)))){
+ for(t in c(2:(length(Data)))){
# Y_t=e^{A\Delta t}Y_{t-\Delta t}+e^{A\left(\Delta t\right)}e\left(\Delta G_{t}\right)^{2}
Process_Y <- cbind(Process_Y, (expm(B*delta)%*%(Process_Y[,t-1]+e*squaredG[t])))
+# Process_Y1 <- merge(Process_Y1, (expm(B*delta)%*%(Process_Y1[,t-1]+e*squaredG[t])))
#Process_Y <- cbind(Process_Y, (Process_Y[,t-1]+delta*B%*%Process_Y[,t-1]+e*squaredG[t]))
# sim[t,3:ncolsim]<-sim[t-1,3:ncolsim]+(AMatrix*Delta)%*%sim[t-1,3:ncolsim]+evect*sim[t-1,2]*incr.L[2,t-1]
# sim[t,2]<-value.a0+tavect%*%sim[t-1,3:ncolsim]
@@ -90,7 +92,13 @@
var_V[t] <- cost + t(a)%*%Process_Y[,t-1]
}
#\Delta L_{t}=\frac{\Delta G_{t}}{\sqrt{V_{t}}}.
-
+
incr.L<- DeltaG/sqrt(var_V)
- return(incr.L)
-}
\ No newline at end of file
+ CogDum <- cbind(var_V,t(Process_Y))
+ DumRet <- as.numeric(Data[1]+cumsum(DeltaG))
+ Cog <- cbind(DumRet,CogDum)
+ colnames(Cog) <- lab
+ Cogarch <- setData(Cog, delta = delta)
+ res<-list(incr.L=incr.L, Cogarch=Cogarch)
+ return(res)
+}
Modified: pkg/yuima/R/qmle.R
===================================================================
--- pkg/yuima/R/qmle.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/qmle.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -122,7 +122,7 @@
}
qmle <- function(yuima, start, method="BFGS", fixed = list(), print=FALSE,
- lower, upper, joint=FALSE, Est.Incr="Carma.IncPar",aggregation=TRUE, threshold=NULL, ...){
+ lower, upper, joint=FALSE, Est.Incr="NoIncr",aggregation=TRUE, threshold=NULL, ...){
if(is(yuima at model, "yuima.carma")){
NoNeg.Noise<-FALSE
cat("\nStarting qmle for carma ... \n")
@@ -144,10 +144,10 @@
res <- NULL
if("grideq" %in% names(as.list(call)[-(1:2)])){
res <- PseudoLogLik.COGARCH(yuima, start, method=method, fixed = list(),
- lower, upper, Est.Incr, call, ...)
+ lower, upper, Est.Incr, call, aggregation = aggregation, ...)
}else{
res <- PseudoLogLik.COGARCH(yuima, start, method=method, fixed = list(),
- lower, upper, Est.Incr, call, grideq = FALSE,...)
+ lower, upper, Est.Incr, call, grideq = FALSE, aggregation = aggregation,...)
}
return(res)
@@ -1007,12 +1007,12 @@
method = method, nobs=yuima.nobs, model=yuima at model)
}
} else {
- if( Est.Incr=="Carma.IncPar" || Est.Incr=="Carma.Inc" ){
+ if( Est.Incr=="IncrPar" || Est.Incr=="Incr" ){
final_res<-new("yuima.carma.qmle", call = call, coef = coef, fullcoef = unlist(mycoef),
vcov = vcov, min = min, details = oout, minuslogl = minusquasilogl,
method = method, nobs=yuima.nobs, logL.Incr = NULL)
}else{
- if(Est.Incr=="Carma.Par"){
+ if(Est.Incr=="NoIncr"){
final_res<-new("mle", call = call, coef = coef, fullcoef = unlist(mycoef),
vcov = vcov, min = min, details = oout, minuslogl = minusquasilogl,
method = method, nobs=yuima.nobs)
Modified: pkg/yuima/R/simulate.R
===================================================================
--- pkg/yuima/R/simulate.R 2016-03-26 23:05:54 UTC (rev 424)
+++ pkg/yuima/R/simulate.R 2016-03-29 11:35:38 UTC (rev 425)
@@ -815,7 +815,7 @@
# Simulate method for an object of class cogarch.gmm.incr
-setMethod("simulate","cogarch.gmm.incr",
+setMethod("simulate","cogarch.est.incr",
function(object, nsim=1, seed=NULL, xinit, ...){
out <-aux.simulategmm(object=object, nsim=nsim, seed=seed, xinit=xinit, ...)
Added: pkg/yuima/man/cogarch.est.incr.rd
===================================================================
--- pkg/yuima/man/cogarch.est.incr.rd (rev 0)
+++ pkg/yuima/man/cogarch.est.incr.rd 2016-03-29 11:35:38 UTC (rev 425)
@@ -0,0 +1,37 @@
+\name{cogarch.est.incr-class}
+\docType{class}
+\alias{cogarch.est.incr-class}
+\alias{plot,cogarch.est.incr,ANY-method}
+\alias{est.cogarch.incr-class}
+\alias{cogarch.est.incr-class}
+\alias{simulate,cogarch.est.incr-method}
+%%\alias{setSampling,yuima.carma-method}
+
+\title{Class for Estimation of COGARCH(p,q) model with underlying increments}
+\description{
+ The \code{cogarch.est.incr} class is a class of the \pkg{yuima} package that extends the \code{\link{cogarch.est-class}} and is filled by the function \code{\link{gmm}} or \code{\link{qmle}}.
+}
+\section{Slots}{
+ \describe{
+ \item{\code{Incr.Lev}:}{is an object of class \code{\link{zoo}} that contains the estimated increments of the noise obtained using \code{\link{cogarchNoise}}.}
+ \item{\code{model}:}{is an object of of \code{\link{yuima-class}}.}
+ \item{\code{logL.Incr}:}{is an object of class \code{numeric} that contains the value of the log-likelihood for estimated Levy increments.}
+ \item{\code{objFun}:}{is an object of class \code{character} that indicates the objective function used in the minimization problem. See the documentation of the function \code{\link{gmm}} or \code{\link{qmle}} for more details.}
+ \item{\code{call}:}{is an object of class \code{language}. }
+ \item{\code{coef}:}{is an object of class \code{numeric} that contains estimated parameters.}
+ \item{\code{fullcoef}:}{is an object of class \code{numeric} that contains estimated and fixed parameters.}
+ \item{\code{vcov}:}{is an object of class \code{matrix}.}
+ \item{\code{min}:}{is an object of class \code{numeric}.}
+ \item{\code{minuslogl}:}{is an object of class \code{function}.}
+ \item{\code{method}:}{is an object of class \code{character}.}
+ }
+}
+\section{Methods}{
+ \describe{
+ \item{simulate}{simulation method. For more information see \code{\link{simulate}}.}
+ \item{plot}{Plot method for estimated increment of the noise.}
+ \item{Methods mle}{All methods for \code{mle-class} are available.}
+ }
+}
+\author{The YUIMA Project Team}
+\keyword{classes}
Added: pkg/yuima/man/cogarch.est.rd
===================================================================
--- pkg/yuima/man/cogarch.est.rd (rev 0)
+++ pkg/yuima/man/cogarch.est.rd 2016-03-29 11:35:38 UTC (rev 425)
@@ -0,0 +1,30 @@
+\name{cogarch.est.-class}
+\docType{class}
+\alias{cogarch.est-class}
+\alias{plot,cogarch.est.,ANY-method}
+%%\alias{setSampling,yuima.carma-method}
+
+\title{Class for Generalized Method of Moments Estimation for COGARCH(p,q) model}
+\description{
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/yuima -r 425
More information about the Yuima-commits
mailing list