To determine p-values and or standard errors for estimates, we combine the estimation method for COGARCH in yuima with the bootstrap methodology. Here you can find an example of the code tha we used in our work “Discrete‐Time Approximation of a Cogarch(p,q) Model and its Estimation”.

library(yuima)
Cog13<-setCogarch(p=1,q=3,
                      measure=list(intensity="1", df=list("dnorm(z, 0, 1)")),
                      measure.type = "CP",XinExpr = TRUE)


par13 <- list(a0=6.77e-04,a1=1.93,b1=10.6 ,b2=35.54,b3=37.5,
               y01=1.903289e-05,y02=0,y03=0)
 
check13 <- Diagnostic.Cogarch(Cog13, param = par13)

sampCP <- setSampling(0, Terminal = 3200, n = 64000)

set.seed(123)
simCP13 <- simulate(object = Cog13, true.parameter = par13,
                    sampling = sampCP, method = "euler")
x11()
plot(simCP13)

system.time(
  res13 <- qmle(yuima = simCP13,start = par13,
                        method="Nelder-Mead", grideq=TRUE,
                        Est.Incr ="Incr")
)
coef(res13)
summary(res13)
Diagnostic.Cogarch(res13)
}

library(doSNOW)
library(foreach)

numbCores <- 4
cl<-makeCluster(numbCores, type="SOCK")
registerDoSNOW(cl)
clusterSetupRNG(cl, seed = c(1,251,501,751))

estCoeff13 <- foreach(i=1:numbCores, .combine = rbind) %dopar%{
 
  library(yuima)

  cog13 <- setCogarch(p = 1,q = 3,
    measure = list(intensity = "1", df = list("dnorm(z, 0, 1)")),
    measure.type = "CP", XinExpr=TRUE)
 
  par13 <- list(a0=6.77e-04,a1=1.93,b1=10.6 ,b2=35.54,b3=37.5,
    y01=1.903289e-05,y02=0,y03=0)
 
 
 
  samp13 <- setSampling(0, 3200, n=64000)
 
  nRep <- 250

  RESforCor <- matrix(NA,nRep,5)
  for(j in c(1:nRep)){
    sim13 <-simulate(object = cog13, true.parameter = par13,
      sampling = samp13, method="euler")
    res13 <- NULL
    res13 <- tryCatch(qmle(yuima = sim13, start = par13, grideq=TRUE,
      method = "SANN", control=list(maxit=15000)), error = function(par13){NULL})
    if(is.null(res13)){
      RESforCor[j,] <- matrix(NaN, 1, 5)
    }else{
      namecoef <- c("a0", "a1", "b1", "b2", "b3")
      cond <- sum((coef(res13)[namecoef]-unlist(par13)[namecoef])^2)==0
      if(cond){
        res13a <- tryCatch(qmle(yuima = sim13, start = par13, grideq=TRUE,
                                method = "Nelder-Mead"), error = function(par13){NULL})
        res13<-res13a
      }
      RESforCor[j,] <- coef(res13)[c("a0", "a1", "b1", "b2", "b3")]
    }
  }
  RESforCor
}

stopCluster(cl)

colnames(estCoeff13) <- c("a0 = 6.77e-04", "a1 = 1.93", "b1 = 10.6", "b2 = 35.54", "b3 = 37.5")
write.csv(x=estCoeff13,file = "cog13_CPSANNACC.csv")
colMeans(estCoeff13)
View(estCoeff13)

sd(estCoeff13[,1])
sd(estCoeff13[,3])
sd(estCoeff13[,4])
sd(estCoeff13[,5])
sd(estCoeff13[,6])
write.csv(x = estCoeff13, file = "COG13_CP.csv")