| Title: | High Performance Algorithms for Vine Copula Modeling |
|---|---|
| Description: | Provides an interface to 'vinecopulib', a C++ library for vine copula modeling. The 'rvinecopulib' package implements the core features of the popular 'VineCopula' package, in particular inference algorithms for both vine copula and bivariate copula models. Advantages over 'VineCopula' are a sleeker and more modern API, improved performances, especially in high dimensions, nonparametric and multi-parameter families, and the ability to model discrete variables. The 'rvinecopulib' package includes 'vinecopulib' as header-only C++ library (currently version 1.0.0). Thus users do not need to install 'vinecopulib' itself in order to use 'rvinecopulib'. Since their initial releases, 'vinecopulib' is licensed under the MIT License, and 'rvinecopulib' is licensed under the GNU GPL version 3. |
| Authors: | Thomas Nagler [aut, cre], Thibault Vatter [aut] |
| Maintainer: | Thomas Nagler <[email protected]> |
| License: | GPL-3 | file LICENSE |
| Version: | 1.0.0.1.0 |
| Built: | 2026-09-18 16:40:12 UTC |
| Source: | https://github.com/vinecopulib/rvinecopulib |
as_margin() validates fitted margins and adapts the legacy
list(distr = ...) representation through stats_margin().
as_margin(margin)as_margin(margin)
margin |
a fitted margin or legacy stats distribution list. |
A validated fitted margin.
as_rvine_structure and as_rvine_matrix are new S3 generics allowing to
coerce objects into R-vine structures and matrices (see rvine_structure()
and rvine_matrix()).
as_rvine_structure(x, ...) as_rvine_matrix(x, ...) ## S3 method for class 'rvine_structure' as_rvine_structure(x, ..., validate = FALSE) ## S3 method for class 'rvine_structure' as_rvine_matrix(x, ..., validate = FALSE) ## S3 method for class 'list' as_rvine_structure(x, ..., is_natural_order = FALSE) ## S3 method for class 'list' as_rvine_matrix(x, ..., is_natural_order = FALSE) ## S3 method for class 'rvine_matrix' as_rvine_structure(x, ..., validate = FALSE) ## S3 method for class 'rvine_matrix' as_rvine_matrix(x, ..., validate = FALSE) ## S3 method for class 'matrix' as_rvine_structure(x, ..., validate = TRUE) ## S3 method for class 'matrix' as_rvine_matrix(x, ..., validate = TRUE)as_rvine_structure(x, ...) as_rvine_matrix(x, ...) ## S3 method for class 'rvine_structure' as_rvine_structure(x, ..., validate = FALSE) ## S3 method for class 'rvine_structure' as_rvine_matrix(x, ..., validate = FALSE) ## S3 method for class 'list' as_rvine_structure(x, ..., is_natural_order = FALSE) ## S3 method for class 'list' as_rvine_matrix(x, ..., is_natural_order = FALSE) ## S3 method for class 'rvine_matrix' as_rvine_structure(x, ..., validate = FALSE) ## S3 method for class 'rvine_matrix' as_rvine_matrix(x, ..., validate = FALSE) ## S3 method for class 'matrix' as_rvine_structure(x, ..., validate = TRUE) ## S3 method for class 'matrix' as_rvine_matrix(x, ..., validate = TRUE)
x |
An object of class |
... |
Other arguments passed on to individual methods. |
validate |
When 'TRUE“, verifies that the input is a valid rvine-structure (see Details). You may want to suppress this when you know that you already have a valid structure and you want to save some time, or to explicitly enable it if you have a structure that you want to re-check. |
is_natural_order |
A flag indicating whether the |
The coercion to rvine_structure and rvine_matrix can be applied to
different kind of objects Currently, rvine_structure, rvine_matrix,
matrix and list are supported.
For as_rvine_structure:
rvine_structure : the main use case is to re-check an object via
validate = TRUE.
rvine_matrix and matrix : allow to coerce matrices into R-vine
structures (see rvine_structure() for more details). The main difference
between rvine_matrix and matrix is the nature of the validity
checks.
list : must contain named elements order and struct_array to be
coerced into an R-vine structure (see rvine_structure() for more details).
For as_rvine_matrix:
rvine_structure : allow to coerce an rvine_structure into an
R-vine matrix (useful e.g. for printing).
rvine_matrix: similar to as_rvine_structure for rvine_structure,
the main use case is to re-check an object via validate = TRUE.
matrix : allow to coerce matrices into R-vine
matrices (mainly by checking that the matrix defines a valid
R-vine, see rvine_matrix() for more details).
list : must contain named elements order and struct_array to be
coerced into an R-vine matrix (see rvine_structure() for more details).
Either an object of class rvine_structure or of class
rvine_matrix (see rvine_structure() or rvine_matrix()).
rvine_structure rvine_matrix
# R-vine structures can be constructed from the order vector and struct_array rvine_structure(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3), 2 )) # ... or a similar list can be coerced into an R-vine structure as_rvine_structure(list(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3), 2 ))) # similarly, standard matrices can be coerced into R-vine structures mat <- matrix(c(4, 3, 2, 1, 4, 3, 2, 0, 4, 3, 0, 0, 4, 0, 0, 0), 4, 4) as_rvine_structure(mat) # or truncate and construct the structure mat[3, 1] <- 0 as_rvine_structure(mat) # throws an error mat[3, 1] <- 5 try(as_rvine_structure(mat))# R-vine structures can be constructed from the order vector and struct_array rvine_structure(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3), 2 )) # ... or a similar list can be coerced into an R-vine structure as_rvine_structure(list(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3), 2 ))) # similarly, standard matrices can be coerced into R-vine structures mat <- matrix(c(4, 3, 2, 1, 4, 3, 2, 0, 4, 3, 0, 0, 4, 0, 0, 0), 4, 4) as_rvine_structure(mat) # or truncate and construct the structure mat[3, 1] <- 0 as_rvine_structure(mat) # throws an error mat[3, 1] <- 5 try(as_rvine_structure(mat))
Convert list to bicop object
as.bicop(object, check = TRUE)as.bicop(object, check = TRUE)
object |
a list containing entries for |
check |
whether to check for validity of the family/parameter specification. |
A bicop object corresponding to the specification in object.
as.bicop(list(family = "gumbel", rotation = 90, parameters = 2, npars = 1))as.bicop(list(family = "gumbel", rotation = 90, parameters = 2, npars = 1))
Fit a bivariate copula model for continuous or discrete data. The family can be selected automatically from a vector of options.
bicop( data, var_types = c("c", "c"), family_set = "all", par_method = "mle", nonpar_method = "quadratic", mult = 1, selcrit = "aic", weights = numeric(), psi0 = 0.9, presel = TRUE, allow_rotations = TRUE, keep_data = FALSE, cores = 1 )bicop( data, var_types = c("c", "c"), family_set = "all", par_method = "mle", nonpar_method = "quadratic", mult = 1, selcrit = "aic", weights = numeric(), psi0 = 0.9, presel = TRUE, allow_rotations = TRUE, keep_data = FALSE, cores = 1 )
data |
a matrix or data.frame with at least two columns, containing the (pseudo-)observations for the two variables (copula data should have approximately uniform margins). More columns are required for discrete models, see Details. |
var_types |
variable types, a length 2 vector; e.g., |
family_set |
a character vector of families; see Details for additional options. |
par_method |
the estimation method for parametric models, either |
nonpar_method |
the estimation method for nonparametric models, either
|
mult |
multiplier for the smoothing parameters of nonparametric families. Values larger than 1 make the estimate more smooth, values less than 1 less smooth. |
selcrit |
criterion for family selection, either |
weights |
optional vector of weights for each observation. |
psi0 |
see |
presel |
whether the family set should be thinned out according to symmetry characteristics of the data. |
allow_rotations |
whether to allow rotations of the copula. |
keep_data |
whether the data should be stored (necessary for using
|
cores |
number of cores to use; if more than 1, estimation for multiple families is done in parallel. |
If there are missing data (i.e., NA entries), incomplete observations are
discarded before fitting the copula.
When at least one variable is discrete, more than two columns are required
for data: the first block contains realizations of
. The second block contains
realizations of . The minus indicates a
left-sided limit of the cdf. For, e.g., an integer-valued variable, it holds
. For continuous variables the left
limit and the cdf itself coincide. Respective columns can be omitted in the
second block.
The family_set argument accepts all families in bicop_dist() plus the
following convenience definitions:
"all" contains all the families,
"parametric" contains the parametric families (all except "tll"),
"nonparametric" contains the nonparametric families ("indep" and
"tll")
"onepar" contains the parametric families with a single parameter,
("gaussian", "clayton", "gumbel", "frank", and "joe"),
"twopar" contains the parametric families with two parameters,
("t", "bb1", "bb6", "bb7", and "bb8"),
"threepar" contains the parametric families with three parameters,
("tawn"),
"elliptical" contains the elliptical families,
"archimedean" contains the archimedean families,
"ev" contains the extreme-value families,
"BB" contains the BB families,
"itau" families for which estimation by Kendall's tau inversion is
available ("indep","gaussian", "t","clayton", "gumbel",
"frank", "joe").
An object inheriting from classes bicop and bicop_dist . In addition to
the entries contained in bicop_dist(), objects from the bicop class
contain:
data (optionally, if keep_data = TRUE was used), the dataset that was
passed to bicop().
controls, a list with the set of fit controls that was passed to
bicop().
loglik the log-likelihood.
nobs, an integer with the number of observations that was used to fit
the model.
bicop_dist(), plot.bicop(), contour.bicop(), dbicop(),
pbicop(), hbicop(), rbicop()
## fitting a continuous model from simulated data u <- rbicop(100, "clayton", 90, 3) fit <- bicop(u, family_set = "par") summary(fit) ## compare fit with true model contour(fit) contour(bicop_dist("clayton", 90, 3), col = 2, add = TRUE) ## fit a model from discrete data x_disc <- qpois(u, 1) # transform to Poisson margins plot(x_disc) udisc <- cbind(ppois(x_disc, 1), ppois(x_disc - 1, 1)) fit_disc <- bicop(udisc, var_types = c("d", "d")) summary(fit_disc)## fitting a continuous model from simulated data u <- rbicop(100, "clayton", 90, 3) fit <- bicop(u, family_set = "par") summary(fit) ## compare fit with true model contour(fit) contour(bicop_dist("clayton", 90, 3), col = 2, add = TRUE) ## fit a model from discrete data x_disc <- qpois(u, 1) # transform to Poisson margins plot(x_disc) udisc <- cbind(ppois(x_disc, 1), ppois(x_disc - 1, 1)) fit_disc <- bicop(udisc, var_types = c("d", "d")) summary(fit_disc)
Computes the four corner tail-dependence coefficients or Blomqvist's beta for a bivariate copula distribution.
tail_dep(object) blomqvist_beta(object)tail_dep(object) blomqvist_beta(object)
object |
a bicop_dist or fitted bicop object. |
tail_dep() returns a 2 by 2 matrix whose rows refer to the lower
and upper tail of the first variable and whose columns refer to the lower
and upper tail of the second variable. blomqvist_beta() returns a
numeric scalar.
cop <- bicop_dist("clayton", 0, 2) tail_dep(cop) blomqvist_beta(cop)cop <- bicop_dist("clayton", 0, 2) tail_dep(cop) blomqvist_beta(cop)
Create custom bivariate copula models by specifying the family, rotation, parameters, and variable types.
bicop_dist( family = "indep", rotation = 0, parameters = numeric(0), var_types = c("c", "c") )bicop_dist( family = "indep", rotation = 0, parameters = numeric(0), var_types = c("c", "c") )
family |
the copula family, a string containing the family name (see Details for all possible families). |
rotation |
the rotation of the copula, one of |
parameters |
a vector or matrix of copula parameters. |
var_types |
variable types, a length 2 vector; e.g., |
| type | name | name in R |
| - | Independence | "indep" |
| Elliptical | Gaussian | "gaussian" |
| " | Student t | "t" |
| Archimedean | Clayton | "clayton" |
| " | Gumbel | "gumbel" |
| " | Frank | "frank" |
| " | Joe | "joe" |
| " | Clayton-Gumbel (BB1) | "bb1" |
| " | Joe-Gumbel (BB6) | "bb6" |
| " | Joe-Clayton (BB7) | "bb7" |
| " | Joe-Frank (BB8) | "bb8" |
| Extreme-value | Tawn | "tawn" |
| Nonparametric | Transformation kernel | "tll" |
An object of class bicop_dist, i.e., a list containing:
family, a character indicating the copula family.
rotation, an integer indicating the rotation (i.e., either 0, 90, 180,
or 270).
parameters, a numeric vector or matrix of parameters.
npars, a numeric with the (effective) number of parameters.
var_types, the variable types.
bicop_dist(), plot.bicop(), contour.bicop(), dbicop(),
pbicop(), hbicop(), rbicop()
## Clayton 90° copula with parameter 3 cop <- bicop_dist("clayton", 90, 3) cop str(cop) ## visualization plot(cop) contour(cop) plot(rbicop(200, cop)) ## BB8 copula model for discrete data cop_disc <- bicop_dist("bb8", 0, c(2, 0.5), var_types = c("d", "d")) cop_disc## Clayton 90° copula with parameter 3 cop <- bicop_dist("clayton", 90, 3) cop str(cop) ## visualization plot(cop) contour(cop) plot(rbicop(200, cop)) ## BB8 copula model for discrete data cop_disc <- bicop_dist("bb8", 0, c(2, 0.5), var_types = c("d", "d")) cop_disc
Density, distribution function, random generation and h-functions (with their inverses) for the bivariate copula distribution.
dbicop( u, family, rotation, parameters, var_types = c("c", "c"), log = FALSE, deriv = NULL, cores = 1 ) pbicop(u, family, rotation, parameters, var_types = c("c", "c")) rbicop(n, family, rotation, parameters, qrng = FALSE) ## S3 method for class 'bicop_dist' scores(u, vinecop, cores = 1, parameters = NULL, ...) ## S3 method for class 'bicop_dist' hessian(u, vinecop, cores = 1, parameters = NULL, ...) hbicop( u, cond_var, family, rotation, parameters, inverse = FALSE, var_types = c("c", "c"), deriv = NULL, cores = 1 )dbicop( u, family, rotation, parameters, var_types = c("c", "c"), log = FALSE, deriv = NULL, cores = 1 ) pbicop(u, family, rotation, parameters, var_types = c("c", "c")) rbicop(n, family, rotation, parameters, qrng = FALSE) ## S3 method for class 'bicop_dist' scores(u, vinecop, cores = 1, parameters = NULL, ...) ## S3 method for class 'bicop_dist' hessian(u, vinecop, cores = 1, parameters = NULL, ...) hbicop( u, cond_var, family, rotation, parameters, inverse = FALSE, var_types = c("c", "c"), deriv = NULL, cores = 1 )
u |
evaluation points, a matrix with at least two columns, see Details. |
family |
the copula family, a string containing the family name (see
|
rotation |
the rotation of the copula, one of |
parameters |
a vector or matrix of copula parameters. For |
var_types |
variable types, a length 2 vector; e.g., |
log |
whether to return the log-density or a derivative of the log-density. |
deriv |
|
cores |
number of cores used when evaluating observation-specific parameters. |
n |
number of observations. If 'length(n) > 1“, the length is taken to be the number required. |
qrng |
if |
vinecop |
a |
... |
unused. |
cond_var |
either |
inverse |
whether to compute the h-function or its inverse. |
See bicop for the various implemented copula families.
The copula density is defined as joint density divided by marginal densities, irrespective of variable types.
H-functions (hbicop()) are conditional distributions derived
from a copula. If is a copula, then
In other words, the H-function number refers to the conditioning variable.
When inverting H-functions, the inverse is then taken with respect to the
other variable, that is v when cond_var = 1 and u when cond_var = 2.
Setting deriv evaluates a selected first- or second-order derivative of
the density, log-density, or h-function. For example, deriv = "u1"
differentiates with respect to the first argument and
deriv = c("u1", "par2") evaluates the corresponding mixed second
derivative. Derivative order is immaterial. Parameter derivatives use the
natural parameters of the (possibly rotated) copula.
Derivatives are available only for continuous parametric copulas. The backend uses analytic formulas where available and finite-difference fallbacks where needed. Derivatives of inverse h-functions are not available.
When at least one variable is discrete, more than two columns are required
for u: the first block contains realizations of
. The second block contains
realizations of . The minus indicates a
left-sided limit of the cdf. For, e.g., an integer-valued variable, it holds
. For continuous variables the left
limit and the cdf itself coincide. Respective columns can be omitted in the
second block.
dbicop() gives the density or log-density, and pbicop() gives the
distribution function.
rbicop() generates random deviates, and hbicop() gives the h-functions
(and their inverses). If deriv is set, dbicop() and hbicop() return
the selected observation-wise derivative. scores() gives the
observation-wise score matrix and hessian() gives the average Hessian
matrix for a bicop_dist object.
The length of the result is determined by n for rbicop(), and
the number of rows in u for the other functions.
The numerical arguments other than n are recycled to the length of the
result.
For vectorized simulation, parameters must have one row per generated
observation and one column per family parameter.
The functions can optionally be used with a bicop_dist object in place
of the family argument, e.g.,
dbicop(c(0.1, 0.5), bicop_dist("indep")) or
hbicop(c(0.1, 0.5), 2, bicop_dist("indep")).
## evaluate the copula density dbicop(c(0.1, 0.2), "clay", 90, 3) dbicop(c(0.1, 0.2), bicop_dist("clay", 90, 3)) ## evaluate the copula cdf pbicop(c(0.1, 0.2), "clay", 90, 3) ## simulate data plot(rbicop(500, "clay", 90, 3)) ## h-functions joe_cop <- bicop_dist("joe", 0, 3) # h_1(0.1, 0.2) hbicop(c(0.1, 0.2), 1, "bb8", 0, c(2, 0.5)) # h_2^{-1}(0.1, 0.2) hbicop(c(0.1, 0.2), 2, joe_cop, inverse = TRUE) ## derivatives dbicop(c(0.2, 0.4), joe_cop, deriv = "u1") dbicop(c(0.2, 0.4), joe_cop, log = TRUE, deriv = "par1") hbicop(c(0.2, 0.4), 1, joe_cop, deriv = c("u1", "par1")) ## mixed discrete and continuous data x <- cbind(rpois(10, 1), rnorm(10, 1)) u <- cbind(ppois(x[, 1], 1), pnorm(x[, 2]), ppois(x[, 1] - 1, 1)) pbicop(u, "clay", 90, 3, var_types = c("d", "c"))## evaluate the copula density dbicop(c(0.1, 0.2), "clay", 90, 3) dbicop(c(0.1, 0.2), bicop_dist("clay", 90, 3)) ## evaluate the copula cdf pbicop(c(0.1, 0.2), "clay", 90, 3) ## simulate data plot(rbicop(500, "clay", 90, 3)) ## h-functions joe_cop <- bicop_dist("joe", 0, 3) # h_1(0.1, 0.2) hbicop(c(0.1, 0.2), 1, "bb8", 0, c(2, 0.5)) # h_2^{-1}(0.1, 0.2) hbicop(c(0.1, 0.2), 2, joe_cop, inverse = TRUE) ## derivatives dbicop(c(0.2, 0.4), joe_cop, deriv = "u1") dbicop(c(0.2, 0.4), joe_cop, log = TRUE, deriv = "par1") hbicop(c(0.2, 0.4), 1, joe_cop, deriv = c("u1", "par1")) ## mixed discrete and continuous data x <- cbind(rpois(10, 1), rnorm(10, 1)) u <- cbind(ppois(x[, 1], 1), pnorm(x[, 2]), ppois(x[, 1] - 1, 1)) pbicop(u, "clay", 90, 3, var_types = c("d", "c"))
Predictions of the density, distribution function, h-functions (with their inverses) for a bivariate copula model.
## S3 method for class 'bicop_dist' predict(object, newdata, what = "pdf", ...) ## S3 method for class 'bicop' fitted(object, what = "pdf", ...)## S3 method for class 'bicop_dist' predict(object, newdata, what = "pdf", ...) ## S3 method for class 'bicop' fitted(object, what = "pdf", ...)
object |
a |
newdata |
points where the fit shall be evaluated. |
what |
what to predict, one of |
... |
unused. |
fitted() can only be called if the model was fit with the
keep_data = TRUE option.
When at least one variable is discrete, more than two columns are required
for newdata: the first block contains realizations of
. The second block contains
realizations of . The minus indicates a
left-sided limit of the cdf. For, e.g., an integer-valued variable, it holds
. For continuous variables the left
limit and the cdf itself coincide. Respective columns can be omitted in the
second block.
fitted() and logLik() have return values similar to dbicop(),
pbicop(), and hbicop().
# Simulate and fit a bivariate copula model u <- rbicop(500, "gauss", 0, 0.5) fit <- bicop(u, family_set = "par", keep_data = TRUE) # Predictions all.equal(predict(fit, u, "hfunc1"), fitted(fit, "hfunc1"), check.environment = FALSE)# Simulate and fit a bivariate copula model u <- rbicop(500, "gauss", 0, 0.5) fit <- bicop(u, family_set = "par", keep_data = TRUE) # Predictions all.equal(predict(fit, u, "hfunc1"), fitted(fit, "hfunc1"), check.environment = FALSE)
The empirical CDF with tail correction, ensuring that its output is never 0 or 1.
emp_cdf(x)emp_cdf(x)
x |
numeric vector of observations |
The corrected empirical CDF is defined as
A function with signature function(x) that returns .
# fit ECDF on simulated data x <- rnorm(100) cdf <- emp_cdf(x) # output is bounded away from 0 and 1 cdf(-50) cdf(50)# fit ECDF on simulated data x <- rnorm(100) cdf <- emp_cdf(x) # output is bounded away from 0 and 1 cdf(-50) cdf(50)
bicop_dist and vinecop_dist objectsExtracts either the structure matrix (for vinecop_dist only), or
pair-copulas, their parameters, Kendall's taus, or families (for bicop_dist
and vinecop_dist).
get_structure(object) get_pair_copula(object, tree = NA, edge = NA) get_parameters(object, tree = NA, edge = NA) get_ktau(object, tree = NA, edge = NA) get_family(object, tree = NA, edge = NA) get_all_pair_copulas(object, trees = NA) get_all_parameters(object, trees = NA) get_all_ktaus(object, trees = NA) get_all_families(object, trees = NA)get_structure(object) get_pair_copula(object, tree = NA, edge = NA) get_parameters(object, tree = NA, edge = NA) get_ktau(object, tree = NA, edge = NA) get_family(object, tree = NA, edge = NA) get_all_pair_copulas(object, trees = NA) get_all_parameters(object, trees = NA) get_all_ktaus(object, trees = NA) get_all_families(object, trees = NA)
object |
a |
tree |
tree index (not required if |
edge |
edge index (not required if |
trees |
the trees to extract from |
#' The get_structure method (for vinecop_dist or vine_dist
objects only) extracts the structure (see rvine_structure for more
details).
The get_matrix method (for vinecop_dist or vine_dist objects only)
extracts the structure matrix (see rvine_structure for more details).
The other get_xyz methods for vinecop_dist or vine_dist objects return
the entries corresponding to the pair-copula indexed by its tree and
edge. When object is of class bicop_dist, tree and edge are not
required.
get_pair_copula() = the pair-copula itself (see bicop).
get_parameters() = the parameters of the pair-copula (i.e., a numeric
scalar, vector, or matrix).
get_family() = a character for the family (see bicop for
implemented families).
get_ktau() = a numeric scalar with the pair-copula Kendall's tau.
The get_all_xyz methods (for vinecop_dist or vine_dist objects only)
return lists of lists, with each element corresponding to a tree in trees,
and then elements of the sublists correspond to edges. The returned lists
have two additional attributes:
"d" = the dimension of the model.
"trees" = the extracted trees.
The structure matrix, or pair-copulas, their parameters, Kendall's taus, or families.
# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine copula model vc <- vinecop_dist(pcs, mat) # get the structure get_structure(vc) all(get_matrix(vc) == mat) # get pair-copulas get_pair_copula(vc, 1, 1) get_all_pair_copulas(vc) all.equal(get_all_pair_copulas(vc), pcs, check.attributes = FALSE, check.environment = FALSE)# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine copula model vc <- vinecop_dist(pcs, mat) # get the structure get_structure(vc) all(get_matrix(vc) == mat) # get pair-copulas get_pair_copula(vc, 1, 1) get_all_pair_copulas(vc) all.equal(get_all_pair_copulas(vc), pcs, check.attributes = FALSE, check.environment = FALSE)
The arguments configure kde1d::kde1d() and travel with the family object,
rather than being interpreted by vine().
kde1d_family(xmin = NaN, xmax = NaN, mult = 1, bw = NA, deg = 2)kde1d_family(xmin = NaN, xmax = NaN, mult = 1, bw = NA, deg = 2)
xmin, xmax, mult, bw, deg
|
arguments passed to |
A margin-family specification supporting all rvinecopulib variable types.
margin_dist() creates a fitted or fixed margin from evaluation functions
and explicit metadata. The functions must each accept one vector argument
and return a numeric vector of the same length.
margin_dist( d, p, q, family = "custom", type = "c", support = c(-Inf, Inf), npars = 0, loglik = NA_real_ )margin_dist( d, p, q, family = "custom", type = "c", support = c(-Inf, Inf), npars = 0, loglik = NA_real_ )
d |
density or probability mass function. |
p |
distribution function. |
q |
quantile function. |
family |
descriptive family name. |
type |
margin type: |
support |
mathematical support as its lower and upper bounds. |
npars |
effective number of parameters. |
loglik |
maximized marginal log-likelihood, or |
A fitted margin implementing the margin protocol.
normal_margin <- margin_dist( d = function(x) dnorm(x, mean = 1, sd = 2), p = function(x) pnorm(x, mean = 1, sd = 2), q = function(p) qnorm(p, mean = 1, sd = 2), family = "norm", support = c(-Inf, Inf), npars = 2 ) dmargin(1, normal_margin) pmargin(1, normal_margin) qmargin(0.5, normal_margin)normal_margin <- margin_dist( d = function(x) dnorm(x, mean = 1, sd = 2), p = function(x) pnorm(x, mean = 1, sd = 2), q = function(p) qnorm(p, mean = 1, sd = 2), family = "norm", support = c(-Inf, Inf), npars = 2 ) dmargin(1, normal_margin) pmargin(1, normal_margin) qmargin(0.5, normal_margin)
margin_family() defines a candidate family for vine(). Its fitting
function must accept x, weights, and type, and return an object
implementing the fitted-margin protocol. Additional
named fitting arguments can be stored in fit_args.
margin_family(fit, family_name = "custom", types = "c", fit_args = list())margin_family(fit, family_name = "custom", types = "c", fit_args = list())
fit |
function accepting |
family_name |
descriptive family name. |
types |
supported variable types. |
fit_args |
named list of additional arguments passed to |
A margin-family specification.
normal_family <- margin_family( fit = function(x, weights, type) { if (length(weights)) { location <- weighted.mean(x, weights) scale <- sqrt(weighted.mean((x - location)^2, weights)) } else { location <- mean(x) scale <- sd(x) } margin_dist( d = function(y) dnorm(y, location, scale), p = function(y) pnorm(y, location, scale), q = function(p) qnorm(p, location, scale), family = "normal", type = type, npars = 2, loglik = sum(dnorm(x, location, scale, log = TRUE)) ) }, family_name = "normal" )normal_family <- margin_family( fit = function(x, weights, type) { if (length(weights)) { location <- weighted.mean(x, weights) scale <- sqrt(weighted.mean((x - location)^2, weights)) } else { location <- mean(x) scale <- sd(x) } margin_dist( d = function(y) dnorm(y, location, scale), p = function(y) pnorm(y, location, scale), q = function(p) qnorm(p, location, scale), family = "normal", type = type, npars = 2, loglik = sum(dnorm(x, location, scale, log = TRUE)) ) }, family_name = "normal" )
Margin families describe how a candidate is fitted. Implementations provide
fit_margin() and margin_info() methods. Their margin_info() result must
contain family_name and types. Every fitter receives the observations,
observation weights, and requested variable type. It must either use the
weights or reject them.
fit_margin(family, x, weights = numeric(), type = "c")fit_margin(family, x, weights = numeric(), type = "c")
family |
a margin-family specification. |
x |
observed values. |
weights |
observation weights, or |
type |
requested variable type. |
fit_margin() returns a fitted margin.
Fitted margins provide distribution evaluation and model metadata through a
small set of S3 generics. Methods for dmargin(), pmargin(), and
qmargin() dispatch on margin, their second argument. margin_info()
dispatches on its only argument.
dmargin(x, margin) pmargin(x, margin) qmargin(p, margin) margin_info(object)dmargin(x, margin) pmargin(x, margin) qmargin(p, margin) margin_info(object)
x |
vector of evaluation points. |
margin |
a fitted marginal distribution object. |
p |
vector of probabilities. |
object |
a fitted margin or margin-family specification. |
A fitted margin class must implement all four generics documented here.
Evaluation methods must be vectorized. margin_info() returns a list with
entries family_name, type, support, npars, and loglik. type is
"c" for a continuous distribution, "d" for an integer-supported
distribution, or "zi" for a continuous distribution with an atom at zero.
support contains the mathematical lower and upper bounds. npars must be
non-negative or NA when unavailable; loglik may be NA for a fixed
margin.
The type determines left limits throughout rvinecopulib. They are F(x) for
continuous margins, F(x - 1) for integer-supported margins, and
F(0) - P(X = 0) at zero for zero-inflated margins. Consequently,
dmargin(0, margin) must return the atom's mass for a zero-inflated margin.
dmargin(), pmargin(), and qmargin() return numeric vectors.
margin_info() returns a list of model information.
Calculates the modified vine copula Bayesian information criterion.
mBICV(object, psi0 = 0.9, newdata = NULL)mBICV(object, psi0 = 0.9, newdata = NULL)
object |
a fitted |
psi0 |
baseline prior probability of a non-independence copula. |
newdata |
optional; a new data set. |
The modified vine copula Bayesian information criterion (mBICv) is defined as
where is the log-likelihood and is the
(effective) number of parameters of the model, is the tree level
is the prior probability of having a non-independence copula and
is the number of non-independence copulas in tree . The
mBICv is a consistent model selection criterion for parametric sparse vine
copula models.
Nagler, T., Bumann, C., Czado, C. (2019). Model selection for sparse high-dimensional vine copulas with application to portfolio risk. Journal of Multivariate Analysis, in press (https://arxiv.org/pdf/1801.09739)
u <- sapply(1:5, function(i) runif(50)) fit <- vinecop(u, family_set = "par", keep_data = TRUE) mBICV(fit, 0.9) # with a 0.9 prior probability of a non-independence copula mBICV(fit, 0.1) # with a 0.1 prior probability of a non-independence copulau <- sapply(1:5, function(i) runif(50)) fit <- vinecop(u, family_set = "par", keep_data = TRUE) mBICV(fit, 0.9) # with a 0.9 prior probability of a non-independence copula mBICV(fit, 0.1) # with a 0.1 prior probability of a non-independence copula
This function provides pair plots for copula data. It shows bivariate contour plots on the lower panel, scatter plots and correlations on the upper panel and histograms on the diagonal panel.
pairs_copula_data(data, main = "", ...)pairs_copula_data(data, main = "", ...)
data |
the data (must lie in the unit hypercube). |
main |
an overall title for the plot. |
... |
other parameters passed to |
u <- replicate(3, runif(100)) pairs_copula_data(u)u <- replicate(3, runif(100)) pairs_copula_data(u)
Conversion between Kendall's tau and parameters
par_to_ktau(family, rotation, parameters) ktau_to_par(family, tau)par_to_ktau(family, rotation, parameters) ktau_to_par(family, tau)
family |
a copula family (see |
rotation |
the rotation of the copula, one of |
parameters |
vector or matrix of copula parameters, not used when
|
tau |
Kendall's |
# the following are equivalent par_to_ktau(bicop_dist("clayton", 0, 3)) par_to_ktau("clayton", 0, 3) ktau_to_par("clayton", 0.5) ktau_to_par(bicop_dist("clayton", 0, 3), 0.5)# the following are equivalent par_to_ktau(bicop_dist("clayton", 0, 3)) par_to_ktau("clayton", 0, 3) ktau_to_par("clayton", 0.5) ktau_to_par(bicop_dist("clayton", 0, 3), 0.5)
Covariance matrices and Wald confidence intervals for the estimated pair-copula parameters.
## S3 method for class 'vinecop_dist' vcov(object, newdata = NULL, step_wise = TRUE, cores = 1, ...) ## S3 method for class 'bicop_dist' vcov(object, newdata = NULL, step_wise = TRUE, cores = 1, ...) ## S3 method for class 'vinecop_dist' confint(object, parm, level = 0.95, ...) ## S3 method for class 'bicop_dist' confint(object, parm, level = 0.95, ...)## S3 method for class 'vinecop_dist' vcov(object, newdata = NULL, step_wise = TRUE, cores = 1, ...) ## S3 method for class 'bicop_dist' vcov(object, newdata = NULL, step_wise = TRUE, cores = 1, ...) ## S3 method for class 'vinecop_dist' confint(object, parm, level = 0.95, ...) ## S3 method for class 'bicop_dist' confint(object, parm, level = 0.95, ...)
object |
an object of class |
newdata |
copula data to evaluate the derivatives at. Defaults to the
data stored in |
step_wise |
if |
cores |
number of cores to use; if larger than one, computations are
done in parallel on |
... |
passed on to |
parm |
a specification of which parameters to give intervals for: a vector of names or of positions. Defaults to all of them. |
level |
the confidence level. |
The estimator solves with
, so it is an M-estimator.
Writing
vcov() returns the sandwich of White (1982).
step_wise selects the estimating equation, and it changes the nature of
.
With step_wise = FALSE, is the gradient of the joint
log-likelihood and is minus its averaged Hessian: a symmetric
matrix, and the classical setting for both formulas.
With step_wise = TRUE (the default), the pseudo-observations entering each
pair copula are held fixed. This is the estimating equation that sequential,
tree-by-tree estimation actually solves, but it is not the gradient of any
objective: is a Jacobian rather than a Hessian, and it is block
triangular rather than symmetric, which is why the sandwich is written with
on the right.
No derivatives are implemented for nonparametric pair copulas, so a model containing one is refused. For models with discrete variables the backend differentiates by central finite differences rather than in closed form; the result is correct but carries the accuracy of a finite-difference approximation.
These functions treat the copula data as observed, so they quantify
uncertainty in the pair-copula parameters given the margins. A model
fitted by vine() estimates its margins too, and the second stage inherits
the first stage's sampling variability: intervals that ignore it are too
short. In a small experiment (three-dimensional Gaussian D-vine,
equicorrelation 0.5, n = 1000, 1000 replications) nominal 95% intervals
attained 0.956 when the true margins were used and 0.929 when they were
estimated by ranks.
Correcting for it needs the influence function of the marginal estimator,
and fit_margin() is a user-supplied black box, so there is nothing to
differentiate. We therefore provide no method for "vine" objects rather
than a silently anticonservative one. To quantify uncertainty for a complete
vine distribution, resample:
boot <- replicate(B, {
idx <- sample(nrow(x), replace = TRUE)
fit <- vine(x[idx, ], margins_controls = mc, copula_controls = cc)
unlist(get_all_parameters(fit$copula))
})
apply(boot, 1, quantile, c(0.025, 0.975))
vcov() returns a covariance matrix with one row and column per
estimated pair-copula parameter, ordered by tree, then edge, then parameter
within the pair copula, with dimnames giving that position. confint()
returns a matrix with columns for the lower and upper endpoints.
White H (1982). "Maximum Likelihood Estimation of Misspecified Models." Econometrica, 50(1), 1–25. doi:10.2307/1912526
Stoeber J, Schepsmeier U (2013). "Estimating Standard Errors in Regular Vine Copula Models." Computational Statistics, 28(6), 2679–2707. doi:10.1007/s00180-013-0423-8
scores(), hessian(), vinecop(), vine()
# bivariate copula ------------------------------------------ u <- rbicop(500, "gaussian", 0, 0.5) fit <- bicop(u, family_set = "parametric", keep_data = TRUE) sqrt(diag(vcov(fit))) confint(fit) # vine copula ----------------------------------------------- u <- pseudo_obs(matrix(rnorm(500 * 3), 500, 3)) fit <- vinecop(u, family_set = "parametric", keep_data = TRUE) vcov(fit) confint(fit, level = 0.9) # derivatives of the joint log-likelihood instead of the step-wise ones sqrt(diag(vcov(fit, step_wise = FALSE)))# bivariate copula ------------------------------------------ u <- rbicop(500, "gaussian", 0, 0.5) fit <- bicop(u, family_set = "parametric", keep_data = TRUE) sqrt(diag(vcov(fit))) confint(fit) # vine copula ----------------------------------------------- u <- pseudo_obs(matrix(rnorm(500 * 3), 500, 3)) fit <- vinecop(u, family_set = "parametric", keep_data = TRUE) vcov(fit) confint(fit, level = 0.9) # derivatives of the joint log-likelihood instead of the step-wise ones sqrt(diag(vcov(fit, step_wise = FALSE)))
bicop_dist and bicop objectsThere are several options for plotting bicop_dist objects. The density of a bivariate copula density can be visualized as surface/perspective or contour plot. Optionally, the density can be coupled with standard normal margins (default for contour plots).
## S3 method for class 'bicop_dist' plot(x, type = "surface", margins, size, ...) ## S3 method for class 'bicop' plot(x, type = "surface", margins, size, ...) ## S3 method for class 'bicop_dist' contour(x, margins = "norm", size = 100L, ...) ## S3 method for class 'bicop' contour(x, margins = "norm", size = 100L, ...)## S3 method for class 'bicop_dist' plot(x, type = "surface", margins, size, ...) ## S3 method for class 'bicop' plot(x, type = "surface", margins, size, ...) ## S3 method for class 'bicop_dist' contour(x, margins = "norm", size = 100L, ...) ## S3 method for class 'bicop' contour(x, margins = "norm", size = 100L, ...)
x |
|
type |
plot type; either |
margins |
options are: |
size |
integer; the plot is based on values on a |
... |
optional arguments passed to |
bicop_dist(), graphics::contour(), lattice::wireframe()
## construct bicop_dist object for a student t copula obj <- bicop_dist(family = "t", rotation = 0, parameters = c(0.7, 4)) ## plots plot(obj) # surface plot of copula density contour(obj) # contour plot with standard normal margins contour(obj, margins = "unif") # contour plot of copula density## construct bicop_dist object for a student t copula obj <- bicop_dist(family = "t", rotation = 0, parameters = c(0.7, 4)) ## plots plot(obj) # surface plot of copula density contour(obj) # contour plot with standard normal margins contour(obj, margins = "unif") # contour plot of copula density
Plot one or all trees of an R-vine structure.
## S3 method for class 'rvine_structure' plot(x, ...) ## S3 method for class 'rvine_matrix' plot(x, ...)## S3 method for class 'rvine_structure' plot(x, ...) ## S3 method for class 'rvine_matrix' plot(x, ...)
x |
an |
... |
passed to |
plot(cvine_structure(1:5)) plot(rvine_structure_sim(5)) mat <- rbind(c(1, 1, 1), c(2, 2, 0), c(3, 0, 0)) plot(rvine_matrix(mat)) plot(rvine_matrix_sim(5))plot(cvine_structure(1:5)) plot(rvine_structure_sim(5)) mat <- rbind(c(1, 1, 1), c(2, 2, 0), c(3, 0, 0)) plot(rvine_matrix(mat)) plot(rvine_matrix_sim(5))
vinecop_dist and vinecop objects.There are two plotting generics for vinecop_dist objects.
plot.vinecop_dist plots one or all trees of a given R-vine copula
model. Edges can be labeled with information about the corresponding
pair-copula. contour.vinecop_dist produces a matrix of contour plots
(using plot.bicop).
## S3 method for class 'vinecop_dist' plot(x, tree = 1, var_names = "ignore", edge_labels = NULL, ...) ## S3 method for class 'vinecop' plot(x, tree = 1, var_names = "ignore", edge_labels = NULL, ...) ## S3 method for class 'vinecop_dist' contour(x, tree = "ALL", cex.nums = 1, ...) ## S3 method for class 'vinecop' contour(x, tree = "ALL", cex.nums = 1, ...)## S3 method for class 'vinecop_dist' plot(x, tree = 1, var_names = "ignore", edge_labels = NULL, ...) ## S3 method for class 'vinecop' plot(x, tree = 1, var_names = "ignore", edge_labels = NULL, ...) ## S3 method for class 'vinecop_dist' contour(x, tree = "ALL", cex.nums = 1, ...) ## S3 method for class 'vinecop' contour(x, tree = "ALL", cex.nums = 1, ...)
x |
|
tree |
|
var_names |
integer; specifies how to make use of variable names:
|
edge_labels |
character; options are:
|
... |
Unused for |
cex.nums |
numeric; expansion factor for font of the numbers. |
If you want the contour boxes to be perfect squares, the plot height should
be 1.25/length(tree)*(d - min(tree)) times the plot width.
The plot() method returns an object that (among other things) contains the
igraph representation of the graph; see Examples.
Thomas Nagler, Thibault Vatter
# set up vine copula model u <- matrix(runif(20 * 10), 20, 10) vc <- vinecop(u, family_set = "indep") # plot plot(vc, tree = c(1, 2)) plot(vc, edge_labels = "pair") # extract igraph representation plt <- plot(vc, edge_labels = "family_tau") igr_obj <- get("g", plt$plot_env)[[1]] igr_obj # print object igraph::E(igr_obj)$name # extract edge labels # set up another vine copula model pcs <- lapply(1:3, function(j) # pair-copulas in tree j lapply(runif(4 - j), function(cor) bicop_dist("gaussian", 0, cor))) mat <- rvine_matrix_sim(4) vc <- vinecop_dist(pcs, mat) # contour plot contour(vc)# set up vine copula model u <- matrix(runif(20 * 10), 20, 10) vc <- vinecop(u, family_set = "indep") # plot plot(vc, tree = c(1, 2)) plot(vc, edge_labels = "pair") # extract igraph representation plt <- plot(vc, edge_labels = "family_tau") igr_obj <- get("g", plt$plot_env)[[1]] igr_obj # print object igraph::E(igr_obj)$name # extract edge labels # set up another vine copula model pcs <- lapply(1:3, function(j) # pair-copulas in tree j lapply(runif(4 - j), function(cor) bicop_dist("gaussian", 0, cor))) mat <- rvine_matrix_sim(4) vc <- vinecop_dist(pcs, mat) # contour plot contour(vc)
Compute the pseudo-observations for the given data matrix.
pseudo_obs(x, ties_method = "average", lower_tail = TRUE)pseudo_obs(x, ties_method = "average", lower_tail = TRUE)
x |
vector or matrix random variates to be converted (column wise) to pseudo-observations. |
ties_method |
similar to |
lower_tail |
|
Given n realizations ,
of a random vector X, the pseudo-observations are defined via
for
and
, where
denotes the rank of among all ,
.
The pseudo-observations can thus also be computed by component-wise applying
the empirical distribution functions to the data and scaling the result by
. This asymptotically negligible scaling factor is used to force
the variates to fall inside the open unit hypercube, for example, to avoid
problems with density evaluation at the boundaries.
When lower_tail = FALSE, then pseudo_obs() simply returns
1 - pseudo_obs().
a vector of matrix of the same dimension as the input containing the pseudo-observations.
# pseudo-observations for a vector pseudo_obs(rnorm(10)) # pseudo-observations for a matrix pseudo_obs(cbind(rnorm(10), rnorm(10)))# pseudo-observations for a vector pseudo_obs(rnorm(10)) # pseudo-observations for a matrix pseudo_obs(cbind(rnorm(10), rnorm(10)))
The Rosenblatt transform maps observations from a model to independent uniform variates. The inverse transform evaluates conditional quantiles and maps independent uniforms to draws from the model; see Details.
rosenblatt( x, model, cores = 1, randomize_discrete = TRUE, conditioning_set = NULL ) inverse_rosenblatt(u, model, cores = 1, conditioning_set = NULL)rosenblatt( x, model, cores = 1, randomize_discrete = TRUE, conditioning_set = NULL ) inverse_rosenblatt(u, model, cores = 1, conditioning_set = NULL)
x |
matrix of evaluation points; must be in |
model |
a model object; classes currently supported are
|
cores |
if |
randomize_discrete |
Whether to randomize the transform for discrete variables; see Details. |
conditioning_set |
optional variable indices or names that define the
conditioning variables. The transform uses an admissible sampling order
whose tail contains exactly this set. The model is not modified. If |
u |
matrix of evaluation points; must be in |
Let be the sampling order of a random vector
with distribution . For continuous
variables, the Rosenblatt transform is defined by
If the model is correct, the components of are independent standard
uniforms. The result is stored in the original variable columns; the
sampling order determines only the sequence of conditional distributions.
The inverse transform applies the corresponding conditional quantiles:
Thus has distribution when contains
independent standard uniforms.
If a variable has atoms, its conditional cdf jumps at the observation. Let
denote the conditional cdf of given the preceding
variables in the sampling order. Following Brockwell
(10.1016/j.spl.2007.02.008), rosenblatt() returns
where is the left limit and the are
independent standard uniforms. This randomization is used by default and
yields uniform components under the fitted model. Set
randomize_discrete = FALSE to return the upper endpoint
instead.
# simulate data with some dependence x <- replicate(3, rnorm(200)) x[, 2:3] <- x[, 2:3] + x[, 1] pairs(x) # estimate a vine distribution model fit <- vine(x, copula_controls = list(family_set = "par")) # transform into independent uniforms u <- rosenblatt(x, fit) pairs(u) # inversion pairs(inverse_rosenblatt(u, fit)) # works similarly for vinecop models vc <- fit$copula rosenblatt(pseudo_obs(x), vc)# simulate data with some dependence x <- replicate(3, rnorm(200)) x[, 2:3] <- x[, 2:3] + x[, 1] pairs(x) # estimate a vine distribution model fit <- vine(x, copula_controls = list(family_set = "par")) # transform into independent uniforms u <- rosenblatt(x, fit) pairs(u) # inversion pairs(inverse_rosenblatt(u, fit)) # works similarly for vinecop models vc <- fit$copula rosenblatt(pseudo_obs(x), vc)
R-vine structures are compressed representations encoding the tree structure
of the vine, i.e. the conditioned/conditioning variables of each edge. The
functions [cvine_structure()] or [dvine_structure()] give a simpler way
to construct C-vines (every tree is a star) and D-vines (every tree is a
path), respectively (see Examples).
rvine_structure(order, struct_array = list(), is_natural_order = FALSE) cvine_structure(order, trunc_lvl = Inf) dvine_structure(order, trunc_lvl = Inf) rvine_matrix(matrix)rvine_structure(order, struct_array = list(), is_natural_order = FALSE) cvine_structure(order, trunc_lvl = Inf) dvine_structure(order, trunc_lvl = Inf) rvine_matrix(matrix)
order |
a vector of positive integers. |
struct_array |
a list of vectors of positive integers. The vectors
represent rows of the r-rvine structure and the number of elements have to
be compatible with the |
is_natural_order |
whether |
trunc_lvl |
the truncation level |
matrix |
an R-vine matrix, see Details. |
The R-vine structure is essentially a lower-triangular matrix/triangular array, with a notation that differs from the one in the VineCopula package. An example array is
4 4 4 4 3 3 3 2 2 1
which encodes the following pair-copulas:
| tree | edge | pair-copulas |
| 0 | 0 | (1, 4) |
| 1 | (2, 4) |
|
| 2 | (3, 4) |
|
| 1 | 0 | (1, 3; 4) |
| 1 | (2, 3; 4) |
|
| 2 | 0 | (1, 2; 3, 4)
|
An R-vine structure can be converted to an R-vine matrix using
as_rvine_matrix(), which encodes the same model with a square matrix filled
with zeros. For instance, the matrix corresponding to the structure above is:
4 4 4 4 3 3 3 0 2 2 0 0 1 0 0 0
Similarly, an R-vine matrix can be converted to an R-vine structure using
as_rvine_structure().
Denoting by M[i, j] the array entry in row i and column j (the
pair-copula index for edge e in tree t of a d dimensional vine is
(M[d + 1 - e, e], M[t, e]; M[t - 1, e], ..., M[1, e]). Less formally,
Start with the counter-diagonal element of column e (first conditioned
variable).
Jump up to the element in row t (second conditioned variable).
Gather all entries further up in column e (conditioning set).
Internally, the diagonal is stored separately from the off-diagonal elements, which are stored as a triangular array. For instance, the off-diagonal elements off the structure above are stored as
4 4 4 3 3 2
for the structure above. The reason is that it allows for parsimonious representations of truncated models. For instance, the 2-truncated model is represented by the same diagonal and the following truncated triangular array:
4 4 4 3 3
A valid R-vine structure or matrix must satisfy several conditions which are
checked when rvine_structure(), rvine_matrix(), or some coercion methods
(see as_rvine_structure() and as_rvine_matrix() are called:
It can only contain numbers between 1 and d (and additionally zeros for R-vine matrices).
The anti-diagonal must contain the numbers 1, ..., d.
The anti-diagonal entry of a column must not be contained in any column further to the right.
The entries of a column must be contained in all columns to the left.
The proximity condition must hold: For all t = 1, ..., d - 2 and e = 1,
..., d - t there must exist an index j > d, such that
(M[t, e], {M[1, e], ..., M[t - 1, e]}) equals either
(M[d + 1 - j, j], {M[1, j], ..., M[t - 1, j]}) or
(M[t - 1, j], {M[d + 1 - j, j], M[1, j], ..., M[t - 2, j]}).
Condition 5 already implies conditions 2-4, but is more difficult to check by hand.
Either an rvine_structure or an rvine_matrix.
as_rvine_structure(), as_rvine_matrix(),
plot.rvine_structure(), plot.rvine_matrix(),
rvine_structure_sim(), rvine_matrix_sim()
# R-vine structures can be constructed from the order vector and struct_array rvine_structure(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3), 2 )) # R-vine matrices can be constructed from standard matrices mat <- matrix(c(4, 3, 2, 1, 4, 3, 2, 0, 4, 3, 0, 0, 4, 0, 0, 0), 4, 4) rvine_matrix(mat) # coerce to R-vine structure str(as_rvine_structure(mat)) # truncate and construct the R-vine matrix mat[3, 1] <- 0 rvine_matrix(mat) # or use directly the R-vine structure constructor rvine_structure(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3) )) # throws an error mat[3, 1] <- 5 try(rvine_matrix(mat)) # C-vine structure cvine <- cvine_structure(1:5) cvine plot(cvine) # D-vine structure dvine <- dvine_structure(c(1, 4, 2, 3, 5)) dvine plot(dvine)# R-vine structures can be constructed from the order vector and struct_array rvine_structure(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3), 2 )) # R-vine matrices can be constructed from standard matrices mat <- matrix(c(4, 3, 2, 1, 4, 3, 2, 0, 4, 3, 0, 0, 4, 0, 0, 0), 4, 4) rvine_matrix(mat) # coerce to R-vine structure str(as_rvine_structure(mat)) # truncate and construct the R-vine matrix mat[3, 1] <- 0 rvine_matrix(mat) # or use directly the R-vine structure constructor rvine_structure(order = 1:4, struct_array = list( c(4, 4, 4), c(3, 3) )) # throws an error mat[3, 1] <- 5 try(rvine_matrix(mat)) # C-vine structure cvine <- cvine_structure(1:5) cvine plot(cvine) # D-vine structure dvine <- dvine_structure(c(1, 4, 2, 3, 5)) dvine plot(dvine)
Simulates from a uniform distribution over all R-vine structures on d
variables. rvine_structure_sim() returns an rvine_structure() object,
rvine_matrix_sim() an rvine_matrix().
rvine_structure_sim(d, natural_order = FALSE) rvine_matrix_sim(d, natural_order = FALSE)rvine_structure_sim(d, natural_order = FALSE) rvine_matrix_sim(d, natural_order = FALSE)
d |
the number of variables |
natural_order |
boolean; whether the structures should be in natural
order (counter-diagonal is |
rvine_structure(), rvine_matrix(),
plot.rvine_structure(), plot.rvine_matrix()
rvine_structure_sim(10) rvine_structure_sim(10, natural_order = TRUE) # counter-diagonal is 1:d rvine_matrix_sim(10)rvine_structure_sim(10) rvine_structure_sim(10, natural_order = TRUE) # counter-diagonal is 1:d rvine_matrix_sim(10)
rvinecopulib provides high-performance tools for constructing, fitting, selecting, simulating, and visualizing bivariate and vine copula models. It also fits complete multivariate distributions by combining marginal models with a vine copula. Continuous, integer-valued discrete, mixed, and zero-inflated variables are supported.
The package is the R interface to the header-only 'vinecopulib' C++ library, which is bundled with rvinecopulib. Users do not need to install the C++ library separately. 'vinecopulib' is licensed under the MIT License and rvinecopulib under the GNU GPL version 3.
rvinecopulib provides:
parametric, nonparametric, rotated, and extreme-value pair copulas;
automatic pair-copula family, vine structure, truncation, and threshold selection;
complete multivariate distributions with nonparametric, parametric, or custom margins;
continuous, integer-valued discrete, mixed, and zero-inflated data;
observation weights and parallel fitting;
conditional simulation and Rosenblatt transforms;
likelihood scores, Hessians, and observation-specific parameters;
custom marginal families and custom tree-selection criteria.
The framework is exposed at three levels. vine() models observations on
their original scale by fitting the marginal distributions and dependence
model together. vinecop() models uniform pseudo-observations when only the
dependence model is required. bicop() provides the corresponding
two-variable workflow.
The constructors vine_dist(), vinecop_dist(), and bicop_dist() create
models from explicitly specified components.
Start with vignette("getting-started"). Dedicated articles cover
bivariate copulas, vine structures and selection, marginal models,
discrete data, conditional simulation, and likelihood derivatives. The
package website contains the
rendered articles and complete function reference.
Thomas Nagler, Thibault Vatter
Useful links:
Report bugs at https://github.com/vinecopulib/rvinecopulib/issues
## complete distribution on the original data scale x <- data.frame(a = rnorm(50), b = rexp(50), n = rpois(50, 2)) fit <- vine( x, var_types = c("c", "c", "d"), copula_controls = list(family_set = "onepar") ) rvine(3, fit) ## dependence model on the copula scale u <- pseudo_obs(as.matrix(USArrests)) copula_fit <- vinecop(u, family_set = "onepar") summary(copula_fit) ## bivariate model uv <- rbicop(100, "gaussian", 0, 0.5) bicop(uv, family_set = "par")## complete distribution on the original data scale x <- data.frame(a = rnorm(50), b = rexp(50), n = rpois(50, 2)) fit <- vine( x, var_types = c("c", "c", "d"), copula_controls = list(family_set = "onepar") ) rvine(3, fit) ## dependence model on the copula scale u <- pseudo_obs(as.matrix(USArrests)) copula_fit <- vinecop(u, family_set = "onepar") summary(copula_fit) ## bivariate model uv <- rbicop(100, "gaussian", 0, 0.5) bicop(uv, family_set = "par")
stats_margin() adapts the univariate distributions documented in
stats::Distributions to the fitted-margin protocol. Distribution
parameters are supplied through ...; the resulting margin reports the
distribution's parameter count but has no fitted log-likelihood.
stats_margin(family, ...)stats_margin(family, ...)
family |
distribution suffix such as |
... |
parameters passed to the distribution's density, CDF, and quantile functions. |
A fixed margin implementing the margin protocol.
normal <- stats_margin("norm", mean = 1, sd = 2) poisson <- stats_margin("pois", lambda = 3) margin_info(normal) margin_info(poisson)normal <- stats_margin("norm", mean = 1, sd = 2) poisson <- stats_margin("pois", lambda = 3) margin_info(normal) margin_info(poisson)
Extracts a truncated sub-vine based on a truncation level supplied by user.
truncate_model(object, ...) ## S3 method for class 'rvine_structure' truncate_model(object, trunc_lvl, ...) ## S3 method for class 'rvine_matrix' truncate_model(object, trunc_lvl, ...) ## S3 method for class 'vinecop_dist' truncate_model(object, trunc_lvl, ...) ## S3 method for class 'vine_dist' truncate_model(object, trunc_lvl, ...)truncate_model(object, ...) ## S3 method for class 'rvine_structure' truncate_model(object, trunc_lvl, ...) ## S3 method for class 'rvine_matrix' truncate_model(object, trunc_lvl, ...) ## S3 method for class 'vinecop_dist' truncate_model(object, trunc_lvl, ...) ## S3 method for class 'vine_dist' truncate_model(object, trunc_lvl, ...)
object |
a model object. |
... |
further arguments passed to specific methods. |
trunc_lvl |
tree level after which the vine copula should be truncated. |
While a vine model for a d dimensional random vector contains at most d-1
nested trees, this function extracts a sub-model based on a given truncation
level.
For instance, truncate_model(object, 1) results in a 1-truncated
vine (i.e., a vine with a single tree). Similarly truncate_model(object, 2)
results in a 2-truncated vine (i.e., a vine with two trees). Note that
truncate_model(truncate_model(object, 1), 2) returns a 1-truncated vine.
# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine structure structure <- as_rvine_structure(mat) # truncate the model truncate_model(structure, 1) # set up vine copula model vc <- vinecop_dist(pcs, mat) # truncate the model truncate_model(vc, 1)# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine structure structure <- as_rvine_structure(mat) # truncate the model truncate_model(structure, 1) # set up vine copula model vc <- vinecop_dist(pcs, mat) # truncate the model truncate_model(vc, 1)
univariateML_family() adapts one of univariateML's available models to the
margin-family protocol.
univariateML_family(family)univariateML_family(family)
family |
a univariateML family name. |
A parametric margin-family specification.
Automated fitting or creation of custom vine copula models
vine( data, margins_controls = list(), copula_controls = list(family_set = "all", structure = NA, par_method = "mle", nonpar_method = "constant", mult = 1, selcrit = "aic", psi0 = 0.9, presel = TRUE, allow_rotations = TRUE, trunc_lvl = Inf, tree_crit = "tau", threshold = 0, keep_data = FALSE, show_trace = FALSE, cores = 1, tree_algorithm = "mst_prim", conditioning_set = NULL), weights = numeric(), keep_data = FALSE, cores = 1, var_types = NULL ) vine_dist(margins, pair_copulas, structure)vine( data, margins_controls = list(), copula_controls = list(family_set = "all", structure = NA, par_method = "mle", nonpar_method = "constant", mult = 1, selcrit = "aic", psi0 = 0.9, presel = TRUE, allow_rotations = TRUE, trunc_lvl = Inf, tree_crit = "tau", threshold = 0, keep_data = FALSE, show_trace = FALSE, cores = 1, tree_algorithm = "mst_prim", conditioning_set = NULL), weights = numeric(), keep_data = FALSE, cores = 1, var_types = NULL ) vine_dist(margins, pair_copulas, structure)
data |
a matrix or data.frame. Discrete variables have to be declared as
|
margins_controls |
a list with arguments to be passed to marginal fitting. The supported entries are
|
copula_controls |
a list with arguments to be passed to |
weights |
optional numeric vector of weights for each observation, used for both marginal and copula estimation. Every margin family receives these weights and is responsible for using or explicitly rejecting them. |
keep_data |
whether the original data should be stored; if you want to
store the pseudo-observations used for fitting the copula, use the
|
cores |
the number of cores to use for parallel computations. Unless
overridden by |
var_types |
optional variable types, one for each variable after factor
expansion: |
margins |
a list containing one marginal distribution per variable.
Each margin can be a |
pair_copulas |
A nested list of 'bicop_dist' objects, where
|
structure |
an |
vine_dist() creates a vine copula by specifying the margins, a nested list
of bicop_dist objects and a quadratic structure matrix.
vine() provides automated fitting for vine copula models. Margin-specific
fitting options are properties of the family specification, not controls
interpreted by vine(). For example, use
kde1d_family(xmin = 0, mult = 2) in family_set. copula_controls is a
list with the same parameters as vinecop() (except for data).
A character margins_controls$family_set is used for every variable. Its
entries can be "kde1d" or names of models available in univariateML. A
list supplies one candidate set per
variable, for example list(c("norm", "cauchy"), c("pois", "geom")).
Each entry may itself be a character vector, a margin-family object, or a
list combining both. An unnamed list is matched by position; a named list
must contain every expanded variable name exactly once and is reordered to
match the data.
The aliases "par" and "parametric" expand to all univariateML
families, "nonpar" and "nonparametric" select "kde1d", and "all"
combines both. Parametric names and aliases require the suggested
univariateML package; the default "kde1d", kde1d_family(), and custom
margin_family() objects do not. Identical candidate specifications are
fitted only once.
Candidates incompatible with a variable's type are removed before fitting.
rvinecopulib fits every remaining candidate and selects the first minimum of
negative log-likelihood, AIC, or BIC according to selcrit. Competing
candidates therefore need a finite margin_info()$loglik value, and AIC or
BIC additionally requires a finite margin_info()$npars. A sole candidate
without a parameter count is retained with a warning, but model AIC and BIC
are then unavailable. An unsupported candidate may fail while another
candidate succeeds; such failures and all candidate warnings are reported.
Protocol violations always stop fitting. The built-in univariateML family
explicitly rejects observation weights.
Objects inheriting from vine_dist for vine_dist(), and
vine and vine_dist for vine().
Objects from the vine_dist class are lists containing:
margins, a list of marginals (see below).
copula, an object of the class vinecop_dist, see vinecop_dist().
For objects from the vine class, copula is also an object of the class
vine, see vinecop(). Additionally, objects from the vine class contain:
margins_controls, a list with the controls used to fit and select the
marginal families.
copula_controls, a list with the set of fit controls that was passed
to vinecop() when estimating the copula.
data (optionally, if keep_data = TRUE was used), the dataset that was
passed to vine().
nobs, an integer containing the number of observations that was used
to fit the model.
Concerning margins:
For objects created with vine_dist(), it simply corresponds to the margins
argument.
For objects created with vine(), it is a list of fitted objects
implementing the margin protocol.
# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine copula model with Gaussian margins vc <- vine_dist(list(stats_margin("norm")), pcs, mat) # show model summary(vc) # simulate some data x <- rvine(50, vc) # estimate a vine copula model fit <- vine(x, copula_controls = list(family_set = "par")) summary(fit) # parametric margin selection (requires the suggested univariateML package) if (requireNamespace("univariateML", quietly = TRUE)) { fit_par_margins <- vine( x, margins_controls = list(family_set = c("norm", "cauchy")) ) } ## model for discrete data x <- as.data.frame(x) x[, 1] <- ordered(round(x[, 1]), levels = seq.int(-5, 5)) fit_disc <- vine(x, copula_controls = list(family_set = "par")) summary(fit_disc)# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine copula model with Gaussian margins vc <- vine_dist(list(stats_margin("norm")), pcs, mat) # show model summary(vc) # simulate some data x <- rvine(50, vc) # estimate a vine copula model fit <- vine(x, copula_controls = list(family_set = "par")) summary(fit) # parametric margin selection (requires the suggested univariateML package) if (requireNamespace("univariateML", quietly = TRUE)) { fit_par_margins <- vine( x, margins_controls = list(family_set = c("norm", "cauchy")) ) } ## model for discrete data x <- as.data.frame(x) x[, 1] <- ordered(round(x[, 1]), levels = seq.int(-5, 5)) fit_disc <- vine(x, copula_controls = list(family_set = "par")) summary(fit_disc)
Density, distribution function and random generation for the vine based distribution.
dvine(x, vine, cores = 1, log = FALSE) pvine(x, vine, n_mc = 10^4, cores = 1) rvine(n, vine, qrng = FALSE, cores = 1, x_cond = NULL, conditioning_set = NULL)dvine(x, vine, cores = 1, log = FALSE) pvine(x, vine, n_mc = 10^4, cores = 1) rvine(n, vine, qrng = FALSE, cores = 1, x_cond = NULL, conditioning_set = NULL)
x |
evaluation points, either a length d vector or a d-column matrix, where d is the number of variables in the vine. |
vine |
an object of class |
cores |
number of cores to use; if larger than one, computations are
done in parallel on |
log |
if |
n_mc |
number of samples used for quasi Monte Carlo integration. |
n |
number of observations. |
qrng |
if |
x_cond |
optional conditioning values for |
conditioning_set |
variable indices or names corresponding to the
columns of |
See vine for the estimation and construction of vine models. Here, the density, distribution function and random generation for the vine distributions are standard.
The functions are based on dvinecop(), pvinecop() and rvinecop() for
vinecop objects. Margins are evaluated through dmargin(), pmargin(),
and qmargin(). Methods are provided for margins fitted by vine() and for
the fixed stats::Distributions specifications accepted by vine_dist().
dvine() gives the density, pvine() gives the distribution function,
and rvine() generates unconditional or conditional random deviates.
The length of the result is determined by n for rvine(), and
the number of rows in u for the other functions.
The vine object is recycled to the length of the
result.
# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # set up vine copula model mat <- rvine_matrix_sim(3) vc <- vine_dist(list(stats_margin("norm")), pcs, mat) # simulate from the model x <- rvine(200, vc) pairs(x) # evaluate the density and cdf dvine(x[1, ], vc) pvine(x[1, ], vc)# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # set up vine copula model mat <- rvine_matrix_sim(3) vc <- vine_dist(list(stats_margin("norm")), pcs, mat) # simulate from the model x <- rvine(200, vc) pairs(x) # evaluate the density and cdf dvine(x[1, ], vc) pvine(x[1, ], vc)
Predictions of the density and distribution function for a vine copula model.
## S3 method for class 'vine' predict(object, newdata, what = "pdf", n_mc = 10^4, cores = 1, ...) ## S3 method for class 'vine' fitted(object, what = "pdf", n_mc = 10^4, cores = 1, ...)## S3 method for class 'vine' predict(object, newdata, what = "pdf", n_mc = 10^4, cores = 1, ...) ## S3 method for class 'vine' fitted(object, what = "pdf", n_mc = 10^4, cores = 1, ...)
object |
a |
newdata |
points where the fit shall be evaluated. |
what |
what to predict, either |
n_mc |
number of samples used for quasi Monte Carlo integration when
|
cores |
number of cores to use; if larger than one, computations are
done in parallel on |
... |
unused. |
fitted() and predict() have return values similar to dvine()
and pvine().
x <- sapply(1:5, function(i) rnorm(50)) fit <- vine(x, copula_controls = list(family_set = "par"), keep_data = TRUE) all.equal(predict(fit, x), fitted(fit), check.environment = FALSE)x <- sapply(1:5, function(i) rnorm(50)) fit <- vine(x, copula_controls = list(family_set = "par"), keep_data = TRUE) all.equal(predict(fit, x), fitted(fit), check.environment = FALSE)
Automated fitting and model selection for vine copula models with continuous or discrete data. Selection of the structure is performed using the algorithm of Dissmann et al. (2013).
vinecop( data, var_types = rep("c", NCOL(data)), family_set = "all", structure = NA, par_method = "mle", nonpar_method = "constant", mult = 1, selcrit = "aic", weights = numeric(), psi0 = 0.9, presel = TRUE, allow_rotations = TRUE, trunc_lvl = Inf, tree_crit = "tau", threshold = 0, keep_data = FALSE, vinecop_object = NULL, show_trace = FALSE, cores = 1, tree_algorithm = "mst_prim", conditioning_set = NULL )vinecop( data, var_types = rep("c", NCOL(data)), family_set = "all", structure = NA, par_method = "mle", nonpar_method = "constant", mult = 1, selcrit = "aic", weights = numeric(), psi0 = 0.9, presel = TRUE, allow_rotations = TRUE, trunc_lvl = Inf, tree_crit = "tau", threshold = 0, keep_data = FALSE, vinecop_object = NULL, show_trace = FALSE, cores = 1, tree_algorithm = "mst_prim", conditioning_set = NULL )
data |
a matrix or data.frame with at least two columns, containing the (pseudo-)observations for the two variables (copula data should have approximately uniform margins). More columns are required for discrete models, see Details. |
var_types |
variable types, a length d vector; e.g., |
family_set |
a character vector of families; see |
structure |
an |
par_method |
the estimation method for parametric models, either |
nonpar_method |
the estimation method for nonparametric models, either
|
mult |
multiplier for the smoothing parameters of nonparametric families. Values larger than 1 make the estimate more smooth, values less than 1 less smooth. |
selcrit |
criterion for family selection, either |
weights |
optional vector of weights for each observation. |
psi0 |
prior probability of a non-independence copula (only used for
|
presel |
whether the family set should be thinned out according to symmetry characteristics of the data. |
allow_rotations |
whether to allow rotations of the copula. |
trunc_lvl |
the truncation level of the vine copula; |
tree_crit |
the criterion for tree selection, one of |
threshold |
for thresholded vine copulas; |
keep_data |
whether the data should be stored (necessary for using
|
vinecop_object |
a |
show_trace |
logical; whether a trace of the fitting progress should be printed. |
cores |
number of cores to use; if more than 1, estimation of pair copulas within a tree is done in parallel. |
tree_algorithm |
The algorithm for building the spanning
tree ( |
conditioning_set |
variable indices or names to place at the end of the
vine order. The resulting model can be sampled conditionally with
|
If there are missing data (i.e., NA entries), incomplete observations are
discarded before fitting a pair-copula. This is done on a pair-by-pair basis
so that the maximal available information is used.
The dependence measures used to select trees (default: Kendall's tau) are corrected for ties (see wdm::wdm).
Let n be the number of observations and d the number of variables.
When at least one variable is discrete, two types of
"observations" are required in data: the first n x d block
contains realizations of . The second n x d
block contains realizations of . The minus indicates a
left-sided limit of the cdf. For, e.g., an integer-valued variable, it holds
. For continuous variables the left
limit and the cdf itself coincide. Respective columns can be omitted in the
second block.
Selection of the structure is performed using the algorithm of
Dissmann, J. F., E. C. Brechmann, C. Czado, and D. Kurowicka (2013).
Selecting and estimating regular vine copulae and application to
financial returns. Computational Statistics & Data Analysis, 59 (1),
52-69.
The dependence measure used to select trees (default: Kendall's tau) is
corrected for ties and can be changed using the tree_crit
argument, which can be set to "tau", "rho", "hoeffd", "mcor",
"cxi", or "joe", or to a custom function.
Both Prim's (default: "mst_prim") and Kruskal's ("mst_kruskal")
algorithms are available through tree_algorithm to set the
maximum spanning tree selection algorithm.
An alternative to the maximum spanning tree selection is to use random
spanning trees, which can be selected using tree_algorithm and
come in two flavors, both using Wilson's algorithm loop erased random walks:
"random_weighted"' generates a random spanning tree with probability proportional to the product of the weights (i.e., the dependence) of the edges in the tree.
"random_unweighted"' generates a random spanning tree uniformly over all spanning trees satisfying the proximity condition.
It is possible to fix the vine structure only in the first trees and select
the remaining ones automatically. To specify only the first k trees, supply
a k-truncated rvine_structure() or rvine_matrix(). All trees up to
trunc_lvl will then be selected automatically.
Objects inheriting from vinecop and vinecop_dist for vinecop(). In
addition to the entries provided by vinecop_dist(), there are:
threshold, the (set or estimated) threshold used for thresholding the
vine.
data (optionally, if keep_data = TRUE was used), the dataset that was
passed to vinecop().
controls, a list with fit controls that was passed to vinecop().
nobs, the number of observations that were used to fit the model.
Dissmann, J. F., E. C. Brechmann, C. Czado, and D. Kurowicka (2013). Selecting and estimating regular vine copulae and application to financial returns. Computational Statistics & Data Analysis, 59 (1), 52-69.
vinecop(), dvinecop(), pvinecop(), rvinecop(),
plot.vinecop(), contour.vinecop()
## simulate dummy data x <- rnorm(30) * matrix(1, 30, 5) + 0.5 * matrix(rnorm(30 * 5), 30, 5) u <- pseudo_obs(x) ## fit and select the model structure, family and parameters fit <- vinecop(u) summary(fit) plot(fit) contour(fit) ## select by BIC from one-parameter families fit <- vinecop(u, family_set = "onepar", selcrit = "bic") summary(fit) ## 1-truncated, Gaussian D-vine fit <- vinecop( u, structure = dvine_structure(1:5), family_set = "gauss", trunc_lvl = 1 ) plot(fit) contour(fit) ## Partial structure selection with only first tree specified structure <- rvine_structure(order = 1:5, list(rep(5, 4))) structure fit <- vinecop(u, structure = structure, family_set = "gauss") plot(fit) ## Model for discrete data x <- qpois(u, 1) # transform to Poisson margins # we require two types of observations (see Details) u_disc <- cbind(ppois(x, 1), ppois(x - 1, 1)) fit <- vinecop(u_disc, var_types = rep("d", 5)) ## Model for mixed data x <- qpois(u[, 1], 1) # transform first variable to Poisson margin # we require two types of observations (see Details) u_disc <- cbind(ppois(x, 1), u[, 2:5], ppois(x - 1, 1)) fit <- vinecop(u_disc, var_types = c("d", rep("c", 4)))## simulate dummy data x <- rnorm(30) * matrix(1, 30, 5) + 0.5 * matrix(rnorm(30 * 5), 30, 5) u <- pseudo_obs(x) ## fit and select the model structure, family and parameters fit <- vinecop(u) summary(fit) plot(fit) contour(fit) ## select by BIC from one-parameter families fit <- vinecop(u, family_set = "onepar", selcrit = "bic") summary(fit) ## 1-truncated, Gaussian D-vine fit <- vinecop( u, structure = dvine_structure(1:5), family_set = "gauss", trunc_lvl = 1 ) plot(fit) contour(fit) ## Partial structure selection with only first tree specified structure <- rvine_structure(order = 1:5, list(rep(5, 4))) structure fit <- vinecop(u, structure = structure, family_set = "gauss") plot(fit) ## Model for discrete data x <- qpois(u, 1) # transform to Poisson margins # we require two types of observations (see Details) u_disc <- cbind(ppois(x, 1), ppois(x - 1, 1)) fit <- vinecop(u_disc, var_types = rep("d", 5)) ## Model for mixed data x <- qpois(u[, 1], 1) # transform first variable to Poisson margin # we require two types of observations (see Details) u_disc <- cbind(ppois(x, 1), u[, 2:5], ppois(x - 1, 1)) fit <- vinecop(u_disc, var_types = c("d", rep("c", 4)))
Create custom vine copula models by specifying the pair-copulas, structure, and variable types.
vinecop_dist(pair_copulas, structure, var_types = rep("c", dim(structure)[1]))vinecop_dist(pair_copulas, structure, var_types = rep("c", dim(structure)[1]))
pair_copulas |
A nested list of ' |
structure |
an |
var_types |
variable types, a length d vector; e.g., |
Object of class vinecop_dist, i.e., a list containing:
pair_copulas, a list of lists. Each element of pair_copulas corresponds
to a tree, which is itself a list of bicop_dist() objects.
structure, a compressed representation of the vine structure, or an
object that can be coerced into one (see rvine_structure() and
as_rvine_structure()).
npars, a numeric with the number of (effective) parameters.
var_types the variable types.
rvine_structure(), rvine_matrix(), vinecop(),
plot.vinecop_dist(), contour.vinecop_dist(), dvinecop(),
pvinecop(), rvinecop()
# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine copula model vc <- vinecop_dist(pcs, mat) # visualization plot(vc) contour(vc) # simulate from the model pairs(rvinecop(200, vc))# specify pair-copulas bicop <- bicop_dist("bb1", 90, c(3, 2)) pcs <- list( list(bicop, bicop), # pair-copulas in first tree list(bicop) # pair-copulas in second tree ) # specify R-vine matrix mat <- matrix(c(1, 2, 3, 1, 2, 0, 1, 0, 0), 3, 3) # set up vine copula model vc <- vinecop_dist(pcs, mat) # visualization plot(vc) contour(vc) # simulate from the model pairs(rvinecop(200, vc))
Density, distribution function and random generation for the vine copula distribution.
dvinecop( u, vinecop, cores = 1, keep_all = FALSE, parameters = NULL, log = FALSE ) scores(u, vinecop, ...) ## S3 method for class 'vinecop_dist' scores(u, vinecop, step_wise = TRUE, cores = 1, parameters = NULL, ...) hessian(u, vinecop, ...) ## S3 method for class 'vinecop_dist' hessian(u, vinecop, step_wise = TRUE, cores = 1, parameters = NULL, ...) pvinecop(u, vinecop, n_mc = 10^4, cores = 1) rvinecop( n, vinecop, qrng = FALSE, cores = 1, u_cond = NULL, conditioning_set = NULL )dvinecop( u, vinecop, cores = 1, keep_all = FALSE, parameters = NULL, log = FALSE ) scores(u, vinecop, ...) ## S3 method for class 'vinecop_dist' scores(u, vinecop, step_wise = TRUE, cores = 1, parameters = NULL, ...) hessian(u, vinecop, ...) ## S3 method for class 'vinecop_dist' hessian(u, vinecop, step_wise = TRUE, cores = 1, parameters = NULL, ...) pvinecop(u, vinecop, n_mc = 10^4, cores = 1) rvinecop( n, vinecop, qrng = FALSE, cores = 1, u_cond = NULL, conditioning_set = NULL )
u |
matrix of evaluation points; must contain at least d columns, where d is the number of variables in the vine. More columns are required for discrete models, see Details. |
vinecop |
an object of class |
cores |
number of cores to use; if larger than one, computations are
done in parallel on |
keep_all |
if |
parameters |
optional observation-specific parameters for |
log |
if |
... |
unused. |
step_wise |
if |
n_mc |
number of samples used for quasi Monte Carlo integration. |
n |
number of observations. |
qrng |
if |
u_cond |
optional conditioning values for |
conditioning_set |
variable indices or names corresponding to the first
block of |
See vinecop() for the estimation and construction of vine copula
models.
The copula density is defined as joint density divided by marginal densities, irrespective of variable types.
When at least one variable is discrete, two types of
"observations" are required in u: the first block
contains realizations of . The second
block contains realizations of . The minus indicates a
left-sided limit of the cdf. For, e.g., an integer-valued variable, it holds
. For continuous variables the left
limit and the cdf itself coincide. Respective columns can be omitted in the
second block.
dvinecop() gives the density, pvinecop() gives the distribution function,
rvinecop() generates unconditional or conditional random deviates,
scores() gives the observation-wise score matrix, and hessian() gives the
average Hessian matrix.
If keep_all = TRUE, dvinecop() returns a list with entries pdf,
logpdf, pdf_edges, hfunc1, hfunc2, hfunc1_sub, and hfunc2_sub.
The _sub
entries contain h-functions evaluated at left-sided limits when the model
contains discrete variables; they are empty for fully continuous models.
The length of the result is determined by n for rvinecop(), and
the number of rows in u for the other functions.
The vinecop object is recycled to the length of the
result.
vinecop_dist(), vinecop(), plot.vinecop(), contour.vinecop()
## simulate dummy data x <- rnorm(30) * matrix(1, 30, 5) + 0.5 * matrix(rnorm(30 * 5), 30, 5) u <- pseudo_obs(x) ## fit a model vc <- vinecop(u, family_set = "clayton") # simulate from the model u <- rvinecop(100, vc) pairs(u) # evaluate the density and cdf dvinecop(u[1, ], vc) pvinecop(u[1, ], vc) # evaluate derivatives of the log-likelihood scores(u, vc) hessian(u, vc) # derivatives can also be computed for the full likelihood scores(u, vc, step_wise = FALSE) hessian(u, vc, step_wise = FALSE) ## Discrete models vc$var_types <- rep("d", 5) # convert model to discrete # with discrete data we need two types of observations (see Details) x <- qpois(u, 1) # transform to Poisson margins u_disc <- cbind(ppois(x, 1), ppois(x - 1, 1)) dvinecop(u_disc[1:5, ], vc) pvinecop(u_disc[1:5, ], vc) # simulated data always has uniform margins pairs(rvinecop(200, vc)) ## Conditional simulation vc_cond <- vinecop( u, family_set = "gaussian", conditioning_set = c(2, 4) ) u_cond <- c(0.25, 0.75) uc <- rvinecop( 100, vc_cond, u_cond = u_cond, conditioning_set = c(2, 4) ) stopifnot( isTRUE(all.equal(uc[, 2], rep(0.25, 100))), isTRUE(all.equal(uc[, 4], rep(0.75, 100))) )## simulate dummy data x <- rnorm(30) * matrix(1, 30, 5) + 0.5 * matrix(rnorm(30 * 5), 30, 5) u <- pseudo_obs(x) ## fit a model vc <- vinecop(u, family_set = "clayton") # simulate from the model u <- rvinecop(100, vc) pairs(u) # evaluate the density and cdf dvinecop(u[1, ], vc) pvinecop(u[1, ], vc) # evaluate derivatives of the log-likelihood scores(u, vc) hessian(u, vc) # derivatives can also be computed for the full likelihood scores(u, vc, step_wise = FALSE) hessian(u, vc, step_wise = FALSE) ## Discrete models vc$var_types <- rep("d", 5) # convert model to discrete # with discrete data we need two types of observations (see Details) x <- qpois(u, 1) # transform to Poisson margins u_disc <- cbind(ppois(x, 1), ppois(x - 1, 1)) dvinecop(u_disc[1:5, ], vc) pvinecop(u_disc[1:5, ], vc) # simulated data always has uniform margins pairs(rvinecop(200, vc)) ## Conditional simulation vc_cond <- vinecop( u, family_set = "gaussian", conditioning_set = c(2, 4) ) u_cond <- c(0.25, 0.75) uc <- rvinecop( 100, vc_cond, u_cond = u_cond, conditioning_set = c(2, 4) ) stopifnot( isTRUE(all.equal(uc[, 2], rep(0.25, 100))), isTRUE(all.equal(uc[, 4], rep(0.75, 100))) )
Predictions of the density and distribution function for a vine copula model.
## S3 method for class 'vinecop' predict(object, newdata, what = "pdf", n_mc = 10^4, cores = 1, ...) ## S3 method for class 'vinecop' fitted(object, what = "pdf", n_mc = 10^4, cores = 1, ...)## S3 method for class 'vinecop' predict(object, newdata, what = "pdf", n_mc = 10^4, cores = 1, ...) ## S3 method for class 'vinecop' fitted(object, what = "pdf", n_mc = 10^4, cores = 1, ...)
object |
a |
newdata |
points where the fit shall be evaluated. |
what |
what to predict, either |
n_mc |
number of samples used for quasi Monte Carlo integration when
|
cores |
number of cores to use; if larger than one, computations are
done in parallel on |
... |
unused. |
fitted() can only be called if the model was fit with the
keep_data = TRUE option.
When at least one variable is discrete, two types of
"observations" are required in newdata: the first block
contains realizations of . The second
block contains realizations of . The minus indicates a
left-sided limit of the cdf. For, e.g., an integer-valued variable, it holds
. For continuous variables the left
limit and the cdf itself coincide. Respective columns can be omitted in the
second block.
fitted() and predict() have return values similar to dvinecop()
and pvinecop().
u <- sapply(1:5, function(i) runif(50)) fit <- vinecop(u, family_set = "par", keep_data = TRUE) all.equal(predict(fit, u), fitted(fit), check.environment = FALSE)u <- sapply(1:5, function(i) runif(50)) fit <- vinecop(u, family_set = "par", keep_data = TRUE) all.equal(predict(fit, u), fitted(fit), check.environment = FALSE)
Marks a numeric vector as continuous with an atom at zero. The class is
preserved when stored in or subset from a data frame and is detected by
vine(). Alternatively, use var_types = "zi".
zero_inflated(x)zero_inflated(x)
x |
a numeric vector. |
x with an additional zero_inflated class.