copula (version 0.999-1)

K: Kendall distribution function for Archimedean copulas

Description

The distribution function of the Kendall distribution of an Archimedean copula is defined as $$K(u) = P(C(U_1,U_2,\dots,U_d) \le u),$$ where $u \in [0,1]$, and the $d$-dimensional $(U_1,U_2,\dots,U_d)$ is distributed according to the copula $C$. Note that the random variable $C(U_1,U_2,\dots,U_d)$ is known as probability integral transform. Its distribution function $K$ is equal to the identity if $d = 1$, but is non-trivial for $d \ge 2$.

pK, qK, dK, and rK provide the distribution function, quantile function, density, and random number generator, respectively, for the Kendall distribution of an Archimedean copula.

Usage

pK(u, cop, d, n.MC=0, log=FALSE)
qK(u, cop, d, n.MC=0,
   method=c("default", "simple", "sort", "discrete", "monoH.FC"),
   u.grid, ...)
dK(u, cop, d, n.MC=0, log=FALSE)
rK(n, cop, d)

Arguments

u

evaluation point(s) (have to be in $[0,1]$).

cop

acopula with specified parameter, or (currently for rK only) a outer_nacopula.

d

dimension (not used when cop is an outer_nacopula).

n.MC

integer, if positive, a Monte Carlo approach is applied with sample size equal to n.MC to evaluate the generator derivatives involved; otherwise (n.MC=0) the exact form

log

logical indicating if the log of $K$ should be returned.

method

string for the method to compute the quantile function of $K$. Currently, one of [object Object],[object Object],[object Object],[object Object],[object Object]

u.grid

(for method="discrete":) the grid on which $K$ is evaluated, a numeric vector.

...

additional arguments passed to uniroot (for method="default", method="simple", method="sort", and method="monoH.FC") or

n

sample size for rK.

Value

• The distribution function, quantile function, density, and random numbers of a Kendall distribution.

Details

For a completely monotone Archimedean generator $\psi$, $$K(u)=\sum_{k=0}^{d-1} \frac{\psi^{(k)}(\psi^{-1}(u))}{k!} (-\psi^{-1}(u))^k,\ u\in[0,1];$$ see Barbe et al. (1996). The corresponding density is $$\frac{(-1)^d\psi^{(d)}(\psi^{-1}(u))}{(d-1)!} (-(\psi^{-1})'(u))(\psi^{-1}(u))^{d-1}$$

References

Barbe, P., Genest, C., Ghoudi, K., and Ré{e}millard, B. (1996), On Kendall's Process, Journal of Multivariate Analysis 58, 197--229.

Hofert, M., Mächler{Maechler}, M., and McNeil, A. J. (2011b), Likelihood inference for Archimedean copulas, submitted.

See Also

gnacopula, htrafo or emde (where K is used).

Examples

tau <- 0.5
(theta <- copGumbel@tauInv(tau)) # 2
d <- 20
(cop <- onacopulaL("Gumbel", list(theta,1:d)))

## compute Kendall distribution function
u <- seq(0,1, length.out = 255)
Ku    <- pK(u, cop=cop@copula, d=d) # exact
Ku.MC <- pK(u, cop=cop@copula, d=d, n.MC=1000) # via Monte Carlo

## build sample from K
set.seed(1)
n <- 200
W <- rK(n, cop)

## plot empirical distribution function based on W
## and the corresponding theoretical Kendall distribution function
## (exact and via Monte Carlo)
plot(ecdf(W), col="blue", xlim=c(0,1), verticals=TRUE,
     main = expression("Empirical"~ F[n]( C(U) ) ~
                       "and its Kendall distribution"~ K(u)),
     do.points=FALSE, asp=1)
abline(0,1, lty=2); abline(h=0:1, v=0:1, lty=3, col="gray")
lines(u, Ku.MC, col="red")
lines(u, Ku, col="black")
legend(.2,.6, expression(F[n],K(u), K[MC](u)),
       col=c("blue","red","black"), lty=1, bty="n",
       xjust = 0.25, yjust = 0.5)

### testing qK
uexpr <- quote( 0:63/63 );  u <- eval(uexpr)
d <- 10
cop <- onacopulaL("Gumbel", list(theta = 2, 1:d))

system.time(Ku1 <- qK(u, cop=cop@copula, d=d, method="simple"))
system.time(Ku2 <- qK(u, cop=cop@copula, d=d, method="sort"))
system.time(Ku3 <- qK(u, cop=cop@copula, d=d, method="discrete",
                      u.grid=0:1e4/1e4))
system.time(Ku4 <- qK(u, cop=cop@copula, d=d, method="monoH.FC",
                      u.grid=0:5e2/5e2))

cols <- adjustcolor(c("gray80", "light blue", "royal blue", "purple3"), 0.6)
matplot(u, cbind(Ku1,Ku2,Ku3,Ku4), type="l", lwd=2*4:1, lty = 1:4, col= cols,
        xlab=substitute(u == U, list(U=uexpr)), ylab=expression({K^{-1}}(u)))
legend("topleft", col=cols, lwd=2*4:1, lty=1:4, bty="n", inset=.03,
       legend= paste0("method= ",
             sQuote(c("simple", "sort", "discrete", "monoH.FC"))))