vinecopulib::Vinecop class

A class for vine copula models.

A vine copula model is characterized by its structure (see RVineStructure objects) and the pair-copulas.

Public static functions

static auto make_pair_copula_store(const size_t d, const size_t trunc_lvl = std::numeric_limits<size_t>::max()) -> std::vector<std::vector<Bicop>>
Initializes object for storing pair copulas.

Constructors, destructors, conversion operators

Vinecop(size_t d) explicit
Instantiates a D-vine with all pair-copulas set to independence.
Vinecop(const RVineStructure& structure, const std::vector<std::vector<Bicop>>& pair_copulas = {}, const std::vector<std::string>& var_types = {}) explicit
Instantiates an arbitrary vine copula model.
Vinecop(const Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>& matrix, const std::vector<std::vector<Bicop>>& pair_copulas = {}, const std::vector<std::string>& var_types = {}) explicit
Instantiates an arbitrary vine copula model.
Vinecop(const Eigen::MatrixXd& data, const RVineStructure& structure = RVineStructure(), const std::vector<std::string>& var_types = {}, const FitControlsVinecop& controls = FitControlsVinecop()) explicit
Instantiates from data.
Vinecop(const Eigen::MatrixXd& data, const Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>& matrix = Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>(), const std::vector<std::string>& var_types = {}, const FitControlsVinecop& controls = FitControlsVinecop()) explicit
Instantiates from data.
Vinecop(const std::string& filename, const bool check = true) explicit
Instantiates from a JSON file.
Vinecop(const nlohmann::json& input, const bool check = true) explicit
Instantiates from a nlohmann::json object.

Public functions

auto to_json() const -> nlohmann::json
Converts the copula into a nlohmann::json object.
void to_file(const std::string& filename) const
Writes the copula object into a JSON file.
void select(const Eigen::MatrixXd& data, const FitControlsVinecop& controls = FitControlsVinecop())
Automatically fits and selects a vine copula model.
auto select_all(const Eigen::MatrixXd& data, const FitControlsVinecop& controls = FitControlsVinecop()) -> DEPRECATED void
Automatically fits and selects a vine copula model.
auto select_families(const Eigen::MatrixXd& data, const FitControlsVinecop& controls = FitControlsVinecop()) -> DEPRECATED void
Automatically selects all pair-copula families and fits all. parameters.
auto pdf(Eigen::MatrixXd u, const size_t num_threads = 1) const -> Eigen::VectorXd
Evaluates the copula density.
auto cdf(const Eigen::MatrixXd& u, const size_t N = 1e4, const size_t num_threads = 1, std::vector<int> seeds = std::vector<int>()) const -> Eigen::VectorXd
Evaluates the copula distribution.
auto simulate(const size_t n, const bool qrng = false, const size_t num_threads = 1, const std::vector<int>& seeds = std::vector<int>()) const -> Eigen::MatrixXd
Simulates from a vine copula model, see inverse_rosenblatt().
auto rosenblatt(const Eigen::MatrixXd& u, const size_t num_threads = 1) const -> Eigen::MatrixXd
Evaluates the Rosenblatt transform for a vine copula model.
auto inverse_rosenblatt(const Eigen::MatrixXd& u, const size_t num_threads = 1) const -> Eigen::MatrixXd
Evaluates the inverse Rosenblatt transform.
auto get_npars() const -> double
Returns sum of the number of parameters for all pair copulas (see. Bicop::get_npars()).
auto loglik(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const size_t num_threads = 1) const -> double
Evaluates the log-likelihood.
auto aic(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const size_t num_threads = 1) const -> double
Evaluates the Akaike information criterion (AIC).
auto bic(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const size_t num_threads = 1) const -> double
Evaluates the Bayesian information criterion (BIC).
auto mbicv(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const double psi0 = 0.9, const size_t num_threads = 1) const -> double
Evaluates the modified Bayesian information criterion for vines (mBICV).
void truncate(size_t trunc_lvl)
Truncates the vine copula model.
auto str() const -> std::string
Summarizes the model into a string (can be used for printing).

Getters and setters

auto get_pair_copula(size_t tree, size_t edge) const -> Bicop
Gets a pair copula.
auto get_family(size_t tree, size_t edge) const -> BicopFamily
Gets the family of a pair copula.
auto get_rotation(size_t tree, size_t edge) const -> int
Gets the rotation of a pair copula.
auto get_parameters(size_t tree, size_t edge) const -> Eigen::MatrixXd
Gets the parameters of a pair copula.
auto get_tau(size_t tree, size_t edge) const -> double
Gets the Kendall's $ tau $ of a pair copula.
auto get_all_pair_copulas() const -> std::vector<std::vector<Bicop>>
Gets all pair copulas.
auto get_all_families() const -> std::vector<std::vector<BicopFamily>>
Gets the families of all pair copulas.
auto get_all_rotations() const -> std::vector<std::vector<int>>
Gets the rotations of all pair copulas.
auto get_all_parameters() const -> std::vector<std::vector<Eigen::MatrixXd>>
Gets the parameters of all pair copulas.
auto get_all_taus() const -> std::vector<std::vector<double>>
Gets the Kendall's $ tau $ s of all pair copulas.
auto get_dim() const -> size_t
Gets the dimension of the vine copula model.
auto get_order() const -> std::vector<size_t>
Gets the order vector of the vine copula model.
auto get_rvine_structure() const -> RVineStructure
Gets the structure matrix of the vine copula model.
auto get_matrix() const -> Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>
Gets the structure matrix of the vine copula model.
auto get_struct_array(bool natural_order = false) const -> TriangularArray<size_t>
Gets the above diagonal coefficients of the vine copula model.
auto get_threshold() const -> double
Gets the threshold.
auto get_loglik() const -> double
Gets the log-likelihood (throws an error if model has not been. fitted to data).
auto get_nobs() const -> size_t
Gets the number of observations used for the fit.
auto get_aic() const -> double
Gets the AIC.
auto get_bic() const -> double
Gets the BIC.
auto get_mbicv(const double psi0 = 0.9) const -> double
Gets the log-likelihood.
void set_all_pair_copulas(const std::vector<std::vector<Bicop>>& pair_copulas)
Sets all pair-copulas.
void set_var_types(const std::vector<std::string>& var_types)
Sets variable types.
auto get_var_types() const -> std::vector<std::string>
Gets the variable types.

Function documentation

static std::vector<std::vector<Bicop>> vinecopulib::Vinecop::make_pair_copula_store(const size_t d, const size_t trunc_lvl = std::numeric_limits<size_t>::max())

Initializes object for storing pair copulas.

Parameters
d Dimension of the vine copula.
trunc_lvl A truncation level (optional).
Returns A nested list such that pc_store[t][e] contains a Bicop. object for the pair copula corresponding to tree t and edge e.

vinecopulib::Vinecop::Vinecop(size_t d) explicit

Instantiates a D-vine with all pair-copulas set to independence.

Parameters
d The dimension (= number of variables) of the model.

vinecopulib::Vinecop::Vinecop(const RVineStructure& structure, const std::vector<std::vector<Bicop>>& pair_copulas = {}, const std::vector<std::string>& var_types = {}) explicit

Instantiates an arbitrary vine copula model.

Parameters
structure An RVineStructure object specifying the vine structure.
pair_copulas Bicop objects specifying the pair-copulas, namely a nested list such that pc_store[t][e] contains a Bicop object for the pair copula corresponding to tree t and edge e.
var_types Strings specifying the types of the variables, e.g., ("c", "d") means first variable continuous, second discrete. If empty, then all variables are set as continuous.

vinecopulib::Vinecop::Vinecop(const Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>& matrix, const std::vector<std::vector<Bicop>>& pair_copulas = {}, const std::vector<std::string>& var_types = {}) explicit

Instantiates an arbitrary vine copula model.

Parameters
matrix An R-vine matrix specifying the vine structure.
pair_copulas Bicop objects specifying the pair-copulas, namely a nested list such that pc_store[t][e] contains a Bicop object for the pair copula corresponding to tree t and edge e.
var_types Strings specifying the types of the variables, e.g., ("c", "d") means first variable continuous, second discrete. If empty, then all variables are set as continuous.

vinecopulib::Vinecop::Vinecop(const Eigen::MatrixXd& data, const RVineStructure& structure = RVineStructure(), const std::vector<std::string>& var_types = {}, const FitControlsVinecop& controls = FitControlsVinecop()) explicit

Instantiates from data.

Parameters
data An $ n \times d $ matrix of observations.
structure An RVineStructure object specifying the vine structure. If empty, then it is selected as part of the fit.
var_types Strings specifying the types of the variables, e.g., ("c", "d") means first variable continuous, second discrete. If empty, then all variables are set as continuous.
controls See FitControlsVinecop().

Equivalent to creating a default Vinecop() and then selecting the model using select().

vinecopulib::Vinecop::Vinecop(const Eigen::MatrixXd& data, const Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>& matrix = Eigen::Matrix<size_t, Eigen::Dynamic, Eigen::Dynamic>(), const std::vector<std::string>& var_types = {}, const FitControlsVinecop& controls = FitControlsVinecop()) explicit

Instantiates from data.

Parameters
data An $ n \times d $ matrix of observations.
matrix Either an empty matrix (default) or an R-vine structure matrix, see select(). If empty, then it is selected as part of the fit.
var_types Strings specifying the types of the variables, e.g., ("c", "d") means first variable continuous, second discrete. If empty, then all variables are set as continuous.
controls See FitControlsVinecop().

Equivalent to creating a default Vinecop() and then selecting the model using select().

vinecopulib::Vinecop::Vinecop(const std::string& filename, const bool check = true) explicit

Instantiates from a JSON file.

Parameters
filename The name of the JSON file to read.
check Whether to check if the "structure" node of the input represents a valid R-vine structure.

The input file contains 2 attributes : "structure" for the vine structure, which itself contains attributes "array" for the structure triangular array and "order" for the order vector, and "pair copulas". "pair copulas" contains a list of attributes for the trees ("tree1", "tree2", etc), each containing a list of attributes for the edges ("pc1", "pc2", etc). See the corresponding method of Bicop objects for the encoding of pair-copulas.

vinecopulib::Vinecop::Vinecop(const nlohmann::json& input, const bool check = true) explicit

Instantiates from a nlohmann::json object.

Parameters
input The nlohmann::json object to convert from (see to_json() for the structure of the input).
check Whether to check if the "structure" node represents a valid R-vine structure.

nlohmann::json vinecopulib::Vinecop::to_json() const

Converts the copula into a nlohmann::json object.

Returns the nlohmann::json object containing the copula.

The nlohmann::json object contains two nodes : "structure" for the vine structure, which itself contains nodes "array" for the structure triangular array and "order" for the order vector, and "pair copulas". The former two encode the R-Vine structure and the latter is a list of child nodes for the trees ("tree1", "tree2", etc), each containing a list of child nodes for the edges ("pc1", "pc2", etc). See Bicop::to_json() for the encoding of pair-copulas.

void vinecopulib::Vinecop::to_file(const std::string& filename) const

Writes the copula object into a JSON file.

Parameters
filename The name of the JSON file to write.

The output file contains 2 attributes : "structure" for the vine structure, which itself contains attributes "array" for the structure triangular array and "order" for the order vector, and "pair copulas". "pair copulas" contains a list of attributes for the trees ("tree1", "tree2", etc), each containing a list of attributes for the edges ("pc1", "pc2", etc). See the corresponding method of Bicop objects for the encoding of pair-copulas.

void vinecopulib::Vinecop::select(const Eigen::MatrixXd& data, const FitControlsVinecop& controls = FitControlsVinecop())

Automatically fits and selects a vine copula model.

Parameters
data $ n \times (d + k) $ or $ n \times 2d $ matrix of observations, where $ k $ is the number of discrete variables.
controls The controls to the algorithm (see FitControlsVinecop).

select() behaves differently depending on its current truncation level and the truncation level specified in the controls, respectively called trunc_lvl and controls.trunc_lvl in what follows. Essentially, controls.trunc_lvl defines the object's truncation level after calling select():

  • If controls.trunc_lvl <= trunc_lvl, the families and parameters for all pairs in trees smaller or equal to controls.trunc_lvl are selected, using the current structure.
  • If controls.trunc_lvl > trunc_lvl, select() behaves as above for all trees that are smaller or equal to trunc_lvl, and then it selects the structure for higher trees along with the families and parameters. This includes the case where trunc_lvl = 0, namely where the structure is fully unspecified.

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 (see the wdm library).

When at least one variable is discrete, two types of "observations" are required: the first $ n \times d $ block contains realizations of $ F_Y(Y), F_X(X) $ ; the second $ n \times d $ block contains realizations of $ F_Y(Y^-), F_X(X^-), ... $ . The minus indicates a left-sided limit of the cdf. For continuous variables the left limit and the cdf itself coincide. For, e.g., an integer-valued variable, it holds $ F_Y(Y^-) = F_Y(Y - 1) $ . Continuous variables in the second block can be omitted.

If there are missing data (i.e., NaN 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.

DEPRECATED void vinecopulib::Vinecop::select_all(const Eigen::MatrixXd& data, const FitControlsVinecop& controls = FitControlsVinecop())

Automatically fits and selects a vine copula model.

Parameters
data $ n \times (d + k) $ or $ n \times 2d $ matrix of observations, where $ k $ is the number of discrete variables.
controls The controls to the algorithm (see FitControlsVinecop).

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.

When at least one variable is discrete, two types of "observations" are required: the first $ n \times d $ block contains realizations of $ F_Y(Y), F_X(X) $ ; the second $ n \times d $ block contains realizations of $ F_Y(Y^-), F_X(X^-), ... $ . The minus indicates a left-sided limit of the cdf. For continuous variables the left limit and the cdf itself coincide. For, e.g., an integer-valued variable, it holds $ //! F_Y(Y^-) = F_Y(Y - 1) $ . Continuous variables in the second block can be omitted.

DEPRECATED void vinecopulib::Vinecop::select_families(const Eigen::MatrixXd& data, const FitControlsVinecop& controls = FitControlsVinecop())

Automatically selects all pair-copula families and fits all. parameters.

Parameters
data $ n \times (d + k) $ or $ n \times 2d $ matrix of observations, where $ k $ is the number of discrete variables.
controls The controls to the algorithm (see FitControlsVinecop).

When at least one variable is discrete, two types of "observations" are required: the first $ n \times d $ block contains realizations of $ F_Y(Y), F_X(X) $ ; the second $ n \times d $ block contains realizations of $ F_Y(Y^-), F_X(X^-), ... $ . The minus indicates a left-sided limit of the cdf. For continuous variables the left limit and the cdf itself coincide. For, e.g., an integer-valued variable, it holds $ //! F_Y(Y^-) = F_Y(Y - 1) $ . Continuous variables in the second block can be omitted.

If there are missing data (i.e., NaN 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.

Eigen::VectorXd vinecopulib::Vinecop::pdf(Eigen::MatrixXd u, const size_t num_threads = 1) const

Evaluates the copula density.

Parameters
u An $ n \times (d + k) $ or $ n \times 2d $ matrix of evaluation points, where $ k $ is the number of discrete variables (see select()).
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The copula density is defined as joint density divided by marginal densities, irrespective of variable types.

Eigen::VectorXd vinecopulib::Vinecop::cdf(const Eigen::MatrixXd& u, const size_t N = 1e4, const size_t num_threads = 1, std::vector<int> seeds = std::vector<int>()) const

Evaluates the copula distribution.

Parameters
u An $ n \times (d + k) $ or $ n \times 2d $ matrix of evaluation points, where $ k $ is the number of discrete variables (see select()).
N Integer for the number of quasi-random numbers to draw to evaluate the distribution (default: 1e4).
num_threads The number of threads to use for computations; if greater than 1, the function will generate n samples concurrently in num_threads batches.
seeds Seeds to scramble the quasi-random numbers; if empty (default), the random number quasi-generator is seeded randomly.

Because no closed-form expression is available, the distribution is estimated numerically using Monte Carlo integration.

Eigen::MatrixXd vinecopulib::Vinecop::simulate(const size_t n, const bool qrng = false, const size_t num_threads = 1, const std::vector<int>& seeds = std::vector<int>()) const

Simulates from a vine copula model, see inverse_rosenblatt().

Parameters
n Number of observations.
qrng Set to true for quasi-random numbers.
num_threads The number of threads to use for computations; if greater than 1, the function will generate n samples concurrently in num_threads batches.
seeds Seeds of the random number generator; if empty (default), the random number generator is seeded randomly.
Returns An $ n \times d $ matrix of samples from the copula model.

Simulated data is always a continous $ n \times d $ matrix.

Eigen::MatrixXd vinecopulib::Vinecop::rosenblatt(const Eigen::MatrixXd& u, const size_t num_threads = 1) const

Evaluates the Rosenblatt transform for a vine copula model.

Parameters
u An $ n \times d $ or $ n \times 2d $ matrix of evaluation points.
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The Rosenblatt transform converts data from this model into independent uniform variates. Only works for continuous data.

The Rosenblatt transform (Rosenblatt, 1952) $ U = T(V) $ of a random vector $ V = (V_1,\ldots,V_d) ~ F $ is defined as

\[ U_1= F(V_1), U_{2} = F(V_{2}|V_1), \ldots, U_d =F(V_d|V_1,\ldots,V_{d-1}), \]

where $ F(v_k|v_1,\ldots,v_{k-1}) $ is the conditional distribution of $ V_k $ given $ V_1 \ldots, V_{k-1}, k = 2,\ldots,d $ . The vector $ U = (U_1, \dots, U_d) $ then contains independent standard uniform variables. The inverse operation

\[V_1 = F^{-1}(U_1), V_{2} = F^{-1}(U_2|U_1), \ldots, V_d =F^{-1}(U_d|U_1,\ldots,U_{d-1}) \]

can be used to simulate from a distribution. For any copula $ F $ , if $ U$ is a vector of independent random variables, $ V = T^{-1}(U) $ has distribution $ F $ .

The formulas above assume a vine copula model with order $ d, \dots, 1 $ . More generally, Vinecop::rosenblatt() returns the variables

\[ U_{M[d - j, j]}= F(V_{M[d - j, j]} | V_{M[d - j - 1, j - 1]}, \dots, V_{M[0, 0]}), \]

where $ M $ is the structure matrix. Similarly, Vinecop::inverse_rosenblatt() returns

\[ V_{M[d - j, j]}= F^{-1}(U_{M[d - j, j]} | U_{M[d - j - 1, j - 1]}, \dots, U_{M[0, 0]}). \]

Eigen::MatrixXd vinecopulib::Vinecop::inverse_rosenblatt(const Eigen::MatrixXd& u, const size_t num_threads = 1) const

Evaluates the inverse Rosenblatt transform.

Parameters
u An $ n \times d $ matrix of evaluation points.
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The inverse Rosenblatt transform can be used for simulation: the function applied to independent uniform variates resembles simulated data from the vine copula model.

If the problem is too large, it is split recursively into halves (w.r.t. $ n $ , the number of observations). "Too large" means that the required memory will exceed 1 GB. An examplary configuration requiring less than 1 GB is $ n = 1000 $ , $d = 200$ .

Only works for continous models.

The Rosenblatt transform (Rosenblatt, 1952) $ U = T(V) $ of a random vector $ V = (V_1,\ldots,V_d) ~ F $ is defined as

\[ U_1= F(V_1), U_{2} = F(V_{2}|V_1), \ldots, U_d =F(V_d|V_1,\ldots,V_{d-1}), \]

where $ F(v_k|v_1,\ldots,v_{k-1}) $ is the conditional distribution of $ V_k $ given $ V_1 \ldots, V_{k-1}, k = 2,\ldots,d $ . The vector $ U = (U_1, \dots, U_d) $ then contains independent standard uniform variables. The inverse operation

\[V_1 = F^{-1}(U_1), V_{2} = F^{-1}(U_2|U_1), \ldots, V_d =F^{-1}(U_d|U_1,\ldots,U_{d-1}) \]

can be used to simulate from a distribution. For any copula $ F $ , if $ U$ is a vector of independent random variables, $ V = T^{-1}(U) $ has distribution $ F $ .

The formulas above assume a vine copula model with order $ d, \dots, 1 $ . More generally, Vinecop::rosenblatt() returns the variables

\[ U_{M[d - j, j]}= F(V_{M[d - j, j]} | V_{M[d - j - 1, j - 1]}, \dots, V_{M[0, 0]}), \]

where $ M $ is the structure matrix. Similarly, Vinecop::inverse_rosenblatt() returns

\[ V_{M[d - j, j]}= F^{-1}(U_{M[d - j, j]} | U_{M[d - j - 1, j - 1]}, \dots, U_{M[0, 0]}). \]

double vinecopulib::Vinecop::loglik(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const size_t num_threads = 1) const

Evaluates the log-likelihood.

Parameters
u An $ n \times (d + k) $ or $ n \times 2d $ matrix of evaluation points, where $ k $ is the number of discrete variables (see select()).
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The log-likelihood is defined as

\[ \mathrm{loglik} = \sum_{i = 1}^n \log c(U_{1, i}, ..., U_{d, i}), \]

where $ c $ is the copula density, see Vinecop::pdf().

double vinecopulib::Vinecop::aic(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const size_t num_threads = 1) const

Evaluates the Akaike information criterion (AIC).

Parameters
u An $ n \times (d + k) $ or $ n \times 2d $ matrix of evaluation points, where $ k $ is the number of discrete variables (see select()).
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The AIC is defined as

\[ \mathrm{AIC} = -2\, \mathrm{loglik} + 2 p, \]

where $ \mathrm{loglik} $ is the log-liklihood (see loglik()) and $ p $ is the (effective) number of parameters of the model. The AIC is a consistent model selection criterion even for nonparametric models.

double vinecopulib::Vinecop::bic(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const size_t num_threads = 1) const

Evaluates the Bayesian information criterion (BIC).

Parameters
u An $ n \times (d + k) $ or $ n \times 2d $ matrix of evaluation points, where $ k $ is the number of discrete variables (see select()).
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The BIC is defined as

\[ \mathrm{BIC} = -2\, \mathrm{loglik} + \log(n) p, \]

where $ \mathrm{loglik} $ is the log-liklihood (see loglik()) and $ p $ is the (effective) number of parameters of the model. The BIC is a consistent model selection criterion for nonparametric models.

double vinecopulib::Vinecop::mbicv(const Eigen::MatrixXd& u = Eigen::MatrixXd(), const double psi0 = 0.9, const size_t num_threads = 1) const

Evaluates the modified Bayesian information criterion for vines (mBICV).

Parameters
u An $ n \times (d + k) $ or $ n \times 2d $ matrix of evaluation points, where $ k $ is the number of discrete variables (see select()).
psi0 Baseline prior probability of a non-independence copula.
num_threads The number of threads to use for computations; if greater than 1, the function will be applied concurrently to num_threads batches of u.

The mBICV is defined as

\[ \mathrm{mBICV} = -2\, \mathrm{loglik} + \log(n) p, - 2 * \sum_{t=1}^(d - 1) \{q_t \log(\psi_0^t) - (d - t - q_t) \log(1 -\psi_0^t)\},\]

where $ \mathrm{loglik} $ is the log-liklihood, $ p $ is the (effective) number of parameters of the model, $ t $ is the tree level, $ \psi_0 $ is the prior probability of having a non-independence copula in the first tree, and $ q_t $ is the number of non-independence copulas in tree $ t $ ; The vBIC is a consistent model selection criterion for parametric sparse vine copula models when $ d = o(\sqrt{n \log n})$ .

void vinecopulib::Vinecop::truncate(size_t trunc_lvl)

Truncates the vine copula model.

Parameters
trunc_lvl The truncation level.

If the model is already truncated at a level less than trunc_lvl, the function does nothing.

Bicop vinecopulib::Vinecop::get_pair_copula(size_t tree, size_t edge) const

Gets a pair copula.

Parameters
tree Tree index (starting with 0).
edge Edge index (starting with 0).

BicopFamily vinecopulib::Vinecop::get_family(size_t tree, size_t edge) const

Gets the family of a pair copula.

Parameters
tree Tree index (starting with 0).
edge Edge index (starting with 0).

int vinecopulib::Vinecop::get_rotation(size_t tree, size_t edge) const

Gets the rotation of a pair copula.

Parameters
tree Tree index (starting with 0).
edge Edge index (starting with 0).

Eigen::MatrixXd vinecopulib::Vinecop::get_parameters(size_t tree, size_t edge) const

Gets the parameters of a pair copula.

Parameters
tree Tree index (starting with 0).
edge Edge index (starting with 0).

double vinecopulib::Vinecop::get_tau(size_t tree, size_t edge) const

Gets the Kendall's $ tau $ of a pair copula.

Parameters
tree Tree index (starting with 0).
edge Edge index (starting with 0).

std::vector<std::vector<Bicop>> vinecopulib::Vinecop::get_all_pair_copulas() const

Gets all pair copulas.

Returns a nested std::vector with entry [t][e] corresponding to edge e in tree t.

std::vector<std::vector<BicopFamily>> vinecopulib::Vinecop::get_all_families() const

Gets the families of all pair copulas.

Returns a nested std::vector with entry [t][e] corresponding to edge e in tree t.

std::vector<std::vector<int>> vinecopulib::Vinecop::get_all_rotations() const

Gets the rotations of all pair copulas.

Returns a nested std::vector with entry [t][e] corresponding to edge e in tree t.

std::vector<std::vector<Eigen::MatrixXd>> vinecopulib::Vinecop::get_all_parameters() const

Gets the parameters of all pair copulas.

Returns a nested std::vector with entry [t][e] corresponding to edge e in tree t.

std::vector<std::vector<double>> vinecopulib::Vinecop::get_all_taus() const

Gets the Kendall's $ tau $ s of all pair copulas.

Returns a nested std::vector with entry [t][e] corresponding to edge e in tree t.

TriangularArray<size_t> vinecopulib::Vinecop::get_struct_array(bool natural_order = false) const

Gets the above diagonal coefficients of the vine copula model.

Parameters
natural_order Whether indices correspond to natural order.

double vinecopulib::Vinecop::get_threshold() const

Gets the threshold.

Usually zero except select_threshold == TRUE in FitControlsVinecop()).

size_t vinecopulib::Vinecop::get_nobs() const

Gets the number of observations used for the fit.

The function throws an error if model has not been fitted to data.

double vinecopulib::Vinecop::get_aic() const

Gets the AIC.

The function throws an error if model has not been fitted to data.

double vinecopulib::Vinecop::get_bic() const

Gets the BIC.

The function throws an error if model has not been fitted to data.

double vinecopulib::Vinecop::get_mbicv(const double psi0 = 0.9) const

Gets the log-likelihood.

The function throws an error if model has not been fitted to data.

void vinecopulib::Vinecop::set_all_pair_copulas(const std::vector<std::vector<Bicop>>& pair_copulas)

Sets all pair-copulas.

Parameters
pair_copulas A vector of pair-copulas that has to be consistent with the current structure (see Vinecop()).

void vinecopulib::Vinecop::set_var_types(const std::vector<std::string>& var_types)

Sets variable types.

Parameters
var_types A vector specifying the types of the variables, e.g., {"c", "d"} means first varible continuous, second discrete.