# NOT RUN {
kfuncCOPlmom(1, cop=P) # 0.5 * 0.5 = 0.25 is expected joint prob. of independence
#[1] 0.2499999 (in agreement with theory)
ls.str(kfuncCOPlmoms(cop=GHcop, para=4.21)) # Gumbel-Hougaard copula
# lambdas : num [1:5] 0.440617 0.169085 0.011228 -0.000797 0.000249
# ratios : num [1:5] NA 0.38375 0.0664 -0.00472 0.00147
# source : chr "kfuncCOPlmoms" # e.g. L-skew = 0.0664
# }
# NOT RUN {
# }
# NOT RUN {
UV <- simCOP(200, cop=PLcop, para=1/pi, graphics=FALSE)
theta <- PLpar(UV[,1], UV[,2]); zs <- seq(0.01,0.99, by=.01) # for later
# Take the sample estimated parameter and convert to joint probabilities Z
# Convert the Z to the Kendall Function estimates again with the sample parameter
Z <- PLcop(UV[,1], UV[,2], para=theta); KF <- kfuncCOP(Z, cop=PLcop, para=theta)
# Compute L-moments of the "Kendall function" and the sample versions
# Note again though that the L-moment are for the distribution of the Z!
KNFlmr <- kfuncCOPlmoms(cop=PLcop, para=theta); SAMlmr <- lmomco::lmoms(Z)
knftxt <- paste0("Kendall L-moments: ",
paste(round(KNFlmr$lambdas, digits=4), collapse=", "))
samtxt <- paste0("Sample L-moments: " ,
paste(round(SAMlmr$lambdas, digits=4), collapse=", "))
plot(Z, KF, xlim=c(0,1), ylim=c(0,1), lwd=0.8, col=1,
xlab="COPULA(u,v) VALUE [JOINT PROBABILITY]",
ylab="KENDALL FUNCTION, AS NONEXCEEDANCE PROBABILITY")
rug(Z, side=1, col=2, lwd=0.2); rug(KF, side=2, col=2, lwd=0.2) # rug plots
lines(zs, kfuncCOP(zs, cop=PLcop, para=1/pi), col=3)
knf_meanZ <- KNFlmr$lambdas[1]; sam_meanZ <- SAMlmr$lambdas[1]
knf_mean <- kfuncCOP(knf_meanZ, cop=PLcop, para=theta) # theo. Kendall function
sam_mean <- kfuncCOP(sam_meanZ, cop=PLcop, para=theta) # sam. est. of Kendall func
points(knf_meanZ, knf_mean, pch=16, col=4, cex=4) # big blue dot
points(sam_meanZ, sam_mean, pch=16, col=5, cex=2) # smaller pale blue dot
lines(zs, zs-zs*log(zs), lty=2, lwd=0.8) # dash ref line for independence
text(0.2, 0.30, knftxt, pos=4, cex=0.8); text(0.2, 0.25, samtxt, pos=4, cex=0.8)
text(0.2, 0.18, paste0("Note the uniform distribution of the ",
"vertical axis rug."), cex=0.7, pos=4) #
# }
Run the code above in your browser using DataLab