Nothing
R/copula-distributions.R
In rgarch: Flexible GARCH modelling in R
#################################################################################
##
## R package rgarch by Alexios Ghalanos Copyright (C) 2008, 2009, 2010, 2011
## This file is part of the R package rgarch.
##
## The R package rgarch is free software: you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
##
## The R package rgarch is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
#################################################################################
# Implements copula methods for p, d, q, r for multivariate Normal amd Student.
# Future expansions should likely include the multivariate skew Normal and skew Student
# of Azzalini.
# General Setup:
# dcopula:
#-----------------------------------------------------------------------------------------------------------
#
# c(u_1, ... , u_n) = multivariate_pdf( margin_quantile_1(u_1), ... , margin_quantile_n(u_n) )
# ----------------------------------------------------------------------------------------------
# SumProduct_{1..n}[ margin_pdf_i( margin_quantile_i(u_i) ) ]
#
# where u_i is the uniform margin obtained by applying the univariate_cdf to the fitted data (each i can have
# a seperate distribution fitted)
#-----------------------------------------------------------------------------------------------------------
# pcopula:
#-----------------------------------------------------------------------------------------------------------
# C(u_1, ..., u_n) = multivariate_cdf( margin_quantile_1(u_1), ... , margin_quantile_n(u_n) )
#
# where u_i is the uniform margin obtained by applying the univariate_cdf to the fitted data (each i can have
# a seperate distribution fitted)
#-----------------------------------------------------------------------------------------------------------
# conversely, a multivariate density function multivariate_pdf() can be decomposed as:
# dinvcopula:
#-----------------------------------------------------------------------------------------------------------
# multivariate_pdf(x_1, ..., x_n) = c(margin_cdf_1(x_1), ... , margin_cdf_n(x_n)) * SumProduct_{1..n}[ univariate_pdf(x_i) ]
#
#-----------------------------------------------------------------------------------------------------------
# pinvcopula:
#-----------------------------------------------------------------------------------------------------------
# multivariate_cdf(x_1, ..., x_n) = c(margin_cdf_1(x_1), ... , margin_cdf_n(x_n))
#
#-----------------------------------------------------------------------------------------------------------
# Density of full Copula GARCH Model
dcgarch = function(mfit, cfit)
{
# The full density of the model is given by:
# log(copula density) + log(marginal GARCH densities)
# where copula density = log( multivariate_f(q_1*(u_1),...,q_n(u_n)) ) ...
# - (log(univariate_f(q_1*(u_1)) + log(univariate_f(q_n*(u_n)))
garchllh = sum(sapply(mfit@fit, FUN = function(x) likelihood(x)))
copllh = cfit$LLH
llh = garchllh + copllh
return(llh)
}
dcopula.gauss = function(U, Sigma, logvalue = FALSE){
m = dim(Sigma)[2]
# U = uniform(0,1)
Z = apply(U, 2, FUN = function(x) qnorm(x))
# Z = transormed to standard normal variates
ans = dmvnorm(Z, mean = rep(0, m), sigma = Sigma, log = TRUE) - apply(dnorm(Z, log = TRUE), 1, "sum")
# .Call("dcopulaGauss", Z = Z, m = matrix(0, nrow = 1, ncol = m), sigma = Sigma, dnormZ = dnorm(Z, log = TRUE),
# PACKAGE = "rgarch")
if( !logvalue ) ans = exp( ans )
return( ans )
}
pcopula.gauss = function(U, Sigma, ...)
{
m = dim(Sigma)[2]
U = matrix(U, ncol = m)
Z = apply(U, 2, FUN = function(x) qnorm(x))
mu = rep(0, m)
ans = apply(Z, 2, FUN = function(x) pmvnorm(lower = rep(-Inf, m), upper = x, mean = mu, sigma = Sigma, ...))
return( ans )
}
rcopula.gauss = function(n, U, Sigma, ...)
{
m = dim(Sigma)[2]
mu = rep(0, m)
ans = pnorm(rmvnorm(n, mean = mu, sigma = Sigma, method = "eigen"))
return( ans )
}
dcopula.student = function(U, Corr, df, logvalue = FALSE)
{
m = dim(Corr)[2]
Z = apply(U, 2, FUN = function(x) qt(p = x , df = df))
mu = rep(0, m)
ans = dmvt(Z, delta = mu, sigma = Corr, df = df, log = TRUE) - apply(Z, 1, FUN = function(x) sum(dt(x, df = df, log = TRUE)))
# .Call("dcopulaStudent", Z = Z, m = matrix(0, nrow = 1, ncol = m), sigma = Corr, df = df, dtZ = dt(Z, df = df, log = TRUE),
# PACKAGE = "rgarch")
if( !logvalue ) ans = exp(ans)
return( ans )
}
pcopula.student = function(U, Corr, df, ...)
{
m = dim(Corr)[2]
U = matrix(U, ncol = m)
Z = apply(U, 2, FUN = function(x) qt(x, df = df))
mu = rep(0, m)
ans = apply(Z, 1, FUN = function(x) pmvt(lower = rep(-Inf, m), upper = x, delta = mu, corr = Corr, df = df, ...))
return( ans )
}
rcopula.student = function(n, U, Corr, df)
{
m = dim(Corr)[2]
mu = rep(0, m)
ans = pt(rmvt(n, delta = mu, sigma = Corr, df = df), df = df)
return ( ans )
}
#####################################################################################
# Copula Simulation
#------------------------------------------------------------------------------------
.sample.copula = function(dist, object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
dist = tolower(dist)
set.seed(rseed)
ans = switch(dist,
mvt0 = .rtcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control),
mvt1 = .rtvtcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control),
mvnorm0 = .rnormcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control),
mvnorm1 = .rtvnormcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control))
return(ans)
}
.rtvtcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
dccOrder = object@mspec$dccOrder
cf = coef(object, type = "st")
shape = cf["shape"]
mo = max( dccOrder )
dcca = cf[1:dccOrder[1]]
dccb = cf[(dccOrder[1]+1):sum(dccOrder[1:2])]
preQ = Rbar
m = dim(Rbar)[1]
z = array(NA, dim = c(n.sim + n.start + mo, m, m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
z[,,i] = rmvt(n = (n.sim + n.start + mo), delta = rep(0, m), sigma = diag(m), df = shape)
}
xseed = rseed+1
if( parallel ){
os = .Platform$OS.type
if(is.null(parallel.control$pkg)){
if( os == "windows" ) parallel.control$pkg = "snowfall" else parallel.control$pkg = "multicore"
if( is.null(parallel.control$cores) ) parallel.control$cores = 2
} else{
mtype = match(tolower(parallel.control$pkg[1]), c("multicore", "snowfall"))
if(is.na(mtype)) stop("\nParallel Package type not recognized in parallel.control\n")
parallel.control$pkg = tolower(parallel.control$pkg[1])
if( os == "windows" && parallel.control$pkg == "multicore" ) stop("\nmulticore not supported on windows O/S\n")
if( is.null(parallel.control$cores) ) parallel.control$cores = 2 else parallel.control$cores = as.integer(parallel.control$cores[1])
}
if( parallel.control$pkg == "multicore" ){
mtmp = mclapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]),
mc.cores = parallel.control$cores)
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pt(sapply(mtmp, FUN = function(x) x$stdresid[,i]), shape), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
} else{
sfInit(parallel = TRUE, cpus = parallel.control$cores)
sfExport("z", "preQ", "Rbar", "dcca", "dccb", "mo", "n.sim", "n.start", "m", "xseed", local = TRUE)
mtmp = sfLapply(as.list(1:m.sim), fun = function(j)
rgarch:::.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
sfStop()
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pt(sapply(mtmp, FUN = function(x) x$stdresid[,i]), shape), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
}
} else{
mtmp = lapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = pt(matrix(sapply(mtmp, FUN = function(x) x$stdresid[,i]), ncol = m.sim), shape)
for(i in 1:m.sim) Usim[,,i] = matrix(sapply(Ures, FUN = function(x) x[-(1:mo),i]), ncol = m)
}
return(list(Usim = Usim, simR = simR))
}
.rtcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
nsim = n.sim + n.start
m = dim(Rbar)[1]
shape = object@mfit$stpars["shape"]
sim = array(data = NA, dim = c(nsim , m , m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
tmp = rmvt(n = nsim, delta = rep(0, m), sigma = Rbar, df = shape)
sim[,,i] = matrix(pt(tmp, df = shape), nrow = nsim, ncol = m)
}
return( sim )
}
.rtvnormcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
dccOrder = object@mspec$dccOrder
cf = coef(object, type = "st")
mo = max( dccOrder )
dcca = cf[1:dccOrder[1]]
dccb = cf[(dccOrder[1]+1):sum(dccOrder[1:2])]
preQ = Rbar
m = dim(Rbar)[1]
z = array(NA, dim = c(n.sim + n.start + mo, m, m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
z[,,i] = rmvnorm(n = (n.sim + n.start + mo), mean = rep(0, m), sigma = diag(m))
}
xseed = rseed+1
if( parallel ){
os = .Platform$OS.type
if(is.null(parallel.control$pkg)){
if( os == "windows" ) parallel.control$pkg = "snowfall" else parallel.control$pkg = "multicore"
if( is.null(parallel.control$cores) ) parallel.control$cores = 2
} else{
mtype = match(tolower(parallel.control$pkg[1]), c("multicore", "snowfall"))
if(is.na(mtype)) stop("\nParallel Package type not recognized in parallel.control\n")
parallel.control$pkg = tolower(parallel.control$pkg[1])
if( os == "windows" && parallel.control$pkg == "multicore" ) stop("\nmulticore not supported on windows O/S\n")
if( is.null(parallel.control$cores) ) parallel.control$cores = 2 else parallel.control$cores = as.integer(parallel.control$cores[1])
}
if( parallel.control$pkg == "multicore" ){
mtmp = mclapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]),
mc.cores = parallel.control$cores)
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pnorm(sapply(mtmp, FUN = function(x) x$stdresid[,i])), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
} else{
sfInit(parallel = TRUE, cpus = parallel.control$cores)
sfExport("z", "preQ", "Rbar", "dcca", "dccb", "mo", "n.sim", "n.start", "m", "xseed", local = TRUE)
mtmp = sfLapply(as.list(1:m.sim), fun = function(j)
rgarch:::.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
sfStop()
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pnorm(sapply(mtmp, FUN = function(x) x$stdresid[,i])), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
}
} else{
mtmp = lapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = pnorm(matrix(sapply(mtmp, FUN = function(x) x$stdresid[,i]), ncol = m.sim))
for(i in 1:m.sim) Usim[,,i] = matrix(sapply(Ures, FUN = function(x) x[-(1:mo),i]), ncol = m)
}
return(list(Usim = Usim, simR = simR))
}
.rnormcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
nsim = n.sim + n.start
m = dim(Rbar)[1]
sim = array(data = NA, dim = c(nsim , m , m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
tmp = rmvnorm(n = nsim, mean = rep(0, m), sigma = Rbar)
sim[,,i] = matrix(pnorm(tmp), nrow = nsim, ncol = m)
}
return( sim )
}
#------------------------------------------------------------------------------------
Try the rgarch package in your browser
Any scripts or data that you put into this service are public.
rgarch documentation built on May 2, 2019, 5:22 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
#################################################################################
##
## R package rgarch by Alexios Ghalanos Copyright (C) 2008, 2009, 2010, 2011
## This file is part of the R package rgarch.
##
## The R package rgarch is free software: you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
##
## The R package rgarch is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
## GNU General Public License for more details.
##
#################################################################################
# Implements copula methods for p, d, q, r for multivariate Normal amd Student.
# Future expansions should likely include the multivariate skew Normal and skew Student
# of Azzalini.
# General Setup:
# dcopula:
#-----------------------------------------------------------------------------------------------------------
#
# c(u_1, ... , u_n) = multivariate_pdf( margin_quantile_1(u_1), ... , margin_quantile_n(u_n) )
# ----------------------------------------------------------------------------------------------
# SumProduct_{1..n}[ margin_pdf_i( margin_quantile_i(u_i) ) ]
#
# where u_i is the uniform margin obtained by applying the univariate_cdf to the fitted data (each i can have
# a seperate distribution fitted)
#-----------------------------------------------------------------------------------------------------------
# pcopula:
#-----------------------------------------------------------------------------------------------------------
# C(u_1, ..., u_n) = multivariate_cdf( margin_quantile_1(u_1), ... , margin_quantile_n(u_n) )
#
# where u_i is the uniform margin obtained by applying the univariate_cdf to the fitted data (each i can have
# a seperate distribution fitted)
#-----------------------------------------------------------------------------------------------------------
# conversely, a multivariate density function multivariate_pdf() can be decomposed as:
# dinvcopula:
#-----------------------------------------------------------------------------------------------------------
# multivariate_pdf(x_1, ..., x_n) = c(margin_cdf_1(x_1), ... , margin_cdf_n(x_n)) * SumProduct_{1..n}[ univariate_pdf(x_i) ]
#
#-----------------------------------------------------------------------------------------------------------
# pinvcopula:
#-----------------------------------------------------------------------------------------------------------
# multivariate_cdf(x_1, ..., x_n) = c(margin_cdf_1(x_1), ... , margin_cdf_n(x_n))
#
#-----------------------------------------------------------------------------------------------------------
# Density of full Copula GARCH Model
dcgarch = function(mfit, cfit)
{
# The full density of the model is given by:
# log(copula density) + log(marginal GARCH densities)
# where copula density = log( multivariate_f(q_1*(u_1),...,q_n(u_n)) ) ...
# - (log(univariate_f(q_1*(u_1)) + log(univariate_f(q_n*(u_n)))
garchllh = sum(sapply(mfit@fit, FUN = function(x) likelihood(x)))
copllh = cfit$LLH
llh = garchllh + copllh
return(llh)
}
dcopula.gauss = function(U, Sigma, logvalue = FALSE){
m = dim(Sigma)[2]
# U = uniform(0,1)
Z = apply(U, 2, FUN = function(x) qnorm(x))
# Z = transormed to standard normal variates
ans = dmvnorm(Z, mean = rep(0, m), sigma = Sigma, log = TRUE) - apply(dnorm(Z, log = TRUE), 1, "sum")
# .Call("dcopulaGauss", Z = Z, m = matrix(0, nrow = 1, ncol = m), sigma = Sigma, dnormZ = dnorm(Z, log = TRUE),
# PACKAGE = "rgarch")
if( !logvalue ) ans = exp( ans )
return( ans )
}
pcopula.gauss = function(U, Sigma, ...)
{
m = dim(Sigma)[2]
U = matrix(U, ncol = m)
Z = apply(U, 2, FUN = function(x) qnorm(x))
mu = rep(0, m)
ans = apply(Z, 2, FUN = function(x) pmvnorm(lower = rep(-Inf, m), upper = x, mean = mu, sigma = Sigma, ...))
return( ans )
}
rcopula.gauss = function(n, U, Sigma, ...)
{
m = dim(Sigma)[2]
mu = rep(0, m)
ans = pnorm(rmvnorm(n, mean = mu, sigma = Sigma, method = "eigen"))
return( ans )
}
dcopula.student = function(U, Corr, df, logvalue = FALSE)
{
m = dim(Corr)[2]
Z = apply(U, 2, FUN = function(x) qt(p = x , df = df))
mu = rep(0, m)
ans = dmvt(Z, delta = mu, sigma = Corr, df = df, log = TRUE) - apply(Z, 1, FUN = function(x) sum(dt(x, df = df, log = TRUE)))
# .Call("dcopulaStudent", Z = Z, m = matrix(0, nrow = 1, ncol = m), sigma = Corr, df = df, dtZ = dt(Z, df = df, log = TRUE),
# PACKAGE = "rgarch")
if( !logvalue ) ans = exp(ans)
return( ans )
}
pcopula.student = function(U, Corr, df, ...)
{
m = dim(Corr)[2]
U = matrix(U, ncol = m)
Z = apply(U, 2, FUN = function(x) qt(x, df = df))
mu = rep(0, m)
ans = apply(Z, 1, FUN = function(x) pmvt(lower = rep(-Inf, m), upper = x, delta = mu, corr = Corr, df = df, ...))
return( ans )
}
rcopula.student = function(n, U, Corr, df)
{
m = dim(Corr)[2]
mu = rep(0, m)
ans = pt(rmvt(n, delta = mu, sigma = Corr, df = df), df = df)
return ( ans )
}
#####################################################################################
# Copula Simulation
#------------------------------------------------------------------------------------
.sample.copula = function(dist, object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
dist = tolower(dist)
set.seed(rseed)
ans = switch(dist,
mvt0 = .rtcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control),
mvt1 = .rtvtcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control),
mvnorm0 = .rnormcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control),
mvnorm1 = .rtvnormcopula(object, Rbar, n.sim, n.start, m.sim, rseed, parallel, parallel.control))
return(ans)
}
.rtvtcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
dccOrder = object@mspec$dccOrder
cf = coef(object, type = "st")
shape = cf["shape"]
mo = max( dccOrder )
dcca = cf[1:dccOrder[1]]
dccb = cf[(dccOrder[1]+1):sum(dccOrder[1:2])]
preQ = Rbar
m = dim(Rbar)[1]
z = array(NA, dim = c(n.sim + n.start + mo, m, m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
z[,,i] = rmvt(n = (n.sim + n.start + mo), delta = rep(0, m), sigma = diag(m), df = shape)
}
xseed = rseed+1
if( parallel ){
os = .Platform$OS.type
if(is.null(parallel.control$pkg)){
if( os == "windows" ) parallel.control$pkg = "snowfall" else parallel.control$pkg = "multicore"
if( is.null(parallel.control$cores) ) parallel.control$cores = 2
} else{
mtype = match(tolower(parallel.control$pkg[1]), c("multicore", "snowfall"))
if(is.na(mtype)) stop("\nParallel Package type not recognized in parallel.control\n")
parallel.control$pkg = tolower(parallel.control$pkg[1])
if( os == "windows" && parallel.control$pkg == "multicore" ) stop("\nmulticore not supported on windows O/S\n")
if( is.null(parallel.control$cores) ) parallel.control$cores = 2 else parallel.control$cores = as.integer(parallel.control$cores[1])
}
if( parallel.control$pkg == "multicore" ){
mtmp = mclapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]),
mc.cores = parallel.control$cores)
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pt(sapply(mtmp, FUN = function(x) x$stdresid[,i]), shape), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
} else{
sfInit(parallel = TRUE, cpus = parallel.control$cores)
sfExport("z", "preQ", "Rbar", "dcca", "dccb", "mo", "n.sim", "n.start", "m", "xseed", local = TRUE)
mtmp = sfLapply(as.list(1:m.sim), fun = function(j)
rgarch:::.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
sfStop()
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pt(sapply(mtmp, FUN = function(x) x$stdresid[,i]), shape), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
}
} else{
mtmp = lapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = pt(matrix(sapply(mtmp, FUN = function(x) x$stdresid[,i]), ncol = m.sim), shape)
for(i in 1:m.sim) Usim[,,i] = matrix(sapply(Ures, FUN = function(x) x[-(1:mo),i]), ncol = m)
}
return(list(Usim = Usim, simR = simR))
}
.rtcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
nsim = n.sim + n.start
m = dim(Rbar)[1]
shape = object@mfit$stpars["shape"]
sim = array(data = NA, dim = c(nsim , m , m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
tmp = rmvt(n = nsim, delta = rep(0, m), sigma = Rbar, df = shape)
sim[,,i] = matrix(pt(tmp, df = shape), nrow = nsim, ncol = m)
}
return( sim )
}
.rtvnormcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
dccOrder = object@mspec$dccOrder
cf = coef(object, type = "st")
mo = max( dccOrder )
dcca = cf[1:dccOrder[1]]
dccb = cf[(dccOrder[1]+1):sum(dccOrder[1:2])]
preQ = Rbar
m = dim(Rbar)[1]
z = array(NA, dim = c(n.sim + n.start + mo, m, m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
z[,,i] = rmvnorm(n = (n.sim + n.start + mo), mean = rep(0, m), sigma = diag(m))
}
xseed = rseed+1
if( parallel ){
os = .Platform$OS.type
if(is.null(parallel.control$pkg)){
if( os == "windows" ) parallel.control$pkg = "snowfall" else parallel.control$pkg = "multicore"
if( is.null(parallel.control$cores) ) parallel.control$cores = 2
} else{
mtype = match(tolower(parallel.control$pkg[1]), c("multicore", "snowfall"))
if(is.na(mtype)) stop("\nParallel Package type not recognized in parallel.control\n")
parallel.control$pkg = tolower(parallel.control$pkg[1])
if( os == "windows" && parallel.control$pkg == "multicore" ) stop("\nmulticore not supported on windows O/S\n")
if( is.null(parallel.control$cores) ) parallel.control$cores = 2 else parallel.control$cores = as.integer(parallel.control$cores[1])
}
if( parallel.control$pkg == "multicore" ){
mtmp = mclapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]),
mc.cores = parallel.control$cores)
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pnorm(sapply(mtmp, FUN = function(x) x$stdresid[,i])), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
} else{
sfInit(parallel = TRUE, cpus = parallel.control$cores)
sfExport("z", "preQ", "Rbar", "dcca", "dccb", "mo", "n.sim", "n.start", "m", "xseed", local = TRUE)
mtmp = sfLapply(as.list(1:m.sim), fun = function(j)
rgarch:::.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
sfStop()
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = matrix(pnorm(sapply(mtmp, FUN = function(x) x$stdresid[,i])), ncol = m.sim)
for(i in 1:m.sim) Usim[,,i] = sapply(Ures, FUN = function(x) x[-(1:mo),i])
}
} else{
mtmp = lapply(as.list(1:m.sim), FUN = function(j)
.dccsimf(Z = z[,,j], Qbar = preQ, Rbar = Rbar, dcca, dccb, mo, n.sim, n.start, m, xseed[j]))
simR = lapply(mtmp, FUN = function(x) if(is.matrix(x$R)) array(x$R, dim = c(m, m, n.sim)) else last(x$R, n.sim))
Ures = vector(mode = "list", length = m)
Usim = array(NA, dim = c(n.sim+n.start, m, m.sim))
for(i in 1:m) Ures[[i]] = pnorm(matrix(sapply(mtmp, FUN = function(x) x$stdresid[,i]), ncol = m.sim))
for(i in 1:m.sim) Usim[,,i] = matrix(sapply(Ures, FUN = function(x) x[-(1:mo),i]), ncol = m)
}
return(list(Usim = Usim, simR = simR))
}
.rnormcopula = function(object, Rbar, n.sim, n.start = 0, m.sim, rseed, parallel = FALSE, parallel.control)
{
nsim = n.sim + n.start
m = dim(Rbar)[1]
sim = array(data = NA, dim = c(nsim , m , m.sim))
for(i in 1:m.sim){
set.seed(rseed[i])
tmp = rmvnorm(n = nsim, mean = rep(0, m), sigma = Rbar)
sim[,,i] = matrix(pnorm(tmp), nrow = nsim, ncol = m)
}
return( sim )
}
#------------------------------------------------------------------------------------
Try the rgarch package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.