mixCopula {copula}

Create Mixture of Copulas

Description

A mixture of m copulas of dimension d with weights w_j, j=1,2,\ldots,m is itself a d-dimensional copula, with cumulative distribution function

C(x) = \sum_{j=1}^m w_j C_j(x),

and (probability) density function

c(x) = \sum_{j=1}^m w_j c_j(x),

where C_j are the CDFs and c_j are the densities of the m component copulas, j=1,2,\ldots,m.

Usage

mixCopula(coplist, w = NULL)

Arguments

coplist

a list of length m (\ge 1) copulas (each inheriting from parCopula), all of the same dimension.

w

numeric vector of length m of non-negative mixture weights, or NULL, which means equal weights.

Details

It easy to see that the tail dependencies lambda() and Spearman's rank correlation rho() can be computed as mixture of the individual measures.

Value

an object of class mixCopula

See Also

khoudrajiCopula, rotCopula also create new copula models from existing ones.

Examples

mC <- mixCopula(list(gumbelCopula(2.5, dim=3),
                     claytonCopula(pi, dim=3),
                     tCopula(0.7, dim=3)),
                c(2,2,4)/8)
mC
stopifnot(dim(mC) == 3)

set.seed(17)
uM <- rCopula(600, mC)
splom2(uM, main = "mixCopula( (gumbel, clayton, t-Cop) )")
d.uM <- dCopula(uM, mC)
p.uM <- pCopula(uM, mC)

## mix a Gumbel with a rotated Gumbel (with equal weights 1/2):
mGG <- mixCopula(list(gumbelCopula(2), rotCopula(gumbelCopula(1.5))))
rho(mGG)   # 0.57886
lambda(mGG)# both lower and upper tail dependency

loglikCopula(c(2.5, pi, rho.1=0.7, df = 4, w = c(2,2,4)/8),
                    u = uM, copula = mC)
## define (profiled) log-likelihood function of two arguments (df, rho) :
ll.df <- Vectorize(function(df, rho)
                   loglikCopula(c(2.5, pi, rho.1=rho, df=df, w = c(2,2,4)/8),
                                uM, mC))
(df. <- 1/rev(seq(1/8, 1/2, length=21)))# in [2, 8] equidistant in 1/. scale

ll. <- ll.df(df., rho = (rh1 <- 0.7))
plot(df., ll., type = "b", main = "loglikCopula((.,.,rho = 0.7, df, ..), u, <mixCopula>)")

if(!exists("Xtras")) Xtras <- copula:::doExtras() ; cat("Xtras: ", Xtras,"\n")
if (Xtras) withAutoprint({
  Rhos <- seq(0.55, 0.70, by = 0.01)
  ll.m <- matrix(NA, nrow=length(df.), ncol=length(Rhos))
  for(k in seq_along(Rhos)) ll.m[,k] <- ll.df(df., rho = Rhos[k])
  tit <- "loglikelihood(<tCop>, true param. for rest)"
  persp         (df., Rhos, ll.m, phi=30, theta = 50, ticktype="detailed", main = tit)
  filled.contour(df., Rhos, ll.m, xlab="df", ylab = "rho", main = tit)
})

## fitCopula() example -----------------------------------------------------

## 1) with "fixed" weights :

(mNt <- mixCopula(list(normalCopula(0.95), tCopula(-0.7)), w = c(1, 2) / 3))
set.seed(1452) ; U <- pobs(rCopula(1000, mNt))
(m1 <- mixCopula(list(normalCopula(), tCopula()), w = mNt@w))

getTheta(m1, freeOnly = TRUE, attr = TRUE)
getTheta(m1, named=TRUE)
isFree(m1) # all of them; --> now fix the weights :
fixedParam(m1) <- fx <- c(FALSE, FALSE, FALSE, TRUE, TRUE)
stopifnot(identical(isFree(m1), !fx))

if(Xtras) withAutoprint({ ## time
  system.time( # ~ 16 sec (nb-mm4) :
    fit <- fitCopula(m1, start = c(0, 0, 10), data = U)
  )
  fit
  summary(fit) #-> incl  'Std.Error' (which seems small for rho1 !)
})

## 2) with "free" weights (possible since copula 1.0-0):

(mNt2 <- mixCopula(list(normalCopula(0.9), tCopula(-0.8)), w = c(1, 3) / 4))
set.seed(1959) ; U2 <- pobs(rCopula(2000, mNt2))
if(Xtras) withAutoprint({ ## time
  m2 <- mixCopula(list(normalCopula(), tCopula()), w = mNt@w)
  system.time( # ~ 13.5 sec (lynne) :
    f2 <- fitCopula(m2, start = c(0, 0, 10, c(1/2, 1/2)), data = U2)
  )
  f2
  summary(f2) # NA for 'Std. Error' as did *not* estimate.variance
  summary(f2, orig=FALSE) # default 'orig=TRUE': w-scale;  whereas
     coef(f2, orig=FALSE) # 'orig=FALSE' => shows 'l-scale' instead
})