Export significant test results (#37)

Make the list that stores hypothesis-test results into an S3 object sig_test_summary.
Add function to export top-tables and significant gene-lists from a sig_test_summary object to consistently-named files.

* reordered package imports alphabetically

* make class for sig test summaries

* save sig-test results using export() generic

* drop `test_type` from filenames

`test_type` (eg, single_contrast, nested_lrt) should be part
of the user-specified `output_dir` for simplicity.

* fix roxygen for export.default

DESCRIPTION

@@ -17,29 +17,30 @@ Suggests:
     covr,
     mockery
 Imports:
-    magrittr,
-    edgeR,
+    biomaRt,
     Biobase,
-    methods,
+    cqn,
+    edgeR,
     dplyr,
-    readr (>= 1.0.0),
+    limma,
+    magrittr,
+    methods,
     purrr,
-    biomaRt,
-    tibble,
-    stringr,
-    tidyr,
+    readr (>= 1.0.0),
     rlang,
-    cqn,
-    limma
+    stringr,
+    tibble,
+    tidyr
 Collate:
     'coverage.R'
     'data_validity.R'
     'dgelist_functions.R'
     'feature_annotation.R'
     'feature_counts.R'
     'feature_counts_io.R'
-    'hypothesis_tests.R'
     'utilities.R'
     'logfile_helpers.R'
     'logfile_cutadapt_io.R'
     'logfile_hisat2_io.R'
+    'sig_test_summary-class.R'
+    'sig_test_summary-methods.R'

NAMESPACE

@@ -3,6 +3,7 @@
 export(append_feature_annotations)
 export(cbind_feature_counts)
 export(cqn_dgelist)
+export(export)
 export(feature_counts_to_dgelist)
 export(filter_by_read_count)
 export(get_gc_percent)
@@ -41,6 +42,7 @@ importFrom(purrr,map_lgl)
 importFrom(readr,parse_number)
 importFrom(readr,read_file)
 importFrom(readr,read_tsv)
+importFrom(readr,write_tsv)
 importFrom(rlang,.data)
 importFrom(stats,setNames)
 importFrom(stringr,str_replace)

R/hypothesis_tests.R → R/sig_test_summary-class.R

@@ -5,7 +5,6 @@
 #' @noRd
 #'
 get_lrt <- function(fit, testable_features, contrast) {
-  message("Untested: `reeq::get_lrt`")
 
   stopifnot(is(fit, "DGEGLM"))
   stopifnot(all(testable_features %in% rownames(fit)))
@@ -19,7 +18,6 @@ get_lrt <- function(fit, testable_features, contrast) {
 #' @noRd
 #'
 get_sig_features <- function(lrt, p_threshold) {
-  message("Untested: `reeq::get_sig_features`")
 
   # TODO: expand to use limma::decideTests.MArrayLM
   stopifnot(is(lrt, "DGEExact") || is(lrt, "DGELRT"))
@@ -45,19 +43,20 @@ get_sig_test_summary <- function(fit,
                                  testable_features = row.names(fit),
                                  contrast,
                                  p_threshold = 0.05) {
-  message("Untested: `reeq::get_sig_test_summary`")
-
   lrt <- get_lrt(fit, testable_features, contrast)
   sig_features <- get_sig_features(lrt, p_threshold)
   top_table <- edgeR::topTags(lrt, n = Inf)
 
-  list(
-    features = testable_features,
-    sig_features = sig_features,
-    lrt = lrt,
-    top_table = top_table,
-    num_features = length(testable_features),
-    num_sig_features = length(sig_features),
-    p_threshold = p_threshold
+  structure(
+    list(
+      features = testable_features,
+      sig_features = sig_features,
+      lrt = lrt,
+      top_table = top_table,
+      num_features = length(testable_features),
+      num_sig_features = length(sig_features),
+      p_threshold = p_threshold
+    ),
+    class = "sig_test_summary"
   )
 }

R/sig_test_summary-methods.R

@@ -0,0 +1,67 @@
+#' Export a data structure to a (set of) file(s)
+#'
+#' @param        x             Some datastructure. Typically a
+#'   `sig_test_summary` as created by `reeq::get_sig_test_summary`.
+#' @param        ...           Further arguments. Currently unused.
+#'
+#' @export
+
+export <- function(x, ...) {
+  UseMethod("export")
+}
+
+#' Default function for exporting a data-structure.
+#'
+#' @inheritParams   export
+#'
+#' @noRd
+#'
+
+export.default <- function(x, ...) {
+  stop("Not implemented")
+}
+
+#' Export a `sig_test_summary` object to files
+#'
+#' Write the (unfiltered) top-table and significant gene lists from a
+#' `sig_test_summary` object to some files. The file format is
+#' `<output_dir>/<test_name>.top_tags.tsv` or
+#' `<output_dir>/<test_name>.sig_features.<p_threshold>.tsv`, where the p-value
+#' threshold is obtained from inside the `sig_test_summary` object and is
+#' formatted to 4 decimal places.
+#'
+#' @inheritParams   export
+#'
+#' @param        output_dir    Which directory should the results be output to
+#' @param        test_name     A name for the hypothesis test that is
+#'   summarised in the dataset `x`.
+#'
+#' @importFrom   readr         write_tsv
+#'
+
+export.sig_test_summary <- function(x,
+                                    output_dir,
+                                    test_name,
+                                    ...) {
+  p_threshold <- sprintf("p%.4f", x$p_threshold)
+
+  # obtain the top-table for all features
+  tt_file <- paste(
+    test_name, "top_tags", "tsv",
+    sep = "."
+  )
+  tt_path <- file.path(output_dir, tt_file)
+  tt_data <- as.data.frame(x$top_table)
+
+  # obtain the list of all significant features
+  gene_file <- paste(
+    test_name, "sig_features", p_threshold, "tsv",
+    sep = "."
+  )
+  gene_path <- file.path(output_dir, gene_file)
+  gene_data <- x$sig_features
+
+  # write out files containing the results
+  readr::write_tsv(tt_data, tt_path)
+  write(gene_data, gene_path, ncolumns = 1)
+}

tests/testthat/test-sig_test_summary-methods.R

@@ -0,0 +1,61 @@
+###############################################################################
+
+context("Tests for methods on `sig_test_summary` objects")
+
+###############################################################################
+
+test_that("export - invalid input", {
+  expect_error(
+    object = export("NOT A sig_test_summary"),
+    info = "`reeq::export` is only defined for `sig_test_summary` objects"
+  )
+})
+
+test_that("export - creates the expected files", {
+  # TODO: explicit constructor for `sig_test_summary` objects
+
+  # `export.sig_test_summary` requires that `p_threshold`, `sig_features` and
+  #   `top_table` are present in the `sig_test_summary` that is passed in.
+
+  input <- structure(
+    list(
+      # 4 significant digits of the p-threshold are present in the output
+      #   "significant features" file
+      p_threshold = 0.001111,
+      top_table = data.frame(x = 1:3, y = 4:6, row.names = LETTERS[1:3]),
+      sig_features = "A"
+    ),
+    class = "sig_test_summary"
+  )
+
+  out_dir <- file.path(tempdir(), "export1")
+  dir.create(out_dir)
+  top_table_file <- file.path(
+    out_dir, "some_contrast.top_tags.tsv"
+  )
+  sig_features_file <- file.path(
+    out_dir, "some_contrast.sig_features.p0.0011.tsv"
+  )
+
+  expect_silent(
+    export(
+      input,
+      output_dir = out_dir,
+      test_name = "some_contrast"
+    )
+  )
+
+  expect_true(
+    file.exists(top_table_file),
+    info = "`export` should make a top-tags file"
+  )
+  expect_true(
+    file.exists(sig_features_file),
+    info = "`export` should make a significant-features file"
+  )
+  expect_equal(
+    object = scan(sig_features_file, what = character()),
+    expected = input$sig_features,
+    info = "contents of significant-features file should match the input"
+  )
+})