PLACKETTpar {copBasic} | R Documentation |
Estimate the Parameter of the Plackett Copula
Description
The parameter \Theta
of the Plackett copula (Nelsen, 2006, pp. 89–92) (PLACKETTcop
or PLcop
) is related to the Spearman Rho (\rho_S \ne 1
, see rhoCOP
)
\rho_S(\Theta) = \frac{\Theta + 1}{\Theta - 1} - \frac{2\Theta\log(\Theta)}{(\Theta - 1)^2}\mbox{.}
Alternatively, the parameter can be estimated using a median-split estimator, which is best shown as an algorithm. First, compute the two medians:
medx <- median(x); medy <- median(y)
Second and third, compute the number of occurrences where both values are less than their medians and express that as a probability:
k <- length(x[x < medx & y < medy]); m <- k / length(x)
Finally, the median-split estimator of \Theta
is computed by
\Theta = \frac{4m^2}{(1-2m)^2}\mbox{.}
Nelsen (2006, p. 92) and Salvadori et al. (2007, p. 247) provide further details. The input values x
and y
are not used for the median-split estimator if Spearman Rho (see rhoCOP
) is provided by rho
.
Usage
PLACKETTpar(x, y, rho=NULL, byrho=FALSE, cor=NULL, ...)
PLpar(x, y, rho=NULL, byrho=FALSE, cor=NULL, ...)
Arguments
x |
Vector of values for random variable |
y |
Vector of values for random variable |
rho |
Spearman Rho and |
byrho |
Should Spearman Rho be used instead of the median-split estimator; |
cor |
A copBasic syntax for “the correlation coefficient” suitable for the copula—a synonym for |
... |
Additional arguments to pass. |
Value
A value for the Plackett copula \Theta
is returned.
Note
Evidently there “does not appear to be a closed form for \tau(\Theta)
” (Fredricks and Nelsen, 2007, p. 2147), but given \rho(\Theta)
, the equivalent \tau(\Theta)
can be computed by either the tauCOP
function or by approximation. One of the Examples sweeps through \rho \mapsto [0,1; \Delta\rho{=}\delta]
, fits the Plackett \theta(\rho)
, and then solves for Kendall Tau \tau(\theta)
using tauCOP
. A polynomial is then fit between \tau
and \rho
to provide rapid conversion between |\rho|
and \tau
, where the residual standard error is 0.0005705, adjusted R-squared is \approx 1
, the maximum residual is \epsilon < 0.006
. Because of symmetry, it is only necessary to fit positive association and reflect the result by the sign of \rho
. This polynomial is from the Examples
is
rho <- 0.920698 "getPLACKETTtau" <- function(rho) { taupoly <- 0.6229945*abs(rho) + 1.1621854*abs(rho)^2 - 10.7424188*abs(rho)^3 + 48.9687845*abs(rho)^4 - 119.0640264*abs(rho)^5 + 160.0438496*abs(rho)^6 - 111.8403591*abs(rho)^7 + 31.8054602*abs(rho)^8 return(sign(rho)*taupoly) } getPLACKETTtau(rho) # 0.7777726
The following code might be useful in some applications for the inversion of the polynomial for the \rho
as a function of \tau
:
"fun" <- function(rho, tau=NULL) {tp <- getPLACKETTtau(rho); return(tau-tp)} tau <- 0.78 rho <- uniroot(fun, interval=c(0, 1), tau=tau)$root # 0.9220636
Author(s)
W.H. Asquith
References
Fredricks, G.A, and Nelsen, R.B., 2007, On the relationship between Spearman's rho and Kendall's tau for pairs of continuous random variables: Journal of Statistical Planning and Inference, v. 137, pp. 2143–2150.
Nelsen, R.B., 2006, An introduction to copulas: New York, Springer, 269 p.
Salvadori, G., De Michele, C., Kottegoda, N.T., and Rosso, R., 2007, Extremes in Nature—An approach using copulas: Springer, 289 p.
See Also
PLACKETTcop
, PLcop
, PLACKETTsim
, rhoCOP
Examples
## Not run:
Q1 <- rnorm(1000); Q2 <- Q1 + rnorm(1000)
PLpar(Q1, Q2); PLpar(Q1, Q2, byrho=TRUE) # two estimates for same data
PLpar(rho= 0.76) # positive association
PLpar(rho=-0.76) # negative association
tauCOP(cop=PLcop, para=PLpar(rho=-0.15, by.rho=TRUE)) #
## End(Not run)
## Not run:
RHOS <- seq(0, 0.990, by=0.002); TAUS <- rep(NA, length(RHOS))
for(i in 1:length(RHOS)) {
#message("Spearman Rho: ", RHOS[i])
theta <- PLACKETTpar(rho=RHOS[i], by.rho=TRUE); tau <- NA
try(tau <- tauCOP(cop=PLACKETTcop, para=theta), silent=TRUE)
TAUS[i] <- ifelse(is.null(tau), NA, tau)
}
LM <- lm(TAUS~ RHOS + I(RHOS^2) + I(RHOS^3) + I(RHOS^4) +
I(RHOS^5) + I(RHOS^6) + I(RHOS^7) + I(RHOS^8) - 1)
plot(RHOS,TAUS, type="l", xlab="abs(Spearman Rho)", ylab="abs(Kendall Tau)")
lines(RHOS,fitted.values(LM), col=3)#
## End(Not run)