,
Christian Genest [aut, ctb],
Alexander J. McNeil [ctb],
Johanna G. Neslehova [ctb]| Title: | Liouville Copulas |
|---|---|
| Description: | Collections of functions allowing random number generations and estimation of 'Liouville' copulas, as described in Belzile and Neslehova (2017) <doi:10.1016/j.jmva.2017.05.008>. |
| Authors: | Leo Belzile [aut, cre] ,
Christian Genest [aut, ctb],
Alexander J. McNeil [ctb],
Johanna G. Neslehova [ctb] |
| Maintainer: | Leo Belzile <[email protected]> |
| License: | GPL-3 |
| Version: | 1.0.7 |
| Built: | 2025-01-05 04:33:03 UTC |
| Source: | https://github.com/lbelzile/lcopula |
lcopula provides functions to evaluate and sample Liouville copulas. The code is adapted from GumbelLiouville.R used in Genest and Neslehova (2013) and from the routines in LiouvilleFunction.R used by McNeil and Neslehova (2010) and coded by Alexander McNeil. The new implementation draws heavily on functions implemented in the copula package and derived from work by Hofert, Maechler and McNeil (2012).
Warning: a different definition is used by the latter and in the copula package for the Clayton copula. The implementation below uses the generator definition found in the first edition of the book Quantitative Risk Management.
The naming conventions are as follows: 'archi' and 'liouv' denote cores, prefixes 'd', 'p', 'r', 's', 'is' denote respectively density, distribution function, random number generation, survival function and inverse survival function. Suffixes are 'm' for marginal, '_m' for simultaneous multiple margins and '_p' for vectorized versions of the functions.
The functions mostly draws from the theory laid out in McNeil and Neslehova (2010) and Belzile (2014), which are the references for definitions and expressions.
The Archimedean families implemented are Clayton, Gumbel, Frank, Ali-Mikhail-Haq (abbreviated AMH) and Joe. Random number generation from the copulae is available via the function rliouv. Method-of-moment estimators can be used for fitting for the Clayton and Gumbel families for bivariate samples, using the function liouv.maxim.mm. More general models can be fitted using liouv.maxim, but optimization may be slow. The function theta.bci allows to construct bootstrap confidence interval for the parameter . Unless the method-of-moment estimator can be used, it is very slow since it relies on full optimization of each replicated sample under the model.
The copula and survival copula domains of attraction are implemented and the corresponding spectral density and Pickands dependence function can be computed and plotted in the bivariate case.
Leo Belzile, using routines from Pr. Alexander J. McNeil, Pr. Johanna G. Neslehova.
The K.plot was adapted from code provided by Pr. Christian Genest.
Maintainer: Leo Belzile <[email protected]>
McNeil A.J. and Neslehova, J.G. (2010) From Archimedean to Liouville Copulas. J. Multivar. Anal., 101(8): 1772–1790.
Belzile L.R. (2014) Extremal and inferential properties of Liouville copulas. Master thesis, McGill.
Hofert, M., Maechler, M., and McNeil, A.J. (2012) Likelihood inference for Archimedean copulas in high dimensions under known margins. J. Multivar. Anal., 110, 133–150.
Genest, C. and Neslehova, J.G. (2013) Assessing and Modeling Asymmetry in Bivariate Continuous Data. In P. Jaworski, F. Durante, and W.K. Hardle (Eds.), Copulae in Mathematical and Quantitative Finance, Lecture Notes in Statistics, 91–114, Springer: Berlin Heidelberg.
Belzile, L.R. and Neslehova, J.G. (2017), Extremal attractors of Liouville copulas, J. Multivar. Anal., 160, 68–92.
copula package
Spectral density of Coles and Tawn extreme value distribution
.ctspecdens(param, dat, transform = TRUE).ctspecdens(param, dat, transform = TRUE)
param |
vector containing the parameters |
dat |
matrix of pseudo-dat, of dimension |
transform |
logical indicating whether parameters are on the log scale. Default to |
the value of the log likelihood for a sample of size n
theta derived from Kendall's The moment estimates are based on inversion of the formula Kendall's for Clayton or Gumbel Liouville copula
.liouv.iTau_s(tau_hat, family, alphavec).liouv.iTau_s(tau_hat, family, alphavec)
tau_hat |
estimated Kendall's |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
Value of theta
The function computes Kendall's for the given model, given alphavec
.liouv.Tau_s(theta, family, alphavec).liouv.Tau_s(theta, family, alphavec)
theta |
parameter of the corresponding Archimedean copula |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
value of
Pickands dependence function as in Belzile (2014), Proposition 40 and Example 4
Returns the Pickands dependence function of the copula domain of attraction
of the survival copula, the scaled Dirichlet extreme value model. Currently only implemented in the bivariate case.
Setting rho=1 yields the same output as the function in the evd package.
.pickands.dir.uni(t, alpha, rho).pickands.dir.uni(t, alpha, rho)
t |
pseudo-angle in (0,1) |
alpha |
vector of Dirichlet allocations. Currently must be of length 2 |
rho |
index of regular variation parameter |
value of Pickands function for the scaled extremal Dirichlet model
Pickands dependence function as in Belzile (2014), Proposition 41 Returns the Pickands dependence function of the copula domain of attraction of the copula. This is only derived and implemented in the bivariate case.
.pickands.fun.uni(t, rho = 0.5, alpha = c(1, 1)).pickands.fun.uni(t, rho = 0.5, alpha = c(1, 1))
t |
pseudo-angle in (0,1) |
rho |
index of regular variation parameter |
alpha |
vector of Dirichlet allocations. Currently must be of length 2 |
value of Pickands function
The function is used internally for optimization.
.pliouv.opt_old(theta, data, family, alphavec, MC.approx = TRUE).pliouv.opt_old(theta, data, family, alphavec, MC.approx = TRUE)
theta |
parameter of the corresponding Archimedean copula |
data |
sample matrix from a Liouville copula |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
MC.approx |
whether to use Monte-Carlo approximation for the inverse survival function (default is |
value of marginal density
The airquality data frame consists of measurements of pollutants in the city of Leeds. In the saying of Boldi and Davison (2007):
"The dataset comprise daily series of monitoring measurements of ozone levels (O3), nitrogen dioxide (NO2), nitrogen oxide (NO) and particulate matter (PM10), in the city centre of Leeds, UK, over 1994–1998. Levels of the gases are measured in parts per billion, and those of PM 10 in micrograms per cubic metre. "
Downloaded from https://uk-air.defra.gov.uk/. (c) Crown 2015 copyright Defra via uk-air.defra.gov.uk, licenced under the Open Government Licence (OGL).
Heffernan, J., Tawn, J. (2004) A conditional approach for multivariate extreme values (with discussion). J. R. Stat. Soc., Ser. B Stat. Methodol. 66(3), 497–546.
Sabourin, A. , Naveau, P., Fougeres, A.-L. (2013) Bayesian Model averaging for Multivariate extremes. Extremes, 16(3), 325–350.
Boldi, M.O., Davison, A.C. (2007) A mixture model for multivariate extremes. J. R. Stat. Soc., Ser. B Stat. Methodol. 69(2), 217–229.
Cooley, D., Davis, R., Naveau, P. (2010) The pairwise beta distribution: A flexible parametric multivariate model for extremes. J. Multivar. Anal. 101(9), 2103–2117.
The danube dataset contains ranks of base flow observations from the Global River Discharge project of the Oak Ridge National Laboratory Distributed Active Archive Center (ORNL DAAC), a NASA data center. The measurements are monthly average flow rate for two stations situated at Scharding (Austria) on the Inn river and at Nagymaros (Hungary) on the Danube.
The data have been pre-processed to remove any time trend. Specifically, Bacigal et al. (2011) extracted the raw data, and obtain the fast Fourier transformed centered observations. The negative spectrum is retained and a linear time series model with 12 seasonal components is fitted. Residuals are then extracted and AR model fitted to the series, the selection being done based on the AIC criterion with imposed maximum order of 3 and the number of autoregressive components may differ for each series.
This data frame contains the following columns:
A numeric vector containing the rank of prewhitened level observations of the Inn river at Nagyramos.
A numeric vector containing the rank of prewhitened level observations of the Donau river at Scharding.
Vorosmarty, C.J., Fekete, B.M., Tucker, B.A. (1998). Global River Discharge, 1807–1991, V. 1.1(RivDIS). Data set. Available from http://doi.org/10.3334/ORNLDAAC/199
Bacigal, T., Jagr, V., Mesiar, R. (2011) Non-exchangeable random variables, Archimax copulas and their fitting to real data. Kybernetika, 47(4), 519–531.
Genest, C. and Neslehova, J. G. (2013) Assessing and Modeling Asymmetry in Bivariate Continuous Data. In P. Jaworski, F. Durante, and W. K. Hardle (Eds.), Copulae in Mathematical and Quantitative Finance, Lecture Notes in Statistics, 91–114, Springer: Berlin Heidelberg.
Spectral density of scaled Dirichlet and negative Dirichlet extreme value distributions
dirspecdens(param, dat, d, transform = TRUE)dirspecdens(param, dat, d, transform = TRUE)
param |
vector containing the parameters |
dat |
matrix of Frechet or Pareto observations, of dimension |
d |
dimension of model |
transform |
logical indicating whether parameters are on the log scale. Default to |
The function is provided as a wrapper and takes parameters on the log scale for ().
the log-likelihood for the n sample
Computes the Liouville EV model or the scaled Dirichlet EV model spectral density
hmvevdliouv(w, alpha, rho, CDA = c("C", "S"), logdensity = FALSE)hmvevdliouv(w, alpha, rho, CDA = c("C", "S"), logdensity = FALSE)
w |
matrix of points in the unit simplex at which to evaluate the density |
alpha |
vector of Dirichlet allocations (strictly positive). |
rho |
parameter of limiting model corresponding to index of regular variation |
CDA |
copula domain of attraction of either the Liouville copula, |
logdensity |
logical; whether to return the log density or not |
a vector with the same number of rows as w.
hmvevdliouv(seq(0.01,0.99,by=0.01), alpha=c(1,2), rho=0.2, CDA="C") hmvevdliouv(seq(0.01,0.99,by=0.01), alpha=c(0.1,2), rho=0.2, CDA="S") hmvevdliouv(seq(0.01,0.99,by=0.01), alpha=c(1,2), rho=0.2, CDA="S")hmvevdliouv(seq(0.01,0.99,by=0.01), alpha=c(1,2), rho=0.2, CDA="C") hmvevdliouv(seq(0.01,0.99,by=0.01), alpha=c(0.1,2), rho=0.2, CDA="S") hmvevdliouv(seq(0.01,0.99,by=0.01), alpha=c(1,2), rho=0.2, CDA="S")
This function is a wrapper around isliouvm; it allows the user to treat all the data matrix simultaneously
by applying different parameters to each margin.
isliouv_m(u, family, alphavec, theta)isliouv_m(u, family, alphavec, theta)
u |
vector of survival probabilities |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
theta |
parameter of the corresponding Archimedean copula |
a vector of same length as u with the quantile at 1-u
u <- rliouv(n = 10, family = "clayton", alphavec <- c(2,3), theta = 2) isliouv_m(u=u, family="clayton", alphavec=c(2,3), theta=2)u <- rliouv(n = 10, family = "clayton", alphavec <- c(2,3), theta = 2) isliouv_m(u=u, family="clayton", alphavec=c(2,3), theta=2)
This function plots the expectation of the order statistics under the null hypothesis of independence against the ordered empirical copula values. The data is transformed to ranks.
K.plot(data, add = F, ...)K.plot(data, add = F, ...)
data |
a |
add |
whether to superimpose lines to an existing graph. Default to |
... |
additional arguments passed to |
The function uses integrate and may fail for large d or large n. If , the fallback is to generate
a corresponding sample of uniform variates and to compare the empirical copula of the sample generated under the null hypothesis with the one
obtained from the sample.
The Kendall plot corresponding to the data at hand
Pr. Christian Genest (the code was adapted for the multivariate case)
Genest & Boies (2003). Detecting Dependence with Kendall Plots, The American Statistician, 57(4), 275–284.
#Independence K.plot(matrix(runif(2000),ncol=2)) #Negative dependence K.plot(rCopula(n=1000,claytonCopula(param=-0.5,dim=2)),add=TRUE,col=2) #Perfect negative dependence K.plot(rCopula(n=1000,claytonCopula(param=-1,dim=2)),add=TRUE,col=6) #Positive dependence K.plot(rCopula(n=1000,claytonCopula(param=iTau(claytonCopula(0.3),0.5),dim=2)),add=TRUE,col=3) #Perfect positive dependence K.plot(rCopula(n=1000,claytonCopula(param=iTau(claytonCopula(0.3),1),dim=2)),add=TRUE,col=4)#Independence K.plot(matrix(runif(2000),ncol=2)) #Negative dependence K.plot(rCopula(n=1000,claytonCopula(param=-0.5,dim=2)),add=TRUE,col=2) #Perfect negative dependence K.plot(rCopula(n=1000,claytonCopula(param=-1,dim=2)),add=TRUE,col=6) #Positive dependence K.plot(rCopula(n=1000,claytonCopula(param=iTau(claytonCopula(0.3),0.5),dim=2)),add=TRUE,col=3) #Perfect positive dependence K.plot(rCopula(n=1000,claytonCopula(param=iTau(claytonCopula(0.3),1),dim=2)),add=TRUE,col=4)
The moment estimates are based on inversion of the formula Kendall's tau for Clayton or Gumbel Liouville copula
liouv.iTau(tau_hat, family, alphavec)liouv.iTau(tau_hat, family, alphavec)
tau_hat |
estimated vector of Kendall's tau values |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
Vector of theta
liouv.iTau(0.5,family="gumbel", c(1,2)) liouv.iTau(0.5,family="clayton", c(3,2))liouv.iTau(0.5,family="gumbel", c(1,2)) liouv.iTau(0.5,family="clayton", c(3,2))
Two methods, either numerical optimization or method-of-moments
liouv.maxim( data, family, interval, boundary = NULL, lattice.mat = NULL, return_all = FALSE, MC.approx = TRUE )liouv.maxim( data, family, interval, boundary = NULL, lattice.mat = NULL, return_all = FALSE, MC.approx = TRUE )
data |
sample matrix from a Liouville copula |
family |
family of the Liouville copula. Either |
interval |
interval over which to look for |
boundary |
vector of endpoints for search of Dirichlet allocation parameters. Either |
lattice.mat |
matrix of tuples of Dirichlet allocation parameters at which to evaluate the likelihood |
return_all |
should all results (as list) or only maximum value be returned. Defaults to |
MC.approx |
whether to use Monte-Carlo approximation for the inverse survival function (default is |
A wrapper to optim using the Nelder-Mead algorithm or using the methods of moments,
to maxime pointwise given every alphavec over a grid.
Returns the maximum for alphavec and theta.
a list with values of theta and Dirichlet parameter along with maximum found. Gives index of maximum amongst models fitted.
## Not run: data <- rliouv(n=100, family="joe", alphavec=c(1,2), theta=2) liouv.maxim(data=data, family="j", interval=c(1.25,3), boundary=c(2,2),return_all=TRUE) lattice.mat <- t(combn(1:3,2)) liouv.maxim(data=data, family="j", interval=c(1.25,3), lattice.mat=lattice.mat, return_all=FALSE) #data <- rliouv(n=1000, family="gumbel", alphavec=c(1,2), theta=2) liouv.maxim.mm(data=data, family="gumbel", boundary=c(3,3),return_all=TRUE) lattice.mat <- t(combn(1:3,2)) liouv.maxim.mm(data=data, family="gumbel", lattice.mat=lattice.mat, return_all=FALSE) ## End(Not run)## Not run: data <- rliouv(n=100, family="joe", alphavec=c(1,2), theta=2) liouv.maxim(data=data, family="j", interval=c(1.25,3), boundary=c(2,2),return_all=TRUE) lattice.mat <- t(combn(1:3,2)) liouv.maxim(data=data, family="j", interval=c(1.25,3), lattice.mat=lattice.mat, return_all=FALSE) #data <- rliouv(n=1000, family="gumbel", alphavec=c(1,2), theta=2) liouv.maxim.mm(data=data, family="gumbel", boundary=c(3,3),return_all=TRUE) lattice.mat <- t(combn(1:3,2)) liouv.maxim.mm(data=data, family="gumbel", lattice.mat=lattice.mat, return_all=FALSE) ## End(Not run)
The function computes Kendall's for the given model, given alphavec
liouv.Tau(theta, family, alphavec)liouv.Tau(theta, family, alphavec)
theta |
parameter of the corresponding Archimedean copula |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
vector of
liouv.Tau(theta=2, family="gumbel", alphavec=c(1,2)) liouv.Tau(theta=1, family="clayton", alphavec=c(2,1))liouv.Tau(theta=2, family="gumbel", alphavec=c(1,2)) liouv.Tau(theta=1, family="clayton", alphavec=c(2,1))
Multivariate density, survival copula and random generation for the Liouville copulas.
rliouv(n = 100, family, alphavec, theta, reverse = FALSE) pliouv(x, theta, family, alphavec) dliouv(x, family, alphavec, theta, is.log = FALSE)rliouv(n = 100, family, alphavec, theta, reverse = FALSE) pliouv(x, theta, family, alphavec) dliouv(x, family, alphavec, theta, is.log = FALSE)
n |
sample size |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) Specifies (implictly) the dimension of sample |
theta |
parameter of the corresponding Archimedean copula |
reverse |
if |
x |
matrix of quantiles from a Liouville copula |
is.log |
if |
rliouv generates draws from the Liouville copula. dliouv evaluates the density of an n by d matrix of observations. pliouv is the (survival) copula associated with the Liouville vector and is as such the multivariate distribution function for uniform observations.
Liouville copulas were introduced in McNeil and Neslehova (2010), generalizing Archimedean copulas. Like the latter, they
are survival copulas, which means that the copula is evaluated using the (multivariate) survival function of Liouville vectors. See also sliouv for the latter.
The Liouville copula is by definition a survival copula. The function thus
maps marginally observations from the unit interval to the positive half-line using the marginal inverse
survival function isliouvm of the Liouville vector, and then evaluating the survival
distribution at the resulting Liouville vector.
either a matrix of dimension n by length(alphavec) with the corresponding quantile, probability, survival probability or sample from the Liouville vector
McNeil A.J. and Neslehova, J.G. (2010) From Archimedean to Liouville Copulas. J. Multivar. Anal., 101(8): 1772–1790.
## Not run: #Multivariate density of Clayton Liouville copula x <- rliouv(n = 100, family = "clayton", alphavec <- c(2,3), theta = 2) dliouv(x=x, family="clayton", alphavec=c(2,3), theta=2, TRUE) #Distribution function, multivariate sample x <- rliouv(n=100, family="frank", theta=1.5, alphavec=c(2,3)) pliouv(theta=1.5, x=x,family="frank", alphavec=c(2,3)) ## End(Not run)## Not run: #Multivariate density of Clayton Liouville copula x <- rliouv(n = 100, family = "clayton", alphavec <- c(2,3), theta = 2) dliouv(x=x, family="clayton", alphavec=c(2,3), theta=2, TRUE) #Distribution function, multivariate sample x <- rliouv(n=100, family="frank", theta=1.5, alphavec=c(2,3)) pliouv(theta=1.5, x=x,family="frank", alphavec=c(2,3)) ## End(Not run)
Marginal density, distribution, survival and inverse survival functions for Liouville copulas or Liouville vectors.
The inverse survival function of Liouville vectors is not available in closed-form and is obtained numerically by root-finding.
As such, Monte-Carlo approximation have been considered for dealing with inference to avoid computational bottlenecks.
Note: the arguments of sliouv are reversed since they are meant to be called inside optim. The functions borrow
functions and their derivatives from the copula-package.
sliouvm(x, family, alpha, theta) pliouvm(x, family, alpha, theta) isliouvm(u, family, alpha, theta) dliouvm(x, family, alpha, theta)sliouvm(x, family, alpha, theta) pliouvm(x, family, alpha, theta) isliouvm(u, family, alpha, theta) dliouvm(x, family, alpha, theta)
x |
vector of quantiles from a Liouville copula (or a Liouville vector for the survival function , with support on the positive real line) |
family |
family of the Liouville copula. Either |
alpha |
integer Dirichlet parameter |
theta |
parameter of the corresponding Archimedean copula |
u |
vector of quantiles or survival probabilities, (pseudo)-uniform variates |
a vector with the corresponding quantile, probability, survival probabilities
## Not run: #Marginal density samp <- rliouv(n = 100, family = "clayton", alphavec <- c(2,3), theta = 2) dliouvm(x=samp[,1], family="clayton", alpha=2, theta=2) sum(log(dliouvm(x=samp[,1], family="clayton", alpha=2, theta=2))) #Marginal distribution and (inverse) survival function x <- rliouv(n = 100, family = "gumbel", alphavec <- c(2,3), theta = 2) pliouvm(x[,1], family="gumbel", alpha=alphavec[1], theta=2) su <- sliouvm(1-x[,1], family="gumbel", alpha=alphavec[1], theta=2) isliouvm(u=su, family="clayton", alpha=2, theta=2) #pliouv is the same as sliouv(isliouvm) ## End(Not run)## Not run: #Marginal density samp <- rliouv(n = 100, family = "clayton", alphavec <- c(2,3), theta = 2) dliouvm(x=samp[,1], family="clayton", alpha=2, theta=2) sum(log(dliouvm(x=samp[,1], family="clayton", alpha=2, theta=2))) #Marginal distribution and (inverse) survival function x <- rliouv(n = 100, family = "gumbel", alphavec <- c(2,3), theta = 2) pliouvm(x[,1], family="gumbel", alpha=alphavec[1], theta=2) su <- sliouvm(1-x[,1], family="gumbel", alpha=alphavec[1], theta=2) isliouvm(u=su, family="clayton", alpha=2, theta=2) #pliouv is the same as sliouv(isliouvm) ## End(Not run)
Spectral density of scaled negative Dirichlet extreme value distribution
negdirspecdens(param, dat, d, transform = TRUE)negdirspecdens(param, dat, d, transform = TRUE)
param |
vector containing the parameters |
dat |
matrix of Frechet or Pareto observations, of dimension |
d |
dimension of model |
transform |
logical indicating whether parameters are on the log scale. Default to |
The function is provided as a wrapper and takes parameters on the log scale for and .
the log-likelihood for the n sample
The nutrient data frame consists of quintuples consisting of four day measurements for intake of calcium, iron, protein, vitamin A and C from women aged 25 to 50 in the United States as part of the "Continuing Survey of Food Intakes of Individuals" program. The processed data has
737 measurements from a cohort study of the United States Department of Agriculture (USDA) and it is available online from the University of Pennsylvania repository.
This data frame contains the following columns:
A numeric vector containing the identification of the participant.
A numeric vector containing measured calcium.
A numeric vector containing iron measurements.
A numeric vector containing protein measurements.
A numeric vector containing vitamin A measurements.
A numeric vector containing vitamin C measurements.
The survey data was processed by Dr. Andrew Wiesner for a course on multivariate statistics at The Pennsylvania State University.
Genest, C., Neslehova, J., Quessy, J.-F. (2012) Tests of symmetry for bivariate copulas. Annals of the Institute of Statistical Mathematics, 64(4), 811–834.
Genest, C. and Neslehova, J. G. (2013) Assessing and Modeling Asymmetry in Bivariate Continuous Data. In P. Jaworski, F. Durante, & W. K. Hardle (Eds.), Copulae in Mathematical and Quantitative Finance, Lecture Notes in Statistics, 91–114, Springer: Berlin Heidelberg.
Li, B., Genton, M. G. (2013) Nonparametric identification of copula structures. Journal of the American Statistical Association, 108(502), 666–675.
Pickands dependence function as in Belzile (2014), Proposition 40 and Example 4 and Belzile (2014), Proposition 41, assuming that the parameter alpha is integer-valued. Returns the Pickands dependence function of the copula domain of attraction (CDA) of the survival copula, the scaled Dirichlet extreme value model, or the CDA of the copula, the Liouville EV model.
pickands.liouv(t, rho = 0.5, alpha = c(1, 1), CDA = c("C", "S"))pickands.liouv(t, rho = 0.5, alpha = c(1, 1), CDA = c("C", "S"))
t |
pseudo-angle in (0,1) |
rho |
index of regular variation parameter |
alpha |
vector of Dirichlet allocations. Currently must be of length 2 |
CDA |
select the extremal attractor of the copula ( |
value of Pickands function for the scaled Dirichlet EV model
pickands.liouv(seq(0,1,by=0.01),1,c(0.1,0.3),CDA="S") pickands.liouv(t = seq(0,1,by=0.01), rho = 0.5, alpha = c(1,3), CDA="C")pickands.liouv(seq(0,1,by=0.01),1,c(0.1,0.3),CDA="S") pickands.liouv(t = seq(0,1,by=0.01), rho = 0.5, alpha = c(1,3), CDA="C")
The function will draw the Pickands dependence function for output in tikz if
the corresponding function is selected.
pickands.plot( rho, alpha, plot.new = TRUE, CDA = c("C", "S"), tikz = FALSE, ... )pickands.plot( rho, alpha, plot.new = TRUE, CDA = c("C", "S"), tikz = FALSE, ... )
rho |
index of regular variation parameter |
alpha |
vector of Dirichlet allocations. Currently must be of length 2 |
plot.new |
boolean indicating whether a new plotting device should be called |
CDA |
whether to plot Pickands function for the extremal model of the
copula ( |
tikz |
boolean specifying whether to prepare plot for |
... |
additional arguments passed to |
a plot of the Pickands dependence function
pickands.plot( rho = 0.9, alpha = c(1,1), col = "slateblue1", CDA = "C") pickands.plot( rho = 0.9, alpha = c(2,3), col = "slateblue2", CDA = "C", plot.new = FALSE) pickands.plot( rho = 0.5, alpha = c(2,3), col = "slateblue3", CDA = "C", plot.new = FALSE) #Parameters for the Pickands function of the scaled Dirichlet need not be integer pickands.plot( rho = 0.9, alpha = c(1,1), CDA = "S") pickands.plot( rho = 0.9, alpha = c(0.2,0.5), col = "darkred", CDA = "S", plot.new = FALSE) pickands.plot( rho = 0.8, alpha = c(1.2,0.1), col = "red", CDA = "S", plot.new = FALSE)pickands.plot( rho = 0.9, alpha = c(1,1), col = "slateblue1", CDA = "C") pickands.plot( rho = 0.9, alpha = c(2,3), col = "slateblue2", CDA = "C", plot.new = FALSE) pickands.plot( rho = 0.5, alpha = c(2,3), col = "slateblue3", CDA = "C", plot.new = FALSE) #Parameters for the Pickands function of the scaled Dirichlet need not be integer pickands.plot( rho = 0.9, alpha = c(1,1), CDA = "S") pickands.plot( rho = 0.9, alpha = c(0.2,0.5), col = "darkred", CDA = "S", plot.new = FALSE) pickands.plot( rho = 0.8, alpha = c(1.2,0.1), col = "red", CDA = "S", plot.new = FALSE)
The function is used internally for optimization.
pliouv.opt(theta, data, family, alphavec, MC.approx = TRUE)pliouv.opt(theta, data, family, alphavec, MC.approx = TRUE)
theta |
parameter of the corresponding Archimedean copula |
data |
sample matrix from a Liouville copula |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
MC.approx |
whether to use Monte-Carlo approximation for the inverse survival function (default is |
value of marginal density
Sampler based on the Marshall-Olkin algorithm
rarchi(n, family, d, theta)rarchi(n, family, d, theta)
n |
sample size |
family |
family of the Archimedean copula. Either |
d |
dimension of sample |
theta |
parameter of the Archimedean copula |
a sample of dimension n by d from the Archimedean copula
#Sample from a Gumbel Archimedean copula rarchi(n = 100, "gumbel", d = 4, theta = 2) #Sample from the independence copula rarchi(n = 100, "gumbel", d = 4, theta = 1)#Sample from a Gumbel Archimedean copula rarchi(n = 100, "gumbel", d = 4, theta = 2) #Sample from the independence copula rarchi(n = 100, "gumbel", d = 4, theta = 1)
sliouv returns the survival function of a Liouville vector. For the survival copula
and the associated probability, see pliouv.
sliouv(theta, x, family, alphavec)sliouv(theta, x, family, alphavec)
theta |
parameter of the corresponding Archimedean copula |
x |
an |
family |
family of the Liouville copula. Either |
alphavec |
|
This function is a wrapper around sliouv; it allows the user to treat all the data matrix simultaneously
by applying different parameters to each margin.
sliouv_m(x, family, alphavec, theta)sliouv_m(x, family, alphavec, theta)
x |
sample from copula |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
theta |
parameter of the corresponding Archimedean copula |
a matrix of same length as x with the survival probabilities
x <- rliouv(n = 100, family = "gumbel", alphavec <- c(2,3), theta = 2) sliouv_m(x, family="gumbel", alphavec=c(2,3), theta=2) all(sliouv_m(x, family="gumbel", alphavec=c(2,3), theta=2)[,1]- sliouvm(x[,1], family="gumbel", alpha=2, theta=2)==0)x <- rliouv(n = 100, family = "gumbel", alphavec <- c(2,3), theta = 2) sliouv_m(x, family="gumbel", alphavec=c(2,3), theta=2) all(sliouv_m(x, family="gumbel", alphavec=c(2,3), theta=2)[,1]- sliouvm(x[,1], family="gumbel", alpha=2, theta=2)==0)
theta for Liouville copulaThe parametric bootstrap provides confidence intervals by repeatedly sampling datasets from the postulated
Liouvilla copula model. If and the model is either gumbel or clayton, the value of
Kendall's is calculated from the sample, and the confidence interval or the quantiles correspond
to the inverse for the bootstrap quantile values of (using monotonicity).
theta.bci( B = 1999, family, alphavec, n, theta.hat, quant = c(0.025, 0.975), silent = FALSE )theta.bci( B = 1999, family, alphavec, n, theta.hat, quant = c(0.025, 0.975), silent = FALSE )
B |
number of bootstrap replicates |
family |
family of the Liouville copula. Either |
alphavec |
vector of Dirichlet allocations (must be a vector of integers) |
n |
sample size |
theta.hat |
estimate of theta |
quant |
if the vector of probability is specified, the function will return the corresponding bootstrap quantiles |
silent |
boolean for output progress. Default is |
Install package wdm to speed up calculation of Kendall's tau.
Since no closed-form formulas exist for the other models or in higher dimension, the method is extremely slow since it relies on maximization of a new sample from the model and look up the corresponding parameters.
a list with a 95
and the bootstrap values of Kendall's tau in boot_tau if and the model is either gumbel or clayton.
Otherwise, the list contains boot_theta.
## Not run: theta.bci(B=99, family="gumbel", alphavec=c(2,3), n=100, theta.hat=2) theta.bci(B=19, family="AMH", alphavec=c(1,2), n=100, theta.hat=0.5, quant=c(0.05,0.95)) theta.bci(B=19, family="frank", alphavec=c(1,2,3), n=100, theta.hat=0.5, quant=c(0.05,0.95)) ## End(Not run)## Not run: theta.bci(B=99, family="gumbel", alphavec=c(2,3), n=100, theta.hat=2) theta.bci(B=19, family="AMH", alphavec=c(1,2), n=100, theta.hat=0.5, quant=c(0.05,0.95)) theta.bci(B=19, family="frank", alphavec=c(1,2,3), n=100, theta.hat=0.5, quant=c(0.05,0.95)) ## End(Not run)