ellipCopula: Construction of Elliptical Copula Class Objects
In copula: Multivariate Dependence with Copulas

ellipCopula

R Documentation

Construction of Elliptical Copula Class Objects

Description

Creating elliptical copula objects with corresponding dimension and parameters, including the dispersion structure P (pronounced “Rho”).

Usage

ellipCopula (family, param, dim = 2, dispstr = "ex", df = 4, ...)

normalCopula(param, dim = 2, dispstr = "ex")

     tCopula(param, dim = 2, dispstr = "ex", df = 4,
             df.fixed = FALSE, df.min = 0.01)

dispstrToep(perm = NULL, check = TRUE)

## S4 method for signature 'matrix,normalCopula'
pCopula(u, copula, algorithm=NULL, keepAttr=FALSE, ...)

## S4 method for signature 'matrix,tCopula'
pCopula(u, copula, algorithm=NULL, keepAttr=FALSE, ...)

Arguments

`family`	a `character` string specifying the family of an elliptical copula. Must be `"normal"` (the default) or `"t"`.
`param`	a `numeric` vector specifying the parameter values; `P2p()` accesses this vector, whereas `p2P()` and `getSigma()` provide the corresponding “P” matrix, see below.
`dim`	the dimension of the copula.
`dispstr`	a string specifying the “dispersion structure”, i.e., type of the symmetric positive definite matrix characterizing the elliptical copula. Currently available structures are `"ex"` for exchangeable, `"ar1"` for `AR(1)`, `"toep"` for Toeplitz (`toeplitz`), and `"un"` for unstructured. The dispersion structure for Toeplitz can (and often should) now be specified by `dispstrToep()`, see there.
`df`	integer value specifying the number of degrees of freedom of the multivariate t distribution used to construct the t copulas.
`df.fixed`	logical specifying if the degrees of freedom `df` will be considered as a parameter (to be estimated) or not. The default, `FALSE`, means that `df` is to be estimated if the object is passed as argument to `fitCopula`.
`df.min`	non-negative number; the strict lower bound for `df`, mainly during fitting when `df.fixed=FALSE`, with `fitCopula`.
`copula`	an R object of class `"Copula"`, in our case inheriting from `"ellipCopula"`.
`u`	a vector of the copula dimension `d` or a matrix with `d` columns, giving the points where the distribution function needs to be evaluated. Note that values outside of the cube `[0,1]^d` are treated equivalently to those on the cube boundary.
`algorithm`	`NULL` or an `"algorithm"` object for package mvtnorm's `pmvt()` or `pmvnorm()` functions, see `algorithms`. Note that for larger dimensions, the monte-carlo based `GenzBretz(..)` must be used, consequently with slightly random results. By default, `algorithm = NULL`, `algorithm` is chosen separately for each row `x <- u[i,]`, for `normalCopula`: via hidden function `pmvnormAlgo(dim, x, ...)` which currently is defined as pmvnormAlgo <- function(dim, x, ...) { if(dim <= 3 && !anyNA(x) && (!any(xI <- x == Inf) \|\| all(xI))) TVPACK(...) else if(dim <= 5) Miwa(...) else GenzBretz(...) } `tCopula`: via (hidden) `pmvtAlgo(dim, x, ...)` which is the same as `pmvnormAlgo()` above, but as `Miwa()` is not applicable, without the `else if(dim <= 5) Miwa(...)` clause.
`keepAttr`	`logical` passed to `pmvnorm` or `pmvt`, respectively.
`...`	for the `pCopula()` methods, optional arguments to the corresponding algorithm.
`perm`	an `integer` vector of length `d =dim`, which must be a permutation of `1:d` specifying the (column) ordering of the variables which has Toeplitz dispersion.
`check`	a `logical` specifying if the validity of `perm` should be checked (not strictly, currently).

Value

For

ellipCopula(), normalCopula(), or tCopula():

an elliptical copula object of class "normalCopula" or "tCopula".

dispstrToep():

the character string "toep", optionally with attribute (see attributes, attr) "perm" with a permutation p of 1:d, such that (the column permuted) U_p, or in the data case the column-permuted matrix U[,p] has as dispersion matrix toeplitz(c(1, par)), with par the respective parameter vector of bivariate “correlations” \rho_{i,j}.

Note that the result of dispstrToep() is currently stored in the dispstr slot of the copula object.

pCopula(u, *):

the numerical vector of copula probabilites of length nrow(u).

Note

ellipCopula() is a wrapper for normalCopula() and tCopula().

The pCopula() methods for the normal- and t-copulas accept optional arguments to be passed to the underlying (numerical integration) algorithms from package mvtnorm's pmvnorm and pmvt, respectively, notably algorithm, see GenzBretz, or abseps which defaults to 0.001.

Examples

normalCopula(c(0.5, 0.6, 0.7), dim = 3, dispstr = "un")
t.cop <- tCopula(c(0.5, 0.3), dim = 3, dispstr = "toep",
                 df = 2, df.fixed = TRUE)
getSigma(t.cop) # P matrix (with diagonal = 1)
stopifnot(all.equal(toeplitz(c(1, .5, .3)), getSigma(t.cop)))

## dispersion "AR1" :
nC.7 <- normalCopula(0.8, dim = 7, dispstr = "ar1")
getSigma(nC.7)
stopifnot(all.equal(toeplitz(.8^(0:6)), getSigma(nC.7)))

## from the wrapper
norm.cop <- ellipCopula("normal", param = c(0.5, 0.6, 0.7),
                        dim = 3, dispstr = "un")
if(require("scatterplot3d") && dev.interactive(orNone=TRUE)) {
  ## 3d scatter plot of 1000 random observations
  scatterplot3d(rCopula(1000, norm.cop))
  scatterplot3d(rCopula(1000, t.cop))
}
set.seed(12)
uN <- rCopula(512, norm.cop)
set.seed(2); pN1 <- pCopula(uN, norm.cop)
set.seed(3); pN2 <- pCopula(uN, norm.cop)
stopifnot(identical(pN1, pN2)) # no longer random for dim = 3
(Xtras <- copula:::doExtras())
if(Xtras) { ## a bit more accurately:
  set.seed(4); pN1. <- pCopula(uN, norm.cop, abseps = 1e-9)
  set.seed(5); pN2. <- pCopula(uN, norm.cop, abseps = 1e-9)
  stopifnot(all.equal(pN1., pN2., 1e-5))# see 3.397e-6
  ## but increasing the required precision (e.g., abseps=1e-15) does *NOT* help
}

## For smaller copula dimension 'd', alternatives are available and
## non-random, see ?GenzBretz from package 'mvtnorm' :
has_mvtn <- "package:mvtnorm" %in% search() #% << (workaround ESS Rd render bug)
if(!has_mvtn)
  require("mvtnorm")# -> GenzBretz(), Miva(), and TVPACK() are available
## Note that Miwa() would become very slow for dimensions 5, 6, ..
set.seed(4); pN1.M <- pCopula(uN, norm.cop, algorithm = Miwa(steps = 512))
set.seed(5); pN2.M <- pCopula(uN, norm.cop, algorithm = Miwa(steps = 512))
stopifnot(all.equal(pN1.M, pN2.M, tol= 1e-15))# *no* randomness
set.seed(4); pN1.T <- pCopula(uN, norm.cop, algorithm = TVPACK(abseps = 1e-10))
set.seed(5); pN2.T <- pCopula(uN, norm.cop, algorithm = TVPACK(abseps = 1e-14))
stopifnot(all.equal(pN1.T, pN2.T, tol= 1e-15))# *no* randomness (but no effect of 'abseps')
if(!has_mvtn)
   detach("package:mvtnorm")# (revert)


## Versions with unspecified parameters:
tCopula()
allEQ <- function(u,v) all.equal(u, v, tolerance=0)
stopifnot(allEQ(ellipCopula("norm"), normalCopula()),
          allEQ(ellipCopula("t"), tCopula()))
tCopula(dim=3)
tCopula(dim=4, df.fixed=TRUE)
tCopula(dim=5, disp = "toep", df.fixed=TRUE)
normalCopula(dim=4, disp = "un")

## Toeplitz after *permutation* dispersions (new in copula 1.1-0) ---------
tpar <- c(7,5,3)/8 # *gives* pos.def.:
(ev <- eigen(toeplitz(c(1, tpar)), symmetric=TRUE, only.values=TRUE)$values)
stopifnot(ev > 0)
N4.   <- ellipCopula("normal", dim=4, param=tpar, dispstr = "toep") #"regular"
## reversed order is "the same" for toeplitz structure:
N4.pr <- ellipCopula("normal", dim=4, param=tpar, dispstr = dispstrToep(4:1))
N4.p1 <- ellipCopula("normal", dim=4, param=tpar, dispstr = dispstrToep(c(4,1:3)))
N4.p2 <- ellipCopula("normal", dim=4, param=tpar, dispstr = dispstrToep(c(4:3,1:2)))
N4.p3 <- ellipCopula("normal", dim=4, param=tpar, dispstr = dispstrToep(c(2,4,1,3)))

(pm <- attr(N4.p3@dispstr, "perm")) # (2 4 1 3)
ip <- c(3,1,4,2) # the *inverse* permutation of (2 4 1 3) = Matrix::invPerm(pm)
(Sp3 <- getSigma(N4.p3)) # <-- "permuted toeplitz"
Sp3[ip, ip] # re-ordered rows & columns => *is* toeplitz :
stopifnot(exprs = {
  ## permutations  pm  and  ip  are inverses:
  pm[ip] == 1:4
  ip[pm] == 1:4
  is.matrix(T4 <- toeplitz(c(1, tpar)))
  identical(getSigma(N4.),   T4)
  identical(getSigma(N4.pr), T4) # 4:1 and 1:4 is "the same" for Rho
  identical(Sp3[ip, ip]  ,   T4)
  identical(Sp3,      T4[pm,pm])
})
## Data generation -- NB: The U matrices are equal only "in distribution":
set.seed(7); U.p3 <- rCopula(1000, N4.p3)
set.seed(7); U.   <- rCopula(1000, N4.)
stopifnot(exprs = {
 all.equal(loglikCopula(tpar, u=U.p3,      copula= N4.p3),
           loglikCopula(tpar, u=U.p3[,ip], copula= N4.) -> LL3)
 all.equal(loglikCopula(tpar, u=U.,      copula= N4.),
           loglikCopula(tpar, u=U.[,pm], copula= N4.p3) -> LL.)
})
c(LL. , LL3)# similar but different
if(Xtras) {
  fm.  <- fitCopula(N4.  , U.  )
  fm.3 <- fitCopula(N4.p3, U.p3)
  summary(fm.3)
  stopifnot(all.equal(coef(fm.), coef(fm.3), tol = 0.01))# similar but different
}

copula documentation built on Sept. 11, 2024, 7:48 p.m.