Adapt to wide format of results table for results on gene subsets

gagneurlab · Apr 23, 2024 · a082e8c · a082e8c
1 parent 59ddb1e
commit a082e8c
Show file tree

Hide file tree

Showing 2 changed files with 61 additions and 4 deletions.
diff --git a/R/AllGenerics.R b/R/AllGenerics.R
@@ -992,6 +992,8 @@ setMethod("results", "FraserDataSet", function(object,
                 aggregate=aggregate, collapse=collapse, geneColumn=geneColumn,
                 subsetName=NULL, additionalColumns=additionalColumns, 
                 BPPARAM=BPPARAM)
+        # to restore previous column order
+        prevColOrder <- colnames(as.data.table(res))
     }
 
     # add results for FDR_subsets if requested
@@ -1010,16 +1012,63 @@ setMethod("results", "FraserDataSet", function(object,
             res <- unlist(GRangesList(unlist(list(res, resls_subsets))))
         } else{
             res <- unlist(GRangesList(unlist(resls_subsets)))
+            # restore previous column order
+            prevColOrder <- colnames(as.data.table(resls_subsets[[1]]))
         }
+    }
+
+    # sort it if existing
+    if(length(res) > 0){
+
+        # dcast to have only one row per gene-sample combination
+        res <- as.data.table(res)
+        if(isTRUE(aggregate)){
+            dcastCols <- c("pValueGene", "padjustGene")
+        } else{
+            dcastCols <- "padjust"
+        }
+        res <- dcast(res, ... ~ FDR_set, value.var=dcastCols)
 
-        # sort it if existing
-        if(length(res) > 0){
-            res <- res[order(res$pValue)]
+        # rename column back to padjust for tw results after dcast
+        if(any(grepl("transcriptome-wide", colnames(res)))){
+            if(isTRUE(aggregate)){
+                setnames(res, "pValueGene_transcriptome-wide", "pValueGene")
+                setnames(res, "padjustGene_transcriptome-wide", "padjustGene")
+            } else{
+                setnames(res, "transcriptome-wide", "padjust")
+            }
+
+        } else{
+            res[, padjust := NA]
             if(isTRUE(aggregate)){
-                res <- res[!is.na(res$pValueGene)]
+                res[, pValueGene := NA]
+                res[, padjustGene := NA]
             }
         }
+
+        # restore previous column order
+        setcolorder(res, prevColOrder[prevColOrder != "FDR_set"])
+
+        # rename padjust columns for other gene sets to padjust_setname
+        for(setName in FDRsets){
+            setnames(res, setName, paste0("padjust_", setName), 
+                        skip_absent=TRUE)
+        }
+
+        # revert back to GRanges
+        res <- makeGRangesFromDataFrame(res, keep.extra.columns=TRUE)
+
+        # order rows based on transcriptome-wide FDR
+        res <- res[order(res$pValue)]
+        if(isTRUE(aggregate) & isTRUE(returnTranscriptomewideResults)){
+            res <- res[!is.na(res$pValueGene)]
+        }
+
+
+    } else{
+        res$FDR_set <- NULL
     }
+
     return(res)
 })
 

diff --git a/tests/testthat/test_plotSampleResults.R b/tests/testthat/test_plotSampleResults.R
@@ -35,6 +35,14 @@ test_that("Results function", {
     expect_equal(res_gene_signif[type == "psi3", .N], 2)
     expect_equal(res_gene_signif[type == "theta", .N], 3)
 
+    # results on subset of genes during FDR
+    geneList <- list('sample1'=c("TIMMDC1"), 'sample2'=c("MCOLN1"))
+    fds <- calculatePadjValues(fds, type="jaccard", 
+                     subsets=list("exampleSubset"=geneList))
+    expect_equal(length(results(fds, all=TRUE, psiType="jaccard")), 
+                prod(dim(fds)))
+    expect_error(results(fds, all=TRUE, psiType="psi5"))
+
 })
 
 test_that("Main plotting function", {