Refactor analysis functions to utilize cheapr for efficiency (#134)

Author: araikes

DESCRIPTION

@@ -37,7 +37,8 @@ Imports:
     pbmcapply,
     rhdf5,
     tibble,
-    tidyr
+    tidyr,
+    cheapr (>= 1.5.0)
 RoxygenNote: 7.3.3
 Roxygen: list(markdown = TRUE)
 Suggests:

R/ModelArray_Constructor.R

@@ -503,7 +503,7 @@ analyseOneElement.lm <- function(i_element,
     if (flag_initiate) {
       return(list(column_names = NaN, list.terms = NaN))
     } else {
-      onerow <- c(i_element - 1, rep(NaN, (num.stat.output - 1)))
+      onerow <- c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L))
       return(onerow)
     }
   }
@@ -545,7 +545,7 @@ analyseOneElement.lm <- function(i_element,
     if (flag_initiate) {
       return(list(column_names = NaN, list.terms = NaN))
     } else {
-      onerow <- c(i_element - 1, rep(NaN, (num.stat.output - 1)))
+      onerow <- c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L))
       return(onerow)
     }
   }
@@ -693,7 +693,7 @@ analyseOneElement.gam <- function(i_element,
         sp.criterion.attr.name = NaN
       ))
     } else {
-      onerow <- c(i_element - 1, rep(NaN, (num.stat.output - 1)))
+      onerow <- c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L))
       return(onerow)
     }
   }
@@ -739,7 +739,7 @@ analyseOneElement.gam <- function(i_element,
         sp.criterion.attr.name = NaN
       ))
     } else {
-      onerow <- c(i_element - 1, rep(NaN, (num.stat.output - 1)))
+      onerow <- c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L))
       return(onerow)
     }
   }
@@ -976,7 +976,7 @@ analyseOneElement.wrap <- function(i_element,
     if (flag_initiate) {
       return(list(column_names = NaN))
     } else {
-      onerow <- c(i_element - 1, rep(NaN, (num.stat.output - 1)))
+      onerow <- c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L))
       return(onerow)
     }
   }
@@ -1019,7 +1019,7 @@ analyseOneElement.wrap <- function(i_element,
     if (flag_initiate) {
       return(list(column_names = NaN))
     } else {
-      onerow <- c(i_element - 1, rep(NaN, (num.stat.output - 1)))
+      onerow <- c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L))
       return(onerow)
     }
   }
@@ -1033,7 +1033,7 @@ analyseOneElement.wrap <- function(i_element,
       if (on_error == "skip" || on_error == "debug") {
         warning(paste0("analyseOneElement.wrap at element ", i_element, ": ", msg))
         if (flag_initiate) return(list(column_names = NaN))
-        return(c(i_element - 1, rep(NaN, (num.stat.output - 1))))
+        return(c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L)))
       }
       stop(msg)
     }
@@ -1057,7 +1057,7 @@ analyseOneElement.wrap <- function(i_element,
     if (on_error == "skip" || on_error == "debug") {
       warning(paste0("analyseOneElement.wrap at element ", i_element, ": ", msg))
       if (flag_initiate) return(list(column_names = NaN))
-      return(c(i_element - 1, rep(NaN, (num.stat.output - 1))))
+      return(c(i_element - 1, cheapr::rep_len_(NaN, num.stat.output - 1L)))
     }
     stop(msg)
   }

R/analyse-helpers.R

@@ -426,7 +426,7 @@
   response_vals <- scalars(ctx$modelarray)[[ctx$scalar]][i_element, ]
 
   # Start the validity mask with the response scalar
-  masks <- list(is.finite(response_vals))
+  valid_mask <- is.finite(response_vals)
 
   # Read and reorder additional attached scalars
   scalar_values <- list()
@@ -443,22 +443,25 @@
     }
 
     scalar_values[[sname]] <- s_vals
-    masks[[length(masks) + 1L]] <- is.finite(s_vals)
+    valid_mask <- valid_mask & is.finite(s_vals)
   }
 
   # Intersection mask across all scalars
-  valid_mask <- Reduce("&", masks)
-  num_valid <- sum(valid_mask)
+  which_valid <- cheapr::which_(valid_mask)
+  num_valid <- length(which_valid)
 
   if (!(num_valid > num.subj.lthr)) {
     return(list(dat = NULL, sufficient = FALSE, num_valid = num_valid))
   }
 
   # Build filtered data.frame
-  dat <- ctx$phenotypes[valid_mask, , drop = FALSE]
-  for (sname in ctx$attached_scalars) {
-    dat[[sname]] <- scalar_values[[sname]][valid_mask]
-  }
+  dat <- cheapr::sset(ctx$phenotypes, which_valid)
+  scalar_col_list <- lapply(
+    stats::setNames(ctx$attached_scalars, ctx$attached_scalars),
+    function(sname) scalar_values[[sname]][which_valid]
+  )
+  scalar_cols <- cheapr::fast_df(.args = scalar_col_list)
+  dat <- cheapr::col_c(dat, scalar_cols)
 
   list(dat = dat, sufficient = TRUE, num_valid = num_valid)
 }

R/analyse.R

@@ -149,7 +149,9 @@ ModelArray.lm <- function(formula, data, phenotypes, scalar = NULL, element.subs
                           ...) {
   # Validation ----
   .validate_modelarray_input(data)
+
   scalar <- .resolve_formula_scalar(formula, data, scalar)
+
   element.subset <- .validate_element_subset(element.subset, data, scalar)
   phenotypes <- .align_phenotypes(data, phenotypes, scalar)
 
@@ -300,9 +302,10 @@ ModelArray.lm <- function(formula, data, phenotypes, scalar = NULL, element.subs
     return(invisible(NULL))
   }
 
-  df_out <- do.call(rbind, fits_all)
-  df_out <- as.data.frame(df_out)
-  colnames(df_out) <- column_names
+  result_mat <- do.call(rbind, fits_all)
+  col_list <- lapply(seq_len(ncol(result_mat)), function(j) result_mat[, j])
+  names(col_list) <- column_names
+  df_out <- cheapr::fast_df(.args = col_list)
 
   df_out <- .correct_pvalues(df_out, list.terms, correct.p.value.terms, var.terms)
   df_out <- .correct_pvalues(df_out, "model", correct.p.value.model, var.model)
@@ -474,7 +477,9 @@ ModelArray.gam <- function(formula, data, phenotypes, scalar = NULL, element.sub
                            ...) {
   # Validation ----
   .validate_modelarray_input(data)
+
   scalar <- .resolve_formula_scalar(formula, data, scalar)
+
   element.subset <- .validate_element_subset(element.subset, data, scalar)
   phenotypes <- .align_phenotypes(data, phenotypes, scalar)
 
@@ -695,9 +700,10 @@ ModelArray.gam <- function(formula, data, phenotypes, scalar = NULL, element.sub
     return(invisible(NULL))
   }
 
-  df_out <- do.call(rbind, fits_all)
-  df_out <- as.data.frame(df_out)
-  colnames(df_out) <- column_names
+  result_mat <- do.call(rbind, fits_all)
+  col_list <- lapply(seq_len(ncol(result_mat)), function(j) result_mat[, j])
+  names(col_list) <- column_names
+  df_out <- cheapr::fast_df(.args = col_list)
 
   # P-value corrections ----
   df_out <- .correct_pvalues(df_out, list.smoothTerms, correct.p.value.smoothTerms, var.smoothTerms)
@@ -756,9 +762,13 @@ ModelArray.gam <- function(formula, data, phenotypes, scalar = NULL, element.sub
         ...
       )
 
-      reduced.model.df_out <- do.call(rbind, reduced.model.fits)
-      reduced.model.df_out <- as.data.frame(reduced.model.df_out)
-      colnames(reduced.model.df_out) <- reduced.model.column_names
+      result_mat_reduced <- do.call(rbind, reduced.model.fits)
+      col_list_reduced <- lapply(
+        seq_len(ncol(result_mat_reduced)),
+        function(j) result_mat_reduced[, j]
+      )
+      names(col_list_reduced) <- reduced.model.column_names
+      reduced.model.df_out <- cheapr::fast_df(.args = col_list_reduced)
 
       ## Compute delta adj R-sq and partial R-sq ----
       delta_col <- paste0(changed.rsq.term.shortFormat, ".delta.adj.rsq")
@@ -1104,8 +1114,10 @@ ModelArray.wrap <- function(FUN, data, phenotypes, scalar, element.subset = NULL
     return(invisible(NULL))
   }
 
-  df_out <- do.call(rbind, fits_all)
-  df_out <- as.data.frame(df_out)
-  colnames(df_out) <- column_names
+  # Better proposed change:
+  result_mat <- do.call(rbind, fits_all)
+  col_list <- lapply(seq_len(ncol(result_mat)), function(j) result_mat[, j])
+  names(col_list) <- column_names
+  df_out <- cheapr::fast_df(.args = col_list)
   df_out
 }