Merge pull request #77 from rgcca-factory/refactor/use_classes_for_bl…

…ocks

Refactor/use classes for blocks

NAMESPACE

@@ -1,5 +1,18 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(block_init,block)
+S3method(block_init,dual_block)
+S3method(block_init,dual_regularized_block)
+S3method(block_init,primal_regularized_block)
+S3method(block_postprocess,block)
+S3method(block_postprocess,sparse_block)
+S3method(block_project,block)
+S3method(block_project,dual_block)
+S3method(block_project,dual_regularized_block)
+S3method(block_project,primal_regularized_block)
+S3method(block_project,sparse_block)
+S3method(block_update,block)
+S3method(block_update,dual_block)
 S3method(plot,rgcca)
 S3method(plot,rgcca_bootstrap)
 S3method(plot,rgcca_cv)

R/block.R

@@ -0,0 +1,66 @@
+### Create classes
+new_block <- function(x, j, na.rm = TRUE, bias = TRUE,
+                      ..., class = character()) {
+  n <- NROW(x)
+  p <- NCOL(x)
+  N <- ifelse(bias, n, n - 1)
+
+  x <- list(
+    x = x,
+    j = j,
+    n = n,
+    p = p,
+    N = N,
+    na.rm = na.rm,
+    a = NULL,
+    Y = NULL,
+    ...
+  )
+  class(x) <- c(class, "block")
+  return(x)
+}
+
+new_dual_block <- function(x, j, na.rm = TRUE, ..., class = character()) {
+  K <- pm(x, t(x), na.rm = na.rm)
+  new_block(
+    x, j, na.rm, alpha = NULL, K = K, ..., class = c(class, "dual_block")
+  )
+}
+
+new_primal_regularized_block <- function(x, j, tau, ...) {
+  new_block(x, j, tau = tau, M = NULL, ..., class = "primal_regularized_block")
+}
+
+new_dual_regularized_block <- function(x, j, tau, ...) {
+  new_dual_block(
+    x, j, tau = tau, M = NULL, ..., class = "dual_regularized_block"
+  )
+}
+
+new_sparse_block <- function(x, j, sparsity, tol = 1e-08, ...) {
+  const <- sqrt(NCOL(x)) * sparsity
+  new_block(
+    x, j, sparsity = sparsity, const = const,
+    tol = tol, ..., class = "sparse_block"
+  )
+}
+
+### Utility method to choose the adequate class
+create_block <- function(x, j, bias, na.rm, tau, sparsity, tol) {
+  if (sparsity < 1) {
+    res <- new_sparse_block(x, j, sparsity, tol, bias = bias, na.rm = na.rm)
+  } else if (NROW(x) > NCOL(x)) {
+    if (tau < 1) {
+      res <- new_primal_regularized_block(x, j, tau, bias = bias, na.rm = na.rm)
+    } else {
+      res <- new_block(x, j, bias = bias, na.rm = na.rm)
+    }
+  } else {
+    if (tau < 1) {
+      res <- new_dual_regularized_block(x, j, tau, bias = bias, na.rm = na.rm)
+    } else {
+      res <- new_dual_block(x, j, bias = bias, na.rm = na.rm)
+    }
+  }
+  return(res)
+}

R/block_init.R

@@ -0,0 +1,41 @@
+#' @importFrom MASS ginv
+
+block_init <- function(x, init = "svd") {
+  UseMethod("block_init")
+}
+
+#' @export
+block_init.block <- function(x, init = "svd") {
+  if (init == "svd") {
+    x$a <- initsvd(x$x, dual = FALSE)
+  } else {
+    x$a <- rnorm(x$p)
+  }
+
+  return(block_project(x))
+}
+
+#' @export
+block_init.dual_block <- function(x, init = "svd") {
+  if (init == "svd") {
+    x$alpha <- initsvd(x$x, dual = TRUE)
+  } else {
+    x$alpha <- rnorm(x$n)
+  }
+
+  return(block_project(x))
+}
+
+#' @export
+block_init.primal_regularized_block <- function(x, init = "svd") {
+  x$M <- ginv(
+    x$tau * diag(x$p) + (1 - x$tau) * pm(t(x$x), x$x, na.rm = x$na.rm) / x$N
+  )
+  NextMethod()
+}
+
+#' @export
+block_init.dual_regularized_block <- function(x, init = "svd") {
+  x$M <- ginv(x$tau * diag(x$n) + (1 - x$tau) * x$K / x$N)
+  NextMethod()
+}

R/block_postprocess.R

@@ -0,0 +1,33 @@
+block_postprocess <- function(x, ctrl) {
+  UseMethod("block_postprocess")
+}
+
+#' @export
+block_postprocess.block <- function(x, ctrl) {
+  if (ctrl && (x$a[1] < 0)) {
+    x$a <- -x$a
+    x$Y <- -x$Y
+  }
+  return(x)
+}
+
+#' @export
+block_postprocess.sparse_block <- function(x, ctrl) {
+  l2_sat <- norm(x$a, "2")
+  if (abs(l2_sat - 1) > x$tol) {
+    if (l2_sat < .Machine$double.eps) {
+      warning(
+        "The l2 norm of the block weight vector #",
+        x$j, " is too small :", l2_sat
+      )
+    } else {
+      warning(
+        "The l2 constraint is not saturated for block #", x$j,
+        ". The intersection of the l1 and l2 spheres is empty for ",
+        "a sparsity parameter equal to ", x$sparsity,
+        ". Try to increase the value of the sparsity parameter."
+      )
+    }
+  }
+  NextMethod()
+}

R/block_project.R

@@ -0,0 +1,56 @@
+block_project <- function(x) {
+  UseMethod("block_project")
+}
+
+#' @export
+block_project.block <- function(x) {
+  if (any(x$a != 0)) {
+    x$a <- x$a / norm(x$a, type = "2")
+  }
+
+  x$Y <- pm(x$x, x$a, na.rm = x$na.rm)
+  return(x)
+}
+
+#' @export
+block_project.dual_block <- function(x) {
+  if (any(x$alpha != 0)) {
+    x$alpha <- x$alpha / drop(sqrt(t(x$alpha) %*% x$K %*% x$alpha))
+  }
+  x$a <- pm(t(x$x), x$alpha, na.rm = x$na.rm)
+
+  x$Y <- pm(x$x, x$a, na.rm = x$na.rm)
+  return(x)
+}
+
+#' @export
+block_project.primal_regularized_block <- function(x) {
+  if (any(x$a != 0)) {
+    x$a <- x$M %*% x$a / drop(sqrt(t(x$a) %*% x$M %*% x$a))
+  }
+
+  x$Y <- pm(x$x, x$a, na.rm = x$na.rm)
+  return(x)
+}
+
+#' @export
+block_project.dual_regularized_block <- function(x) {
+  if (any(x$alpha != 0)) {
+    x$alpha <- x$M %*% x$alpha / drop(sqrt(
+      t(x$alpha) %*% x$M %*% x$K %*% x$alpha
+    ))
+  }
+  x$a <- pm(t(x$x), x$alpha, na.rm = x$na.rm)
+
+  x$Y <- pm(x$x, x$a, na.rm = x$na.rm)
+  return(x)
+}
+
+#' @export
+block_project.sparse_block <- function(x) {
+  if (any(x$a != 0)) {
+    x$a <- soft_threshold(x$a, x$const)
+  }
+  x$Y <- pm(x$x, x$a, na.rm = x$na.rm)
+  return(x)
+}

R/block_update.R

@@ -0,0 +1,15 @@
+block_update <- function(x, grad) {
+  UseMethod("block_update")
+}
+
+#' @export
+block_update.block <- function(x, grad) {
+  x$a <- pm(t(x$x), grad, na.rm = x$na.rm)
+  return(block_project(x))
+}
+
+#' @export
+block_update.dual_block <- function(x, grad) {
+  x$alpha <- grad
+  return(block_project(x))
+}

R/rgcca.R

@@ -486,7 +486,7 @@ rgcca <- function(blocks, connection = NULL, tau = 1, ncomp = 1,
   gcca_args[["blocks"]] <- blocks
   gcca_args[["disjunction"]] <- opt$disjunction
   gcca_args[[opt$param]] <- rgcca_args[[opt$param]]
-  func_out <- do.call(opt$gcca, gcca_args)
+  func_out <- do.call(rgcca_outer_loop, gcca_args)
 
   ### Format the output
   func_out <- format_output(func_out, rgcca_args, opt, blocks)

R/rgcca_cv.r

@@ -159,6 +159,7 @@
 #'  plot(cv_out, display_order = FALSE)
 #' }
 rgcca_cv <- function(blocks,
+                     connection = NULL,
                      method = "rgcca",
                      response = NULL,
                      par_type = "tau",
@@ -228,7 +229,6 @@ rgcca_cv <- function(blocks,
     par_type <- "sparsity"
   } else if (par_type == "sparsity") {
     rgcca_args$method <- "sgcca"
-    opt$gcca <- sgcca
     opt$param <- "sparsity"
   }