Compute the measure of permutation asymmetry, which can be thought of as bivariate skewness, named for the copBasic package as Nu-Skew \(\nu_\mathbf{C}\) of a copula according to Joe (2014, p. 66) by
$$\nu_\mathbf{C} = 3\mathrm{E}[UV^2 - U^2V] = 6\int\!\!\int_{\mathcal{I}^2} (v-u)\mathbf{C}(u,v)\, \mathrm{d}u\mathrm{d}v\mbox{.}$$
This definition is effectively the type="nu" for the function for which the multiplier \(6\) has been converted to \(96\) as explained in the Note.
Numerical results indicate \(\nu_\mathbf{W} \approx 0\) (W), \(\nu_\mathbf{\Pi} = 0\) (P), \(\nu_\mathbf{M} \approx 0\) (M), \(\nu_\mathbf{PL} \approx 0\) for all \(\Theta\) (PLcop), and the \(\nu^\star_\mathbf{GH} = 0\) (GHcop); copulas with mirror symmetry across the equal value line have \(\nu_\mathbf{C} = 0\).
Asymmetric copulas do exist. For example, consider an asymmetric Gumbel--Hougaard \(\mathbf{GH}\) copula with \(\Theta_p = (5,0.8,p)\):
optimize(function(p) { nuskewCOP(cop=GHcop, para=c(5,0.8, p)) },
c(0,0.99) )$minimum
UV <- simCOP(n=10000, cop=GHcop, c(5,0.8, 0.2836485)) # inspect the graphics
48*mean(UV$U*$V^2 - UV$U^2*UV$V) # -0.2847953 (not the 3rd parameter)
The minimization yields \(\nu_{\mathbf{GH}(5, 0.8, 0.2836485)} = -0.2796104\), which is close the expectation computed where \(48 = 96/2\).
A complementary definition is supported, triggered by type="nustar", and is computed by
$$\nu^\star_\mathbf{C} = 12\int\!\!\int_{\mathcal{I}^2} (v+u)\mathbf{C}(u,v)\, \mathrm{d}u\mathrm{d}v - 4\mbox{,}$$
which has been for the copBasic package, \(\nu^\star_\mathbf{C}\) is named as Nu-Star, which the \(12\) and the \(-4\) have been chosen so that numerical results indicate \(\nu^\star_\mathbf{W} = -1\) (W), \(\nu^\star_\mathbf{\Pi} = 0\) (P), and \(\nu^\star_\mathbf{M} = +1\) (M).
Lastly, the uvlmoms function provides for a quantile-based measure of bivariate skewness based on the difference \(U - V\) that also is discussed by Joe (2014, p. 66).
joeskewCOP(cop=NULL, para=NULL, type=c("nu", "nustar", "nuskew"),
as.sample=FALSE, brute=FALSE, delta=0.002, ...)nuskewCOP(cop=NULL, para=NULL, as.sample=FALSE, brute=FALSE, delta=0.002, ...)
nustarCOP(cop=NULL, para=NULL, as.sample=FALSE, brute=FALSE, delta=0.002, ...)
The value for \(\nu_\mathbf{C}\) or \(\nu^\star_\mathbf{C}\) is returned.
A copula function;
Vector of parameters or other data structure, if needed, to pass to the copula;
The type of metric to compute (nu and nuskew are synonymous for \(\nu_\mathbf{C}\) and nustar is for \(\nu^\star_\mathbf{C}\));
Should brute force be used instead of two nested integrate() functions to perform the double integration;
The \(\mathrm{d}u\) and \(\mathrm{d}v\) for the brute force integration using brute;
A logical controlling whether an optional R data.frame in para is used to compute the sample \(\hat\nu\) or \(\hat\nu^\star\) (see Note). If set to -1, then the message concerning CPU effort will be surpressed; and
Additional arguments to pass.
W.H. Asquith
The implementation of joeskewCOP for copBasic provides the second metric of asymmetry, but why? Consider the results that follow:
joeskewCOP(cop=GHcop, para=c(5, 0.8, 0.2836485), type="nu")
# -0.2796104
joeskewCOP(cop=GHcop, para=c(5, 0.2836485, 0.8), type="nu")
# +0.2796103
joeskewCOP(cop=GHcop, para=c(5, 0.8, 0.2836485), type="nu")
# 0.3571276
joeskewCOP(cop=GHcop, para=c(5, 0.2836485, 0.8), type="nu")
# 0.3571279
tauCOP( cop=GHcop, para=c(5, 0.2836485, 0.8))
# 0.2443377
The demonstration shows---at least for the symmetry (switchability) of the 2nd and 3rd parameters (\(\pi_2\) and \(\pi_3\)) of the asymmetric \(\mathbf{GH}\) copula---that the first definition \(\nu\) is magnitude symmetric but carries a sign change. The demonstration shows magnitude and sign stability for \(\nu^\star\), and ends with Kendall Tau (tauCOP). Collectively, Kendall Tau (or the other symmetric measures of association, e.g. blomCOP, footCOP, giniCOP, hoefCOP, rhoCOP, wolfCOP) when combined with \(\nu\) and \(\nu^\star\) might provide a framework for parameter optimization of the asymmetric \(\mathbf{GH}\) copula (see below).
The asymmetric \(\mathbf{GH}_{(5, 0.2836485, 0.8)}\) is not radial (isCOP.radsym) or permutation (isCOP.permsym), but if \(\pi_2 = \pi_3\) then the resulting \(\mathbf{GH}\) copula is not radially symmetric but is permutation symmetric:
isCOP.radsym( cop=GHcop, para=c(5, 0.2836485, 0.8)) # FALSE
isCOP.permsym(cop=GHcop, para=c(5, 0.2836485, 0.8)) # FALSE
isCOP.radsym( cop=GHcop, para=c(5, 0.8, 0.8)) # FALSE
isCOP.permsym(cop=GHcop, para=c(5, 0.8, 0.8)) # TRUE
The use of \(\nu_\mathbf{C}\) and \(\nu^\star_\mathbf{C}\) with a measure of association is just suggested above for parameter optimization. Suppose we have \(\mathbf{GH}_{(5,0.5,0.7)}\) with Spearman Rho \(\rho = 0.4888\), \(\nu = 0.001475\), and \(\nu^\star = 0.04223\), and the asymmetric \(\mathbf{GH}\) coupla is to be fit. Parameter estimation for the asymmetric \(\mathbf{GH}\) is accomplished by
"fitGHcop" <- function(hats, assocfunc=rhoCOP, init=NA, eps=1E-4, ...) {
H <- GHcop # shorthand for the copula
"objfunc" <- function(par) {
par[1] <- ifelse(par[1] < 1, return(Inf), exp(par[1])) # edge check
par[2:3] <- pnorm(par[2:3]) # detransform
hp <- c(assocfunc(H, par), nuskewCOP(H, par), nustarCOP(H, par))
return(sum((hats-hp)^2))
}
# Theta=1 and Pi2 = Pi3 = 1/2 # as default initial estimates
if(is.na(init)) init <- c(1, rep(1/2, times=2))
opt <- optim(init, objfunc, ...); par <- opt$par
para <- c( exp(par[1]), pnorm(par[2:3]) )
names(para) <- c("Theta", "Pi2", "Pi3")
fit <- c(assocfunc(H, para), nuskewCOP(H, para), nustarCOP(H, para))
txt <- c("AssocMeasure", "NuSkew", "NuStar")
names(fit) <- txt; names(hats) <- txt
if(opt$value > eps) warning("inspect the fit")
return(list(para=para, fit=fit, given=hats, optim=opt))
}
father <- c(5,.5,.7)
densityCOPplot(cop=GHcop, para=father, contour.col=8)
fRho <- rhoCOP( cop=GHcop, father)
fNu <- nuskewCOP(cop=GHcop, father)
fStar <- nustarCOP(cop=GHcop, father) child <- fitGHcop(c(fRho, fNu, fStar))$para
densityCOPplot(cop=GHcop, para=child, ploton=FALSE)
cRho <- rhoCOP( cop=GHcop, child)
cNu <- nuskewCOP(cop=GHcop, child)
cStar <- nustarCOP(cop=GHcop, child)
message("Father stats: ", paste(fRho, fNu, fStar, sep=", "))
message("Child stats: ", paste(cRho, cNu, cStar, sep=", "))
message("Father para: ", paste(father, collapse=", "))
message("Child para: ", paste(child, collapse=", "))
The initial parameter estimate has the value \(\Theta = 1\), which is independence for the one parameter \(\mathbf{GH}\). The two other parameters are set as \(\pi_2 = \pi_3 = 1/2\) to be in the mid-point of their domain. The transformations using the log() \(\leftrightarrow\) exp() and qnorm() \(\leftrightarrow\) pnorm() functions in R are used to keep the optimization in the viable parameter domain. The results produce a fitted copula of \(\mathbf{GH}_{(4.907, 0.5006, 0.7014)}\). This fit aligns well with the parent, and the three statistics are essentially matched during the fitting.
The \(\nu^\star_\mathbf{C}\) can be similar to rhoCOP, but differences do exist. In the presence of radial symmetry, (\(\nu_\mathbf{C} == 0\)), the \(\nu^\star_\mathbf{C}\) is nearly equal to Spearman Rho for some copulas. Let us test further:
p <- 10^seq(0,2,by=.01)
s <- sapply(p, function(t) nustarCOP(cop=GHcop, para=c(t)))
r <- sapply(p, function(t) rhoCOP(cop=GHcop, para=c(t)))
plot(p,s, log="x", type="l", col=3, lwd=3); lines(p,r)
Now let us add some asymmetry
s <- sapply(p, function(t) nustarCOP(cop=GHcop, para=c(t, 0.25, 0.75)))
r <- sapply(p, function(t) rhoCOP(cop=GHcop, para=c(t, 0.25, 0.75)))
plot(p,s, log="x", type="l", col=3, lwd=3); lines(p,r)
Now let us choose a different (the Clayton) copula
s <- sapply(p, function(t) nustarCOP(cop=CLcop, para=c(t)))
r <- sapply(p, function(t) rhoCOP(cop=CLcop, para=c(t)))
plot(p,s, log="x", type="l", col=3, lwd=3); lines(p,r)
Joe, H., 2014, Dependence modeling with copulas: Boca Raton, CRC Press, 462 p.
uvskew, blomCOP, footCOP, giniCOP,
hoefCOP, rhoCOP, tauCOP, wolfCOP