Release 0.7.3.1.0

.github/workflows/R-CMD-check.yaml

@@ -2,9 +2,7 @@
 # Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help
 on:
   push:
-    branches: '*'
   pull_request:
-    branches: '*'
 
 name: R-CMD-check
 

.github/workflows/format-suggest.yaml

@@ -0,0 +1,30 @@
+# Workflow derived from https://github.com/posit-dev/setup-air/tree/main/examples
+on:
+  push:
+  pull_request:
+
+name: format-suggest.yaml
+
+permissions: read-all
+
+jobs:
+  format-suggest:
+    name: format-suggest
+    runs-on: ubuntu-latest
+    permissions:
+      pull-requests: write
+    steps:
+      - uses: actions/checkout@v4
+
+      - name: Install
+        uses: posit-dev/setup-air@v1
+
+      - name: Format
+        run: air format .
+
+      - name: Suggest
+        uses: reviewdog/action-suggester@v1
+        with:
+          level: error
+          fail_level: error
+          tool_name: air

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: rvinecopulib
 Type: Package
 Title: High Performance Algorithms for Vine Copula Modeling
-Version: 0.7.2.1.0
+Version: 0.7.3.1.0
 Authors@R: c(
         person("Thomas", "Nagler",, "[email protected]", role = c("aut", "cre")),
         person("Thibault", "Vatter",, "[email protected]", role = c("aut"))

NEWS.md

@@ -1,3 +1,16 @@
+# rvinecopulib 0.7.3.1.0
+
+### NEW FEATURES
+
+* Allow for random spanning trees as alternatives to the MST-based structure selection using
+  `tree_algorithm` in `vine` and `vinecop`. Options are `"mst_prim"`, `"mst_kruskal"`,
+  `"random_weighted"` or `"random_unweighted"` (#307).
+
+### BUG FIXES
+
+* Decouple edge insertion from criterion computation fix randomness
+  issues in structure selection when using multiple threads ([#640](https://github.com/vinecopulib/vinecopulib/pull/640))
+
 # rvinecopulib 0.7.2.1.0
 
 BUG FIX

R/RcppExports.R

@@ -89,12 +89,12 @@ vinecop_cdf_cpp <- function(u, vinecop_r, N, cores, seeds) {
     .Call(`_rvinecopulib_vinecop_cdf_cpp`, u, vinecop_r, N, cores, seeds)
 }
 
-vinecop_select_cpp <- function(data, structure, family_set, par_method, nonpar_method, mult, truncation_level, tree_criterion, threshold, selection_criterion, weights, psi0, select_truncation_level, select_threshold, preselect_families, select_families, allow_rotations, show_trace, num_threads, var_types) {
-    .Call(`_rvinecopulib_vinecop_select_cpp`, data, structure, family_set, par_method, nonpar_method, mult, truncation_level, tree_criterion, threshold, selection_criterion, weights, psi0, select_truncation_level, select_threshold, preselect_families, select_families, allow_rotations, show_trace, num_threads, var_types)
+vinecop_select_cpp <- function(data, structure, family_set, par_method, nonpar_method, mult, truncation_level, tree_criterion, threshold, selection_criterion, weights, psi0, select_truncation_level, select_threshold, preselect_families, select_families, allow_rotations, show_trace, num_threads, var_types, tree_algorithm, seeds) {
+    .Call(`_rvinecopulib_vinecop_select_cpp`, data, structure, family_set, par_method, nonpar_method, mult, truncation_level, tree_criterion, threshold, selection_criterion, weights, psi0, select_truncation_level, select_threshold, preselect_families, select_families, allow_rotations, show_trace, num_threads, var_types, tree_algorithm, seeds)
 }
 
-vinecop_fit_cpp <- function(data, vinecop_r, par_method, nonpar_method, mult, weights, show_trace, num_threads) {
-    .Call(`_rvinecopulib_vinecop_fit_cpp`, data, vinecop_r, par_method, nonpar_method, mult, weights, show_trace, num_threads)
+vinecop_fit_cpp <- function(data, vinecop_r, par_method, nonpar_method, mult, weights, show_trace, num_threads, tree_algorithm, seeds) {
+    .Call(`_rvinecopulib_vinecop_fit_cpp`, data, vinecop_r, par_method, nonpar_method, mult, weights, show_trace, num_threads, tree_algorithm, seeds)
 }
 
 fit_margins_cpp <- function(data, xmin, xmax, type, mult, bw, deg, weights, num_threads) {

R/as_rvine_structure.R

@@ -118,7 +118,7 @@ as_rvine_matrix.rvine_structure <- function(x, ..., validate = FALSE) {
   matrix <- matrix[d:1, ]
   for (i in seq_len(min(trunc_lvl, d - 1))) {
     newrow <- order[x[["struct_array"]][[i]]]
-    matrix[i, 1:length(newrow)] <- newrow
+    matrix[i, seq_along(newrow)] <- newrow
   }
 
   class(matrix) <- c("rvine_matrix", class(matrix))
@@ -182,13 +182,14 @@ as_rvine_structure.rvine_matrix <- function(x, ..., validate = FALSE) {
   struct_array <- lapply(1:(d - 1), function(i) order[x[i, 1:(d - i)]])
 
   # create and return x
-  structure(list(
-    order = diag(x[d:1, ]),
-    struct_array = struct_array,
-    d = d,
-    trunc_lvl = dim(x)[2]
-  ),
-  class = c("rvine_structure", "list")
+  structure(
+    list(
+      order = diag(x[d:1, ]),
+      struct_array = struct_array,
+      d = d,
+      trunc_lvl = dim(x)[2]
+    ),
+    class = c("rvine_structure", "list")
   )
 }
 

R/bicop.R

@@ -121,23 +121,37 @@
 #' fit_disc <- bicop(udisc, var_types = c("d", "d"))
 #' summary(fit_disc)
 #' @export
-bicop <- function(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 <- function(
+  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
+) {
   assert_that(
     is.character(family_set),
     is.string(par_method),
     is.string(nonpar_method),
-    is.number(mult), mult > 0,
+    is.number(mult),
+    mult > 0,
     is.string(selcrit),
     is.numeric(weights),
-    is.number(psi0), psi0 > 0, psi0 < 1,
+    is.number(psi0),
+    psi0 > 0,
+    psi0 < 1,
     is.flag(presel),
     is.flag(allow_rotations),
     is.flag(keep_data),
-    is.number(cores), cores > 0,
+    is.number(cores),
+    cores > 0,
     correct_var_types(var_types)
   )
 
@@ -273,8 +287,12 @@ as.bicop <- function(object, check = TRUE) {
 #' cop_disc <- bicop_dist("bb8", 0, c(2, 0.5), var_types = c("d", "d"))
 #' cop_disc
 #'
-bicop_dist <- function(family = "indep", rotation = 0, parameters = numeric(0),
-                       var_types = c("c", "c")) {
+bicop_dist <- function(
+  family = "indep",
+  rotation = 0,
+  parameters = numeric(0),
+  var_types = c("c", "c")
+) {
   assert_that(is.string(family), is.number(rotation), is.numeric(parameters))
   assert_that(correct_var_types(var_types))
   if (family %in% setdiff(family_set_nonparametric, "indep")) {
@@ -283,7 +301,9 @@ bicop_dist <- function(family = "indep", rotation = 0, parameters = numeric(0),
       stopifnot(dim(parameters) == c(30, 30))
       margin_integrals <- c(rowMeans(parameters), colMeans(parameters))
       if (any(abs(margin_integrals - 1) > 0.5))
-        warning("margins implied by 'parameters' deviate strongly from the standard uniform distribution.")
+        warning(
+          "margins implied by 'parameters' deviate strongly from the standard uniform distribution."
+        )
     } else {
       parameters <- matrix(1, 30, 30)
     }

R/bicop_methods.R

@@ -125,8 +125,15 @@ rbicop <- function(n, family, rotation, parameters, qrng = FALSE) {
 #'    variable, `cond_var = 2` on the second.
 #' @param inverse whether to compute the h-function or its inverse.
 #' @export
-hbicop <- function(u, cond_var, family, rotation, parameters, inverse = FALSE,
-                   var_types = c("c", "c")) {
+hbicop <- function(
+  u,
+  cond_var,
+  family,
+  rotation,
+  parameters,
+  inverse = FALSE,
+  var_types = c("c", "c")
+) {
   assert_that(in_set(cond_var, 1:2), is.flag(inverse))
   bicop <- args2bicop(family, rotation, parameters, var_types)
   u <- if_vec_to_matrix(u)
@@ -267,16 +274,20 @@ print.bicop_dist <- function(x, ...) {
   if (x$family %in% setdiff(family_set_nonparametric, "indep")) {
     x$parameters <- paste0(round(x$npars, 2), sep = " d.f.")
   }
-  cat("Bivariate copula ('bicop_dist'): ",
-    "family = ", x$family,
-    ", rotation = ", x$rotation,
-    ", parameters = ", ifelse(length(x$parameters) > 1,
-      paste(round(x$parameters, 2),
-        collapse = ", "
-      ),
+  cat(
+    "Bivariate copula ('bicop_dist'): ",
+    "family = ",
+    x$family,
+    ", rotation = ",
+    x$rotation,
+    ", parameters = ",
+    ifelse(
+      length(x$parameters) > 1,
+      paste(round(x$parameters, 2), collapse = ", "),
       x$parameters
     ),
-    ", var_types = ", paste(x$var_types, collapse = ","),
+    ", var_types = ",
+    paste(x$var_types, collapse = ","),
     sep = ""
   )
   cat("\n")
@@ -296,11 +307,16 @@ print.bicop <- function(x, ...) {
   } else {
     pars_formatted <- paste(round(x$parameters, 2), collapse = ", ")
   }
-  cat("Bivariate copula fit ('bicop'): ",
-    "family = ", x$family,
-    ", rotation = ", x$rotation,
-    ", parameters = ", pars_formatted,
-    ", var_types = ", paste(x$var_types, collapse = ","),
+  cat(
+    "Bivariate copula fit ('bicop'): ",
+    "family = ",
+    x$family,
+    ", rotation = ",
+    x$rotation,
+    ", parameters = ",
+    pars_formatted,
+    ", var_types = ",
+    paste(x$var_types, collapse = ","),
     "\n",
     sep = ""
   )

R/bicop_plot.R

@@ -64,7 +64,8 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
   }
 
   if (is.null(list(...)$xlim) & is.null(list(...)$ylim)) {
-    xylim <- switch(margins,
+    xylim <- switch(
+      margins,
       "unif" = c(1e-2, 1 - 1e-2),
       "norm" = c(-3, 3),
       "exp" = c(0, 6),
@@ -76,7 +77,8 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
 
   ## prepare for plotting with selected margins
   if (margins == "unif") {
-    points <- switch(type,
+    points <- switch(
+      type,
       "contour" = seq(1e-5, 1 - 1e-5, length.out = size),
       "surface" = 1:size / (size + 1)
     )
@@ -140,12 +142,19 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
     dens <- pmin(dens, 6)
   }
 
-  jet.colors <- colorRampPalette(c(
-    "#00007F", "blue", "#007FFF", "cyan",
-    "#7FFF7F", "yellow", "#FF7F00", "red",
-    "#7F0000"
-  ),
-  bias = 2
+  jet.colors <- colorRampPalette(
+    c(
+      "#00007F",
+      "blue",
+      "#007FFF",
+      "cyan",
+      "#7FFF7F",
+      "yellow",
+      "#FF7F00",
+      "red",
+      "#7F0000"
+    ),
+    bias = 2
   )
   ## actual plotting
   if (type == "contour") {
@@ -157,13 +166,15 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
       levels = levels,
       xlim = xlim,
       ylim = ylim,
-      xlab = switch(margins,
+      xlab = switch(
+        margins,
         "unif" = expression(u[1]),
         "norm" = expression(z[1]),
         "exp" = expression(e[1]),
         "flexp" = expression(e[1])
       ),
-      ylab = switch(margins,
+      ylab = switch(
+        margins,
         "unif" = expression(u[2]),
         "norm" = expression(z[2]),
         "exp" = expression(e[2]),
@@ -179,7 +190,8 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
       x = c ~ u * v,
       data = lst,
       scales = list(arrows = FALSE),
-      drape = TRUE, colorkey = FALSE,
+      drape = TRUE,
+      colorkey = FALSE,
       screen = list(z = 25, x = -55),
       shade = FALSE,
       aspect = c(1, 1),
@@ -188,20 +200,23 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
       par.settings = list(axis.line = list(col = "transparent")),
       col = NA,
       col.regions = jet.colors(100),
-      xlab = switch(margins,
+      xlab = switch(
+        margins,
         "unif" = expression(u[1]),
         "norm" = expression(z[1]),
         "exp" = expression(e[1]),
         "flexp" = expression(e[1])
       ),
-      ylab = switch(margins,
+      ylab = switch(
+        margins,
         "unif" = expression(u[2]),
         "norm" = expression(z[2]),
         "exp" = expression(e[2]),
         "flexp" = expression(e[2])
       ),
       zlab = "",
-      zlim = switch(margins,
+      zlim = switch(
+        margins,
         "unif" = c(0, max(3, 1.1 * max(lst$c))),
         "norm" = c(0, max(0.4, 1.1 * max(lst$c))),
         "exp" = c(0, max(1, 1.1 * max(lst$c))),
@@ -211,7 +226,8 @@ plot.bicop_dist <- function(x, type = "surface", margins, size, ...) {
   }
 
   pars <- modifyList(pars, list(...))
-  plot_obj <- switch(type,
+  plot_obj <- switch(
+    type,
     contour = do.call(contour, pars),
     surface = print(do.call(wireframe, pars))
   )