From f24ffb082aa43525765ea37333b3500d0302ba3f Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 9 Aug 2021 18:37:03 -0700 Subject: [PATCH 01/28] Initial news --- NEWS.md | 40 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 40 insertions(+) create mode 100644 NEWS.md diff --git a/NEWS.md b/NEWS.md new file mode 100644 index 0000000..2d24ec1 --- /dev/null +++ b/NEWS.md @@ -0,0 +1,40 @@ +# GeomxTools 1.1.1 + +* + + + +# GeomxTools 0.99.4 - concomittant development branch version (includes changes beyond 1.0.0) + +## New features: +* +* Segment `setSegmentQCFlags()` and probe `setBioProbeQC()` QC function +* + + +## Revision: +* Change read + +# GeomxTools 1.1.0 + +* No changes from 1.0.0 + +# GeomxTools 1.0.0 + +* Initial release version, includes load with automatic aggregation of counts to target level. + +# GeomxTools 0.99.3 + +* Reverse change in R version requirement + +# GeomxTools 0.99.2 + +* Update R version requirement + +# GeomxTools 0.99.1 + +* Update R version requirement + +# GeomxTools 0.99.0 + +* Package template creation \ No newline at end of file From a812f7cc7d38a8ac2695cb6375d69f8b016599ba Mon Sep 17 00:00:00 2001 From: Tyler Hether Date: Wed, 11 Aug 2021 11:47:50 -0700 Subject: [PATCH 02/28] rbind to bind_rows --- R/readNanoStringGeoMxSet.R | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/R/readNanoStringGeoMxSet.R b/R/readNanoStringGeoMxSet.R index 2a2dfd1..2dae98f 100644 --- a/R/readNanoStringGeoMxSet.R +++ b/R/readNanoStringGeoMxSet.R @@ -180,12 +180,13 @@ function(dccFiles, # Create protocolData protocol <- - do.call(rbind, + do.call(dplyr::bind_rows, lapply(names(data), function(i) { cbind(data[[i]][["Header"]], data[[i]][["Scan_Attributes"]], data[[i]][["NGS_Processing_Attributes"]]) })) + protocol <- data.frame(protocol, pheno@data[, which(colnames(pheno@data) %in% protocolDataColNames)], check.names = FALSE) From 141ff42e73744b1bc09b6bdbadb9e556e202987f Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Fri, 13 Aug 2021 12:47:41 -0700 Subject: [PATCH 03/28] Bump version --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 118c7cc..6370ba5 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -2,7 +2,7 @@ Package: GeomxTools Title: NanoString GeoMx Tools Description: Tools for NanoString Technologies GeoMx Technology. Package provides functions for reading in DCC and PKC files based on an ExpressionSet derived object. Normalization and QC functions are also included. -Version: 1.1.1 +Version: 1.1.2 Encoding: UTF-8 Authors@R: c(person("Nicole", "Ortogero", email = "nortogero@nanostring.com", role = c("cre", "aut")), person("Zhi", "Yang", email = "zyang@nanostring.com", role = c("aut"))) From fa2dc7d587bc151b1708db0b107ef89d8966c70b Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Fri, 13 Aug 2021 15:10:51 -0700 Subject: [PATCH 04/28] Fix single panel bug in bioprobeqc --- R/NanoStringGeoMxSet-qc.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 7c4ffa9..4499216 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -329,7 +329,7 @@ setGrubbsFlags <- function(object, assayDataElement(modObject, elt="log10") <- logtBase(exprs(modObject), base=10) targetOutliers <- - esBy(modObject, GROUP="TargetName", FUN=function(y) { + esBy(modObject, GROUP="TargetName", simplify=FALSE, FUN=function(y) { currExprs <- assayDataElement(y, elt="log10") currResult <- apply(currExprs, 2, function(z) { if (max(z) < logtBase(minCount, base=10)) { From 721e9e80c208734766010573f1da3ff82df59f0d Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Fri, 13 Aug 2021 15:50:28 -0700 Subject: [PATCH 05/28] Suppress grubbs warnings and output less in vignette --- R/NanoStringGeoMxSet-qc.R | 5 ++++- .../Developer_Introduction_to_the_NanoStringGeoMxSet.Rmd | 2 +- 2 files changed, 5 insertions(+), 2 deletions(-) diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 4499216..3a09bf9 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -339,7 +339,10 @@ setGrubbsFlags <- function(object, return(NA) } grubbsTest <- - outliers::grubbs.test(z, two.sided=TRUE, type=10) + suppressWarnings( + outliers::grubbs.test(z, + two.sided=TRUE, + type=10)) if(grubbsTest$p.value < alphaCutoff & !is.null(names(grubbsTest$p.value))) { if (grepl("lowest", tolower(grubbsTest$alternative))) { diff --git a/vignettes/Developer_Introduction_to_the_NanoStringGeoMxSet.Rmd b/vignettes/Developer_Introduction_to_the_NanoStringGeoMxSet.Rmd index aa37e2b..e7bac36 100644 --- a/vignettes/Developer_Introduction_to_the_NanoStringGeoMxSet.Rmd +++ b/vignettes/Developer_Introduction_to_the_NanoStringGeoMxSet.Rmd @@ -220,7 +220,7 @@ demoData <- shiftCountsOne(demoData, useDALogic=TRUE) demoData <- setSegmentQCFlags(demoData) head(protocolData(demoData)[["QCFlags"]]) demoData <- setBioProbeQCFlags(demoData) -head(featureData(demoData)[["QCFlags"]]) +featureData(demoData)[["QCFlags"]][1:5, 1:4] ``` Probes and Samples that were flagged can be removed from analysis by subsetting. From fdd82ab36305d0fa7110f73f1dd42d4b2dd62163 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 16 Aug 2021 17:34:01 -0700 Subject: [PATCH 06/28] Handle objects with only single probes, speed up agg, add warnings for no agg or probe qc. --- R/NanoStringGeoMxSet-aggregate.R | 32 ++++++++++++++++++++++++++++++-- R/NanoStringGeoMxSet-class.R | 1 + R/NanoStringGeoMxSet-qc.R | 4 +++- 3 files changed, 34 insertions(+), 3 deletions(-) diff --git a/R/NanoStringGeoMxSet-aggregate.R b/R/NanoStringGeoMxSet-aggregate.R index 2f476ef..91d8de8 100644 --- a/R/NanoStringGeoMxSet-aggregate.R +++ b/R/NanoStringGeoMxSet-aggregate.R @@ -16,9 +16,36 @@ #' @export #' aggregateCounts <- function(object, FUN=ngeoMean) { - object <- summarizeNegatives(object) - targetCounts <- do.call(rbind, esBy(object, GROUP = "TargetName", + if(featureType(object) == "Target") { + stop("GeoMxSet object feature type is already target-level. ", + "No further aggregation can be performed.") + } + + # Skip targets with single probe + multiProbeTable <- with(object, table(TargetName)) > 1L + multiProbeTargs <- + names(multiProbeTable)[which(multiProbeTable, arr.ind=TRUE)] + if (length(multiProbeTargs) > 0) { + multiObject <- + object[fData(object)[["TargetName"]] %in% multiProbeTargs, ] + object <- summarizeNegatives(object) + } else { + warning("Object has no multiprobe targets. ", + "No aggregation was performed.") + featureNames(object) <- fData(object)[["TargetName"]] + featureType(object) <- "Target" + return(object) + } + + targetCounts <- do.call(rbind, esBy(multiObject, GROUP = "TargetName", FUN=function(x) {esApply(x, 2, FUN)}, simplify=FALSE)) + singleProbeObject <- subset(object, + subset=!TargetName %in% multiProbeTargs) + singleProbeCounts <- exprs(singleProbeObject) + rownames(singleProbeCounts) <- fData(singleProbeObject)[["TargetName"]] + targetCounts <- rbind(targetCounts, singleProbeCounts) + targetCounts <- targetCounts[unique(fData(object)[["TargetName"]]), ] + targetFeats <- featureData(object)@data targetFeats <- targetFeats[!duplicated(targetFeats[["TargetName"]]), ] @@ -29,6 +56,7 @@ aggregateCounts <- function(object, FUN=ngeoMean) { targetFeats <- AnnotatedDataFrame(targetFeats[rownames(targetCounts), ], dimLabels = c("featureNames", "featureColumns")) + targetObject <- NanoStringGeoMxSet( assayData = targetCounts, phenoData = phenoData(object), diff --git a/R/NanoStringGeoMxSet-class.R b/R/NanoStringGeoMxSet-class.R index 74e709f..842e838 100644 --- a/R/NanoStringGeoMxSet-class.R +++ b/R/NanoStringGeoMxSet-class.R @@ -21,6 +21,7 @@ function(object) { methods::callNextMethod(object) cat("feature: ") cat(featureType(object)) + cat("\n") }) # Constructors diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 3a09bf9..8cfaf04 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -320,7 +320,9 @@ setGrubbsFlags <- function(object, multiObject <- object[fData(object)[["TargetName"]] %in% multiProbeTargs, ] } else { - multiObject <- object + warning("Object has no targets with at least 3 probes.", + "No outlier testing can be performed.") + return(object) } grubbsResults <- From d37e9c293dcfc1a04b881a67e6e78b180d874926 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 17 Aug 2021 17:16:02 -0700 Subject: [PATCH 07/28] probe ratio optimization --- R/NanoStringGeoMxSet-qc.R | 53 ++++++++++++++++++++++++++------------- R/utils.R | 9 +++---- 2 files changed, 39 insertions(+), 23 deletions(-) diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 8cfaf04..b68a721 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -283,24 +283,41 @@ setBioProbeQCFlags <- function(object, setProbeRatioFlags <- function(object, cutoff=DEFAULTS[["minProbeRatio"]]) { - rawTargetCounts <- collapseCounts(object) - rawTargetCounts[["Mean"]] <- - apply(rawTargetCounts[, sampleNames(object)], - MARGIN=1, FUN=ngeoMean) - probeMeans <- apply(assayDataElement(object, elt="exprs"), + # Skip targets with single probes + multiProbeTable <- with(object, table(TargetName)) > 1L + multiProbeTargs <- + names(multiProbeTable)[which(multiProbeTable, arr.ind=TRUE)] + if (length(multiProbeTargs) > 0) { + multiObject <- + object[fData(object)[["TargetName"]] %in% multiProbeTargs, ] + } else { + warning("Object has no multi-probe targets. ", + "No ratio testing can be performed.") + return(object) + } + + probeCounts <- data.table(cbind(fData(multiObject)[, c("TargetName", + "Module")], + list(RTS_ID=featureNames(multiObject)), + exprs(multiObject))) + collapsedCounts <- as.data.frame(probeCounts[, lapply(.SD, ngeoMean), + .SDcols=!c("RTS_ID"), + by=c("TargetName", "Module")]) + targetMeans <- + apply(collapsedCounts[, + colnames(collapsedCounts) %in% sampleNames(multiObject)], + 1, ngeoMean) + names(targetMeans) <- collapsedCounts[["TargetName"]] + + probeMeans <- apply(assayDataElement(multiObject, elt="exprs"), MARGIN=1, FUN=ngeoMean) - probeRatios <- - lapply(names(probeMeans), function(x) { - currTarget <- fData(object)[x, "TargetName"] - currModule <- fData(object)[x, "Module"] - currTargetMean <- - rawTargetCounts[rawTargetCounts[["TargetName"]] == currTarget & - rawTargetCounts[["Module"]] == currModule, - "Mean"] - ratios <- probeMeans[[x]] / currTargetMean - }) - probeRatios <- data.frame("ProbeRatio"=unlist(probeRatios), - row.names=names(probeMeans)) + multiProbeRatios <- probeMeans / targetMeans[probeCounts$TargetName] + + probeRatios <- data.frame("ProbeRatio"=rep(1, dim(object)[1L]), + row.names=featureNames(object)) + probeRatios[names(multiProbeRatios), "ProbeRatio"] <- + multiProbeRatios + featureData(object)[["ProbeRatio"]] <- probeRatios probeRatioFlags <- probeRatios <= cutoff colnames(probeRatioFlags) <- "LowProbeRatio" @@ -320,7 +337,7 @@ setGrubbsFlags <- function(object, multiObject <- object[fData(object)[["TargetName"]] %in% multiProbeTargs, ] } else { - warning("Object has no targets with at least 3 probes.", + warning("Object has no targets with at least 3 probes. ", "No outlier testing can be performed.") return(object) } diff --git a/R/utils.R b/R/utils.R index b482d05..4fc0998 100644 --- a/R/utils.R +++ b/R/utils.R @@ -118,11 +118,10 @@ countsShiftedByOne <- function(object) { return(experimentData(object)@other$shiftedByOne) } -#NEO rewrite in data.tables for speed collapseCounts <- function(object) { - probeCounts <- cbind(fData(object)[, c("TargetName", "Module")], - assayDataElement(object, elt="exprs")) - collapsedCounts <- aggregate(formula=. ~ TargetName + Module, - data=probeCounts, FUN=ngeoMean) + probeCounts <- data.table(cbind(fData(object)[, c("TargetName", "Module")], + assayDataElement(object, elt="exprs"))) + collapsedCounts <- probeCounts[, lapply(.SD, ngeoMean), + by=c("TargetName", "Module")] return(collapsedCounts) } \ No newline at end of file From 0c001045d5a4fd83434e38f534889dc28d6824a4 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 17 Aug 2021 17:28:50 -0700 Subject: [PATCH 08/28] Remove failed probe ratio probes prior to Grubbs --- R/NanoStringGeoMxSet-qc.R | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index b68a721..8543a59 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -342,6 +342,14 @@ setGrubbsFlags <- function(object, return(object) } + if (!"ProbeRatio" %in% fvarLabels(object)) { + warning("Probe ratio QC has not yet been performed. ", + "Suggests running probe ratio QC then re-running Grubbs QC.") + } else { + multiObject <- + subset(multiObject, subset=QCFlags[["LowProbeRatio"]] == FALSE) + } + grubbsResults <- lapply(unique(fData(multiObject)[["Module"]]), FUN=function(x) { modObject <- subset(multiObject, subset=Module == x) From f75379f793565f7a4b6e56c462cdb20db9917e0b Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 18 Aug 2021 10:19:11 -0700 Subject: [PATCH 09/28] Fix subsetting in qc --- R/NanoStringGeoMxSet-qc.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 8543a59..ae141c3 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -347,7 +347,7 @@ setGrubbsFlags <- function(object, "Suggests running probe ratio QC then re-running Grubbs QC.") } else { multiObject <- - subset(multiObject, subset=QCFlags[["LowProbeRatio"]] == FALSE) + multiObject[!fData(multiObject)[["QCFlags"]][, "LowProbeRatio"], ] } grubbsResults <- From 3e95823d494761efe3a491bda485934a9a8c6b59 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 18 Aug 2021 14:18:00 -0700 Subject: [PATCH 10/28] Update aggregate unit test --- R/NanoStringGeoMxSet-aggregate.R | 7 +- tests/testthat/test_GeoMxSet_probeCollapse.R | 137 +++++++++---------- 2 files changed, 70 insertions(+), 74 deletions(-) diff --git a/R/NanoStringGeoMxSet-aggregate.R b/R/NanoStringGeoMxSet-aggregate.R index 91d8de8..7e4b3d1 100644 --- a/R/NanoStringGeoMxSet-aggregate.R +++ b/R/NanoStringGeoMxSet-aggregate.R @@ -28,7 +28,12 @@ aggregateCounts <- function(object, FUN=ngeoMean) { if (length(multiProbeTargs) > 0) { multiObject <- object[fData(object)[["TargetName"]] %in% multiProbeTargs, ] - object <- summarizeNegatives(object) + if ("Negative" %in% unique(fData(object)[["CodeClass"]])) { + object <- summarizeNegatives(object) + } else { + warning("Object has no negatives. ", + "No summary statistics for negatives will be calculated.") + } } else { warning("Object has no multiprobe targets. ", "No aggregation was performed.") diff --git a/tests/testthat/test_GeoMxSet_probeCollapse.R b/tests/testthat/test_GeoMxSet_probeCollapse.R index 1485576..eb93e38 100644 --- a/tests/testthat/test_GeoMxSet_probeCollapse.R +++ b/tests/testthat/test_GeoMxSet_probeCollapse.R @@ -1,81 +1,72 @@ -# sum(matches) == nrow(testData@assayData$exprs) should result in TRUE - library(GeomxTools) library(testthat) library(EnvStats) -datadir <- system.file("extdata", "DSP_NGS_Example_Data", - package="GeomxTools") -DCCFiles <- dir(datadir, pattern=".dcc$", full.names=TRUE) -PKCFiles <- unzip(zipfile = file.path(datadir, "/pkcs.zip")) -SampleAnnotationFile <- file.path(datadir, "annotations.xlsx") - -testData <- - suppressWarnings(readNanoStringGeoMxSet(dccFiles = DCCFiles, # QuickBase: readNanoStringGeomxSet, need to change it. - pkcFiles = PKCFiles, - phenoDataFile = SampleAnnotationFile, - phenoDataSheet = "CW005", - phenoDataDccColName = "Sample_ID", - protocolDataColNames = c("aoi", - "cell_line", - "roi_rep", - "pool_rep", - "slide_rep"), - experimentDataColNames = c("panel"))) - - -test_Data <- aggregateCounts(testData) - -DCCFiles <- DCCFiles[!basename(DCCFiles) %in% unique(sData(testData)$NTC_ID)] - -PKC <- readPKCFile(PKCFiles) -PKC$RTS_ID <- gsub("RNA", "RTS00", PKC$RTS_ID) - -numDCC <- 10 - -#random subset of 10 DCC files -DCCFiles <- DCCFiles[sample(1:length(DCCFiles), numDCC)] - -matches <- NULL -dcc <- 1 -while(dcc <= length(DCCFiles) & all(matches == TRUE)){ - DCC_file <- DCCFiles[dcc] - DCC <- suppressWarnings(readDccFile(DCC_file)) - DCC_file <- basename(DCC_file) - rownames(DCC$Code_Summary) <- gsub("RNA", "RTS00", rownames(DCC$Code_Summary)) - - for(i in unique(fData(test_Data)[["TargetName"]])){ - probes <- PKC$RTS_ID[which(PKC$Target == i)] - if(sum(is.na(DCC$Code_Summary[probes,"Count"])) > 0){ - #NAs are changed to 0 counts - DCC$Code_Summary[probes, "Count"][is.na(DCC$Code_Summary[probes, "Count"])] <- 0 - - #0 count doesn't meet minimum threshold so all counts are increased by threshold - #0.5 is default threshold from thresholdValues() - DCC$Code_Summary[probes, "Count"] <- DCC$Code_Summary[probes, "Count"]+0.5 - } - - matches <- c(matches, EnvStats::geoMean(DCC$Code_Summary[probes,"Count"]) == - test_Data@assayData$exprs[i,DCC_file]) - } - - dcc <- dcc + 1 -} - - - - -# req 1: test that the number of collapsed probe is correct:------ -testthat::test_that("test that the number of collapsed probe is correct", { - expect_true(sum(matches) == nrow(test_Data@assayData$exprs)*numDCC) +testData <- + readRDS(file= system.file("extdata", "DSP_NGS_Example_Data", + "demoData.rds", package = "GeomxTools")) +testData <- (shiftCountsOne(testData, elt="exprs", useDALogic=TRUE)) + +subData <- + subset(testData, subset=Module == "VnV_GeoMx_Hs_CTA_v1.2", + select=sampleNames(testData) %in% + sample(sampleNames(testData), + 10, replace=FALSE)) +subData <- + subset(subData, + subset=TargetName %in% + sample(unique(fData(subData)[["TargetName"]]), + 10, replace=FALSE)) + +subAggd <- aggregateCounts(subData) + +testthat::test_that("Feature type changed after aggregation", { + expect_true(featureType(subData) == "Probe") + expect_true(featureType(subAggd) == "Target") }) +testthat::test_that("Aggregated object has target dimension and annotations", { + expect_true(all(fData(subData)[["TargetName"]] %in% featureNames(subAggd))) + expect_true(length(unique(fData(subData)[["TargetName"]])) == + dim(subAggd)[[1L]]) + expect_true(dim(subData)[[2L]] == dim(subAggd)[[2L]]) + + targetLabels <- + fvarLabels(subData)[!fvarLabels(subData) %in% + c("RTS_ID", "QCFlags", "ProbeID")] + expect_true(all(targetLabels %in% fvarLabels(subAggd))) + expect_true(all(svarLabels(subData) %in% svarLabels(subAggd))) +}) +testthat::test_that("Target expression matrix contains aggregated counts", { + expect_true(all(colnames(exprs(subData)) == colnames(exprs(subAggd)))) + sameTargs <- intersect(fData(subData)[["TargetName"]], + rownames(exprs(subAggd))) + expect_true(length(sameTargs) == + length(unique(fData(subData)[["TargetName"]]))) + subList <- lapply(featureNames(subAggd), function(currTarg) { + aggdCounts <- exprs(subAggd)[currTarg, sampleNames(subData)] + targData <- exprs(subset(subData, subset=TargetName == currTarg)) + targMeans <- apply(targData, 2L, function(sampCount) { + ifelse(length(sampCount) > 1, ngeoMean(sampCount), sampCount)}) + return(all(aggdCounts == targMeans)) + }) + expect_true(all(unlist(subList))) +}) -# a = PKC$RTS_ID[which(PKC$Target == "ANGPT2")] -# EnvStats::geoMean(DCC$Code_Summary[a,"Count"]) -# EnvStats::geoMean(DCC$Code_Summary[a,"Count"],na.rm = T) -# test_Data@assayData$exprs["ANGPT2",DCC_file] -# DCC$Code_Summary[a,"Count"] -# EnvStats::geoMean(c(9,1,2,1)) -# EnvStats::geoMean(c(9,1,2,1,1)) \ No newline at end of file +testthat::test_that("Other aggregation functions work", { + subSum <- aggregateCounts(subData, FUN=sum) + expect_true(all(colnames(exprs(subData)) == colnames(exprs(subSum)))) + sameTargs <- intersect(fData(subData)[["TargetName"]], + rownames(exprs(subSum))) + expect_true(length(sameTargs) == + length(unique(fData(subData)[["TargetName"]]))) + subList <- lapply(featureNames(subSum), function(currTarg) { + aggdCounts <- exprs(subSum)[currTarg, sampleNames(subData)] + targData <- exprs(subset(subData, subset=TargetName == currTarg)) + targSums <- apply(targData, 2L, function(sampCount) { + ifelse(length(sampCount) > 1, sum(sampCount), sampCount)}) + return(all(aggdCounts == targSums)) + }) + expect_true(all(unlist(subList))) +}) \ No newline at end of file From db80ec85d794f23f50e9778cab7141c4bb73ab39 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 18 Aug 2021 16:48:50 -0700 Subject: [PATCH 11/28] Fix normalization test --- tests/testthat/test_normalization.R | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/tests/testthat/test_normalization.R b/tests/testthat/test_normalization.R index bd27540..c3bed0b 100644 --- a/tests/testthat/test_normalization.R +++ b/tests/testthat/test_normalization.R @@ -5,6 +5,7 @@ library(EnvStats) # load data demoData <- readRDS(file= system.file("extdata","DSP_NGS_Example_Data", "demoData.rds", package = "GeomxTools")) +demoData <- shiftCountsOne(demoData) demoData <- normalize(demoData , data_type="RNA", norm_method="quant", desiredQuantile = .9, toElt = "q_norm") ## test with user parameter demoData <- normalize(demoData, norm_method="quant") ## test defaults quantile = .75 @@ -25,7 +26,7 @@ test_that("quantile norm factors are present", { # Compute normalization factors from demoData # compute 90% quantile count for samples and divide by geomean thresh <- assayDataApply(demoData, 2, quantile, probs = .9) -expectedOutputData <- thresh/geoMean(thresh) +expectedOutputData <- thresh/ngeoMean(thresh) expectedOutputData <- unname(expectedOutputData) # Extract normalized data from object @@ -37,7 +38,7 @@ test_that("quantile norm factors are correct", { #### req 2b verify calculation of q75 norm factors # compute 75% quantile norm factors for samples and divide by geomean thresh <- assayDataApply(demoData, 2, quantile, probs = .75) -expectedOutputData <- thresh/geoMean(thresh) +expectedOutputData <- thresh/ngeoMean(thresh) expectedOutputData <- unname(expectedOutputData) # Extract normalized data from object @@ -51,7 +52,7 @@ test_that("quantile norm factors for 75% are correct", { # compute normalize count for samples and divide by geomean thresh <- assayDataApply(demoData, 2, quantile, probs = .75) norm_quant <- function(x){ - x <- x/(thresh/geoMean(thresh)) + x <- x/(thresh/ngeoMean(thresh)) } expectedOutputData <- t(assayDataApply(demoData, 1, norm_quant)) @@ -71,7 +72,7 @@ sub_target_demoData <- normalize(sub_target_demoData , data_type="RNA", norm_met # compute negative norm factors for samples and divide by geomean negfactors <- assayDataApply(negativeControlSubset(sub_target_demoData), 2, ngeoMean) norm_neg <- function(x){ - x <- x/(negfactors/geoMean(negfactors))} + x <- x/(negfactors/ngeoMean(negfactors))} expectedOutputData <- t(assayDataApply(sub_target_demoData, 1, norm_neg)) # Extract normalized data from object @@ -100,7 +101,7 @@ test_that("negative norm factors for multipanel are correct", { }) # check the normalized values for multipanel -if (length(unique(fData(target_demoData)[["Module"]]) > 1)){ +if (length(unique(fData(target_demoData)[["Module"]])) > 1){ # compute negative norm factors for samples and divide by geomean # subset data per module pool <- as.list(unique(fData(target_demoData)[["Module"]])) @@ -108,12 +109,11 @@ if (length(unique(fData(target_demoData)[["Module"]]) > 1)){ poolSubSet <- subset(target_demoData, subset = Module == x) negfactors <- assayDataApply(negativeControlSubset(poolSubSet), 2, ngeoMean) norm_neg <- function(x){ - x <- x/(negfactors/geoMean(negfactors)) + x <- x/(negfactors/ngeoMean(negfactors)) } expectedNormMat_i <- t(assayDataApply(poolSubSet, 1, norm_neg)) }) expectedOutputData <- do.call(rbind, expectedNormMat) - expectedOutputData <- expectedOutputData[order(row.names(expectedOutputData)), ] } # Extract normalized data from object @@ -131,7 +131,7 @@ target_demoData <- normalize(target_demoData , data_type="RNA", norm_method="hk" hksubset <- target_demoData[which(featureData(target_demoData)[["TargetName"]] %in% housekeepers),] hkfactors <- assayDataApply(hksubset, 2, ngeoMean) norm_hk <- function(x){ - x <- x/(hkfactors/geoMean(hkfactors)) + x <- x/(hkfactors/ngeoMean(hkfactors)) } expectedOutputData <- t(assayDataApply(target_demoData, 1, norm_hk)) From 5f98b10a608efebe7bf525b7ef608d68f6067e07 Mon Sep 17 00:00:00 2001 From: David Henderson Date: Thu, 19 Aug 2021 11:30:22 -0700 Subject: [PATCH 12/28] Changing from hard coded columns in lsm to selecting by name --- R/NanoStringGeoMxSet-de.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/NanoStringGeoMxSet-de.R b/R/NanoStringGeoMxSet-de.R index ed2451e..2a49806 100644 --- a/R/NanoStringGeoMxSet-de.R +++ b/R/NanoStringGeoMxSet-de.R @@ -100,7 +100,7 @@ mixedModelDE <- function(object, elt = "exprs", modelFormula = NULL, lsm <- lmerTest::ls_means(lmOut, which = groupVar, pairwise = TRUE) } lmOut <- matrix(anova(lmOut)[groupVar, "Pr(>F)"], ncol = 1, dimnames = list(groupVar, "Pr(>F)")) - lsmOut <- matrix(cbind(lsm[,1], lsm[,2]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "PR(>|t|)"))) + lsmOut <- matrix(cbind(lsm[,"Estimate"], lsm[,"Pr(>|t|)"]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "Pr(>|t|)"))) return(list(anova = lmOut, lsmeans = lsmOut)) } mixedOut <- assayDataApply(object, 1, deFunc, groupVar, pDat, formula(paste("expr", as.character(modelFormula)[2], sep = " ~ ")), pairwise, elt = elt) From 33111073997947e7e592f076f49795e5b7a1c134 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 23 Aug 2021 14:05:34 -0700 Subject: [PATCH 13/28] Keep object structure and NAs when thresholding values --- R/utils.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/utils.R b/R/utils.R index 4fc0998..e44aaa5 100644 --- a/R/utils.R +++ b/R/utils.R @@ -57,7 +57,7 @@ thresholdValues <- function(x, thresh=0.5) { thresh <- 0.5 } if (min(x, na.rm = TRUE) < thresh) { - x <- x[!is.na(x)] + thresh + x <- x + thresh } return(x) } From e215ff4f126b43ee2bb72a2ebfa7618756f3f6ad Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 23 Aug 2021 15:38:16 -0700 Subject: [PATCH 14/28] Updated to current commit --- NEWS.md | 54 +++++++++++++++++++++++++++++++++--------------------- 1 file changed, 33 insertions(+), 21 deletions(-) diff --git a/NEWS.md b/NEWS.md index 2d24ec1..04da62a 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,19 +1,40 @@ +# GeomxTools 1.1.2 + +## Revisions: +* Allow users to use more than one DCC version +* Speed improvements in `setBioProbeQC()` and `aggregateCounts()` + +## Bug fixes: +* Fix error in `setBioProbeQC()` with single panel objects +* Fix mixed model output to reference correct p-value +* Fix thresholding in utility functions to keep format matrix format inputs + # GeomxTools 1.1.1 -* +## New features: +* Differential expression with linear mixed model method `mixedModelDE()` +## Revisions: +* Handle multi-panel normalization +## Bug fixes: +* Fix skipping of vectors that don't meet Grubbs requirements +* Fix build warning from knitr update -# GeomxTools 0.99.4 - concomittant development branch version (includes changes beyond 1.0.0) +# GeomxTools 0.99.4 - concomittant development branch version +# Includes changes beyond 1.0.0 ## New features: -* -* Segment `setSegmentQCFlags()` and probe `setBioProbeQC()` QC function -* - +* New slot FeatureType to indicate if data is probe- or target-level +* Segment QC `setSegmentQCFlags()` and probe QC `setBioProbeQC()` +* Count aggregation method `aggregateCounts()` +* Common GeoMx normalizations `normalize()` +* Log and count thresholding methods added to utils -## Revision: -* Change read +## Revisions: +* Updated `readDccFile()` to expand dcc file versions accepted +* Allow user to without auto-aggregating counts to target-level +* Probe annotations attached to featureData with readPKCFile # GeomxTools 1.1.0 @@ -21,20 +42,11 @@ # GeomxTools 1.0.0 -* Initial release version, includes load with automatic aggregation of counts to target level. - -# GeomxTools 0.99.3 - -* Reverse change in R version requirement - -# GeomxTools 0.99.2 - -* Update R version requirement - -# GeomxTools 0.99.1 +* Initial release, includes load with automatic aggregation to target-level -* Update R version requirement +## User notes: +* This version was included in Bioconductor release 3.13 # GeomxTools 0.99.0 -* Package template creation \ No newline at end of file +* Package template creation From 893061a8c744c6e3c4b7d48b9d974ecdb0ce5393 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 23 Aug 2021 18:39:03 -0700 Subject: [PATCH 15/28] Fix p-val in multicore de and update test to check p-val --- R/NanoStringGeoMxSet-de.R | 2 +- tests/testthat/test_de.R | 188 +++++++++++++++++++++----------------- 2 files changed, 104 insertions(+), 86 deletions(-) diff --git a/R/NanoStringGeoMxSet-de.R b/R/NanoStringGeoMxSet-de.R index 2a49806..204183e 100644 --- a/R/NanoStringGeoMxSet-de.R +++ b/R/NanoStringGeoMxSet-de.R @@ -71,7 +71,7 @@ mixedModelDE <- function(object, elt = "exprs", modelFormula = NULL, lsm <- lmerTest::ls_means(lmOut, which = groupVar, pairwise = TRUE) } lmOut <- matrix(anova(lmOut)[groupVar, "Pr(>F)"], ncol = 1, dimnames = list(groupVar, "Pr(>F)")) - lsmOut <- matrix(cbind(lsm[,1], lsm[,2]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "PR(>|t|)"))) + lsmOut <- matrix(cbind(lsm[,"Estimate"], lsm[,"Pr(>|t|)"]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "PR(>|t|)"))) return(list(anova = lmOut, lsmeans = lsmOut)) } diff --git a/tests/testthat/test_de.R b/tests/testthat/test_de.R index 6651a3a..f756a83 100644 --- a/tests/testthat/test_de.R +++ b/tests/testthat/test_de.R @@ -2,7 +2,6 @@ library(lmerTest) library(stringr) library(testthat) - datadir <- system.file("extdata", "DSP_NGS_Example_Data", package = "GeomxTools") demoData <- readRDS(file.path(datadir, "/demoData.rds")) target_demoData <- aggregateCounts(demoData) @@ -13,100 +12,119 @@ pData(target_demoData)[["slide"]]<- factor(pData(target_demoData)[["slide name"] protocolData(target_demoData)[["cell_line"]] <- factor(protocolData(target_demoData)[["cell_line"]]) protocolData(target_demoData)[["pool_rep"]] <- factor(protocolData(target_demoData)[["pool_rep"]]) -test_that("test that the outputs of methods are same", { -# Run multiple cores, fixed effect and random intercept -mixedOutmc <- mixedModelDE(target_demoData, - elt = "exprs_norm", - modelFormula = ~ pool_rep + (1 | slide), - groupVar = "pool_rep", - nCores = 12, - multiCore = TRUE, - pAdjust = NULL -) - -# Run parallel -mixedOutp <- mixedModelDE(target_demoData, - elt = "exprs_norm", - modelFormula = ~ pool_rep + (1 | slide), - groupVar = "pool_rep", - nCores = 12, - multiCore = FALSE, - pAdjust = NULL -) - # Run Serial mixedOuts <- mixedModelDE(target_demoData, - elt = "exprs_norm", - modelFormula = ~ pool_rep + (1 | slide), - groupVar = "pool_rep", - nCores = 1, - multiCore = FALSE, - pAdjust = NULL + elt = "exprs_norm", + modelFormula = ~ pool_rep + (1 | slide), + groupVar = "pool_rep", + nCores = 1, + multiCore = FALSE, + pAdjust = NULL ) - -# req 1: test that the outputs from mc and parallel methods is the same:------ - expect_equal(mixedOutmc, mixedOutp) - - -# req 2: test that the outputs from mc and serial methods is the same:------ - expect_equal(mixedOutmc, mixedOuts) +test_that("test output p-value is as expected", { +# Randomly test one target manually +testTargetExprs <- + assayDataElement(target_demoData, elt="exprs_norm") +testTargetExprs <- + t(testTargetExprs[sample(nrow(testTargetExprs), 1), , drop=FALSE]) +sampleExprs <- data.frame(testTargetExprs, + sData(target_demoData)[, c("pool_rep", "slide")]) +target <- names(sampleExprs)[1L] +testMod <- + suppressMessages(lmerTest::lmer(paste0(target, " ~ pool_rep + (1 | slide)"), + sampleExprs)) +manualResult <- lmerTest::ls_means(testMod, which="pool_rep", pairwise=TRUE) +expect_equal(mixedOuts["lsmeans", target][[1]][, "PR(>|t|)"], + manualResult[, "Pr(>|t|)"]) }) -test_that("mixed model function works for random slope and random intercept", { - -# req 3. test that function works for model with random slope and random intercept -design(target_demoData) <- ~ pool_rep + (1 + slide| slide) -# Run multiple cores -mixedOutmc <- mixedModelDE(target_demoData, - elt = "exprs_norm", - modelFormula = NULL, - groupVar = "pool_rep", - nCores = 12, - multiCore = TRUE, - pAdjust = NULL -) +if (Sys.info()['sysname'] != "Windows") { -# Run parallel -mixedOutp <- mixedModelDE(target_demoData, - elt = "exprs_norm", - modelFormula = NULL, - groupVar = "pool_rep", - nCores = 12, - multiCore = FALSE, - pAdjust = NULL -) + test_that("test that the outputs of methods are same", { + # Run multiple cores, fixed effect and random intercept + mixedOutmc <- mixedModelDE(target_demoData, + elt = "exprs_norm", + modelFormula = ~ pool_rep + (1 | slide), + groupVar = "pool_rep", + nCores = 12, + multiCore = TRUE, + pAdjust = NULL + ) - expect_equal(mixedOutmc, mixedOutp) -}) - - -# req 4: test that function outputs an error when model terms are not in sData:------ -# Run parallel -test_that("test that error when model terms are not in sdata", { - expect_error( + # Run parallel mixedOutp <- mixedModelDE(target_demoData, - elt = "exprs_norm", - modelFormula = ~ cell_line + (1 | fakevar) , - groupVar = "cell_line", - nCores = 12, - multiCore = FALSE, - pAdjust = NULL) - ) -}) - - -# req 5: test that function outputs an error when groupVar is not in model formula:------ -# Run parallel -test_that("test that error when group var not in model", { - expect_error( + elt = "exprs_norm", + modelFormula = ~ pool_rep + (1 | slide), + groupVar = "pool_rep", + nCores = 12, + multiCore = FALSE, + pAdjust = NULL + ) + + # req 1: test that the outputs from mc and parallel methods is the same:------ + expect_equal(mixedOutmc, mixedOutp) + + + # req 2: test that the outputs from mc and serial methods is the same:------ + expect_equal(mixedOutmc, mixedOuts) + }) + + test_that("mixed model function works for random slope and random intercept", { + + # req 3. test that function works for model with random slope and random intercept + design(target_demoData) <- ~ pool_rep + (1 + slide| slide) + # Run multiple cores + mixedOutmc <- mixedModelDE(target_demoData, + elt = "exprs_norm", + modelFormula = NULL, + groupVar = "pool_rep", + nCores = 12, + multiCore = TRUE, + pAdjust = NULL + ) + + # Run parallel mixedOutp <- mixedModelDE(target_demoData, elt = "exprs_norm", - modelFormula = ~ cell_line + (1 | slide) , - groupVar = "segment", + modelFormula = NULL, + groupVar = "pool_rep", nCores = 12, multiCore = FALSE, - pAdjust = NULL) - ) -}) - + pAdjust = NULL + ) + + expect_equal(mixedOutmc, mixedOutp) + }) + + + # req 4: test that function outputs an error when model terms are not in sData:------ + # Run parallel + test_that("test that error when model terms are not in sdata", { + expect_error( + mixedOutp <- mixedModelDE(target_demoData, + elt = "exprs_norm", + modelFormula = ~ cell_line + (1 | fakevar) , + groupVar = "cell_line", + nCores = 12, + multiCore = FALSE, + pAdjust = NULL) + ) + }) + + + # req 5: test that function outputs an error when groupVar is not in model formula:------ + # Run parallel + test_that("test that error when group var not in model", { + expect_error( + mixedOutp <- mixedModelDE(target_demoData, + elt = "exprs_norm", + modelFormula = ~ cell_line + (1 | slide) , + groupVar = "segment", + nCores = 12, + multiCore = FALSE, + pAdjust = NULL) + ) + }) + +} \ No newline at end of file From dbd691c5df90db29eb92cb4e2e636c1f2add4867 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 23 Aug 2021 19:05:40 -0700 Subject: [PATCH 16/28] Fix capitalization --- R/NanoStringGeoMxSet-de.R | 2 +- tests/testthat/test_de.R | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/R/NanoStringGeoMxSet-de.R b/R/NanoStringGeoMxSet-de.R index 204183e..f19ff61 100644 --- a/R/NanoStringGeoMxSet-de.R +++ b/R/NanoStringGeoMxSet-de.R @@ -71,7 +71,7 @@ mixedModelDE <- function(object, elt = "exprs", modelFormula = NULL, lsm <- lmerTest::ls_means(lmOut, which = groupVar, pairwise = TRUE) } lmOut <- matrix(anova(lmOut)[groupVar, "Pr(>F)"], ncol = 1, dimnames = list(groupVar, "Pr(>F)")) - lsmOut <- matrix(cbind(lsm[,"Estimate"], lsm[,"Pr(>|t|)"]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "PR(>|t|)"))) + lsmOut <- matrix(cbind(lsm[,"Estimate"], lsm[,"Pr(>|t|)"]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "Pr(>|t|)"))) return(list(anova = lmOut, lsmeans = lsmOut)) } diff --git a/tests/testthat/test_de.R b/tests/testthat/test_de.R index f756a83..e27229d 100644 --- a/tests/testthat/test_de.R +++ b/tests/testthat/test_de.R @@ -35,7 +35,7 @@ testMod <- suppressMessages(lmerTest::lmer(paste0(target, " ~ pool_rep + (1 | slide)"), sampleExprs)) manualResult <- lmerTest::ls_means(testMod, which="pool_rep", pairwise=TRUE) -expect_equal(mixedOuts["lsmeans", target][[1]][, "PR(>|t|)"], +expect_equal(mixedOuts["lsmeans", target][[1]][, "Pr(>|t|)"], manualResult[, "Pr(>|t|)"]) }) From c5bb3696339e8f170e3f6e0fb2183d79111d2e77 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 24 Aug 2021 13:41:53 -0700 Subject: [PATCH 17/28] Remove duplicated neg mean --- R/NanoStringGeoMxSet-qc.R | 1 - tests/testthat/test_QC.R | 7 ++++--- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index ae141c3..007f1b5 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -168,7 +168,6 @@ setLowNegFlags <- function(object, cutoff=DEFAULTS[["minNegativeCount"]]) { FUN=function( x ) { assayDataApply( x, MARGIN = 2, FUN=ngeoMean, elt="exprs" ) }) - protocolData(object)[["NegGeoMean"]] <- negativeGeoMeans lowNegs <- data.frame("LowNegatives"=apply(negativeGeoMeans < cutoff, 1, sum) > 0) object <- appendSampleFlags(object, lowNegs) diff --git a/tests/testthat/test_QC.R b/tests/testthat/test_QC.R index 49b6785..cd81b45 100644 --- a/tests/testthat/test_QC.R +++ b/tests/testthat/test_QC.R @@ -61,14 +61,15 @@ testData <- setBackgroundQCFlags(testData, qcCutoffs=list(minNegativeCount=10, maxNTCCount=60)) - +testData <- summarizeNegatives(testData) prData <- protocolData(testData) TechBgQC <- as.data.frame(prData[["QCFlags"]]) # req 1: test that the number of Low Negatives is correct:------ testthat::test_that("test that the number of Low Negatives is correct", { - expect_true(sum(TechBgQC$LowNegatives) == - sum(apply(prData@data$NegGeoMean < 10, 1, sum) > 0)) + numberFail <- (pData(testData)[, "NegGeoMean_VnV_GeoMx_Hs_CTA_v1.2"] < 10) + + (pData(testData)[, "NegGeoMean_Six-gene_test_v1_v1.1"] < 10) + expect_true(sum(TechBgQC$LowNegatives) == sum(numberFail > 0)) }) From 69caee052e5d7bb2cf05bde407413d54daad1033 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 24 Aug 2021 16:56:46 -0700 Subject: [PATCH 18/28] Ignore note since no negs in test --- tests/testthat/test_GeoMxSet_probeCollapse.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/testthat/test_GeoMxSet_probeCollapse.R b/tests/testthat/test_GeoMxSet_probeCollapse.R index eb93e38..a4b2871 100644 --- a/tests/testthat/test_GeoMxSet_probeCollapse.R +++ b/tests/testthat/test_GeoMxSet_probeCollapse.R @@ -18,7 +18,7 @@ subData <- sample(unique(fData(subData)[["TargetName"]]), 10, replace=FALSE)) -subAggd <- aggregateCounts(subData) +subAggd <- suppressWarnings(aggregateCounts(subData)) testthat::test_that("Feature type changed after aggregation", { expect_true(featureType(subData) == "Probe") @@ -55,7 +55,7 @@ testthat::test_that("Target expression matrix contains aggregated counts", { }) testthat::test_that("Other aggregation functions work", { - subSum <- aggregateCounts(subData, FUN=sum) + subSum <- suppressWarnings(aggregateCounts(subData, FUN=sum)) expect_true(all(colnames(exprs(subData)) == colnames(exprs(subSum)))) sameTargs <- intersect(fData(subData)[["TargetName"]], rownames(exprs(subSum))) From d5645c4205e00686c9e0b94490529c27ca8c4f75 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 24 Aug 2021 17:33:40 -0700 Subject: [PATCH 19/28] Avoid import warning --- NAMESPACE | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 3cac3e1..121a0f2 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -2,9 +2,8 @@ ### Imports import(S4Vectors) import(Biobase) -import(data.table) importClassesFrom(NanoStringNCTools, SignatureSet) -importClassesFrom(S4Vectors,DataFrame) +importClassesFrom(S4Vectors, DataFrame) importFrom(NanoStringNCTools, SignatureSet) importFrom(rjson, fromJSON) importFrom(readxl, read_xlsx) @@ -20,6 +19,7 @@ importFrom(methods, is) importFrom(methods, validObject) importFrom(utils, write.table) importFrom(outliers, grubbs.test) +importFrom(data.table, data.table, .SD) importFrom(BiocGenerics, design) importFrom(BiocGenerics, "design<-") From a267de44580615582aea8eddf77803c1985bd3fc Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 24 Aug 2021 17:51:00 -0700 Subject: [PATCH 20/28] import from nctools --- NAMESPACE | 1 + 1 file changed, 1 insertion(+) diff --git a/NAMESPACE b/NAMESPACE index 121a0f2..68737fd 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -5,6 +5,7 @@ import(Biobase) importClassesFrom(NanoStringNCTools, SignatureSet) importClassesFrom(S4Vectors, DataFrame) importFrom(NanoStringNCTools, SignatureSet) +importFrom(NanoStringNCTools, normalize) importFrom(rjson, fromJSON) importFrom(readxl, read_xlsx) importFrom(EnvStats, geoMean) From 197efbc2a195e087c84b2a9e44ad8beb1488a327 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Tue, 24 Aug 2021 18:08:01 -0700 Subject: [PATCH 21/28] Add imports for DE --- DESCRIPTION | 2 +- NAMESPACE | 6 ++++++ 2 files changed, 7 insertions(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 6370ba5..5cc9c8f 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -8,7 +8,7 @@ Authors@R: c(person("Nicole", "Ortogero", email = "nortogero@nanostring.com", ro person("Zhi", "Yang", email = "zyang@nanostring.com", role = c("aut"))) Depends: R (>= 3.6), NanoStringNCTools Imports: Biobase, S4Vectors, rjson, readxl, EnvStats, reshape2, methods, - utils, stats, data.table, outliers, BiocGenerics + utils, stats, data.table, outliers, lmerTest, parallel, BiocGenerics Suggests: rmarkdown, knitr, diff --git a/NAMESPACE b/NAMESPACE index 68737fd..9b1c679 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -21,6 +21,12 @@ importFrom(methods, validObject) importFrom(utils, write.table) importFrom(outliers, grubbs.test) importFrom(data.table, data.table, .SD) +importFrom(lmerTest, lmer) +importFrom(lmerTest, ls_means) +importFrom(parallel, mclapply) +importFrom(parallel, parLapply) +importFrom(parallel, makeCluster) +importFrom(parallel, stopCluster) importFrom(BiocGenerics, design) importFrom(BiocGenerics, "design<-") From 28e37e9777da8575703b16cf311140af91c728df Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 25 Aug 2021 12:18:38 -0700 Subject: [PATCH 22/28] Fix documentations iss104 --- NAMESPACE | 8 +++--- R/NanoStringGeoMxSet-de.R | 1 + R/NanoStringGeoMxSet-normalize.R | 7 +++-- R/NanoStringGeoMxSet-qc.R | 7 ++++- R/NanoStringGeoMxSet-signatures.R | 6 ++-- R/NanoStringGeoMxSet-validity.R | 2 +- R/utils.R | 2 ++ man/NanoStringGeoMxSet-class.Rd | 16 +++++++++++ man/normalize-NanoStringGeoMxSet-method.Rd | 5 +++- man/readNanoStringGeoMxSet.Rd | 4 +-- man/setQCFlags-NanoStringGeoMxSet-method.Rd | 31 --------------------- 11 files changed, 44 insertions(+), 45 deletions(-) delete mode 100644 man/setQCFlags-NanoStringGeoMxSet-method.Rd diff --git a/NAMESPACE b/NAMESPACE index 9b1c679..6d8c404 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -2,10 +2,8 @@ ### Imports import(S4Vectors) import(Biobase) -importClassesFrom(NanoStringNCTools, SignatureSet) +import(NanoStringNCTools) importClassesFrom(S4Vectors, DataFrame) -importFrom(NanoStringNCTools, SignatureSet) -importFrom(NanoStringNCTools, normalize) importFrom(rjson, fromJSON) importFrom(readxl, read_xlsx) importFrom(EnvStats, geoMean) @@ -13,7 +11,9 @@ importFrom(reshape2, dcast) importFrom(utils, read.csv) importFrom(stats, as.formula) importFrom(stats, quantile) -importFrom(stats, weights) +importFrom(stats, anova) +importFrom(stats, formula) +importFrom(stats, p.adjust) importFrom(methods, callGeneric) importFrom(methods, callNextMethod) importFrom(methods, is) diff --git a/R/NanoStringGeoMxSet-de.R b/R/NanoStringGeoMxSet-de.R index f19ff61..3791ffb 100644 --- a/R/NanoStringGeoMxSet-de.R +++ b/R/NanoStringGeoMxSet-de.R @@ -10,6 +10,7 @@ #' @param nCores = 1, number of cores to use, set to 1 if running in serial mode #' @param multiCore = TRUE, set to TRUE to use multiCore, FALSE to run in cluster mode #' @param pAdjust = "BY" method for p-value adjustment +#' @param pairwise boolean to calculate least-square means pairwise differences #' #' @return mixed model output list #' diff --git a/R/NanoStringGeoMxSet-normalize.R b/R/NanoStringGeoMxSet-normalize.R index a442029..23f40e1 100644 --- a/R/NanoStringGeoMxSet-normalize.R +++ b/R/NanoStringGeoMxSet-normalize.R @@ -9,6 +9,8 @@ HOUSEKEEPERS <- c( #' @param data_type the data type of the object. Values maybe RNA, protein. #' @param fromElt name of the assayDataElement to normalize #' @param toElt name of the assayDataElement to store normalized values +#' @param housekeepers optional vector of housekeeper target names +#' @param \code{ldots} optional arguments #' @return a NanoStringGeoMxSet object with normalized counts and normalized factors #' @examples #' datadir <- system.file("extdata", "DSP_NGS_Example_Data", @@ -149,7 +151,7 @@ hkNorm <- function(object, data_type, toElt, fromElt, housekeepers) { # subtract background subtractBackground <- function(object, data_type, toElt, fromElt) { if (!featureType(object) == "Target") { - negsubset <- subset(object, subset = Codeclass %in% c("Negative01", "Negative")) + negsubset <- subset(object, subset = CodeClass %in% c("Negative01", "Negative")) negs <- apply(exprs(negsubset), 2, function(x) ngeoMean(x)) assayDataElement(object, toElt) <- t(assayDataApply(object, MARGIN = 1L, FUN = `-`, t(negs), elt = fromElt)) @@ -181,8 +183,9 @@ setGeneric("checkQCFlags", #' checkQCFlags -#' @param NanoStringGeoMxSet +#' @param object name of the NanoStringGeoMxSet object to check the QC Flags #' @param removeLowLocalOutliers logical, if TRUE it sets outlier counts to zero, default is FALSE, +#' @param \code{ldots} optional arguments #' @return NanoStringGeoMxSet #' @export #' diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 007f1b5..ea5cfb5 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -531,6 +531,7 @@ appendFeatureFlags <- function(object, currFlags) { #' @return the object that QC was performed on #' #' @examples +#' setQCFlags(object) #' setMethod("setQCFlags", signature(object="NanoStringGeoMxSet"), @@ -544,7 +545,7 @@ setMethod("setQCFlags", } else if (featureType(object) == "Target") { # object <- setTargetFlags(object=object, qcCutoffs=qcCutoffs) } else { - valid(Object(x)) + valid(object) } return(object) }) @@ -562,6 +563,10 @@ setTargetFlags <- function(object, qcCutoffs=DEFAULTS) { return(object) } +setHighCountFlags <- function(object, qcCutoffs=DEFAULTS) { + object +} + setLOQFlags <- function(object, cutoff=DEFAULTS[["loqCutoff"]]) { if (featureType(object) == "Target") { negativeObject <- negativeControlSubset(object) diff --git a/R/NanoStringGeoMxSet-signatures.R b/R/NanoStringGeoMxSet-signatures.R index 2bd1a0d..83396e9 100644 --- a/R/NanoStringGeoMxSet-signatures.R +++ b/R/NanoStringGeoMxSet-signatures.R @@ -23,13 +23,13 @@ setMethod("signatureScores", "NanoStringGeoMxSet", exprs <- t(assayDataElement2(object, elt)) colnames(exprs) <- featureData(object)[["TargetName"]] sigFuncList <- signatureFuncs( object ) - linWeights <- stats::weights( signatures( object ) )[names( sigFuncList )[which( sigFuncList %in% "default" )]] + linWeights <- weights( signatures( object ) )[names( sigFuncList )[which( sigFuncList %in% "default" )]] nonLinFuncs <- sigFuncList[which( !( sigFuncList %in% "default" ) )] scores <- .sigCalc(exprs, linWeights) while (length(idx <- which(rowSums(is.na(scores)) > 0L))) { subscores <- .sigCalc(cbind(exprs, t(scores[-idx, , drop = FALSE])), - stats::weights(signatures(object))[idx]) + weights(signatures(object))[idx]) if (all(is.na(subscores))) break else @@ -39,7 +39,7 @@ setMethod("signatureScores", "NanoStringGeoMxSet", if (ncol(nonLinScores) > 0) { scores <- rbind( scores , nonLinScores ) } - return( scores[names( stats::weights( signatures( object ) ) ),] ) + return( scores[names( weights( signatures( object ) ) ),] ) }) setGeneric( "signatureGroups" , signature = "object" , diff --git a/R/NanoStringGeoMxSet-validity.R b/R/NanoStringGeoMxSet-validity.R index 9041fa3..2bd1016 100644 --- a/R/NanoStringGeoMxSet-validity.R +++ b/R/NanoStringGeoMxSet-validity.R @@ -72,7 +72,7 @@ function(object) if (any(numTargets == 0L)) { msg <- c(msg, "'signatures' vectors must be non-empty") } else { - targets <- names(unlist(unname(stats::weights(signatures(object))))) + targets <- names(unlist(unname(weights(signatures(object))))) if (is.null(targets) || any(nchar(targets) == 0L)) { msg <- c(msg, "'signatures' vectors must be named") } else if(!all(unique(targets) %in% diff --git a/R/utils.R b/R/utils.R index e44aaa5..e5d1b98 100644 --- a/R/utils.R +++ b/R/utils.R @@ -67,6 +67,8 @@ thresholdValues <- function(x, thresh=0.5) { #' @param object name of the NanoStringGeoMxSet object #' @param elt expression matrix element in \code{assayDataElement} #' to shift all counts by +#' @param useDALogic boolean to use the same logic in DA (impute 0s to 1s) +#' or set to FALSE to shift all counts by 1 #' #' @return object of NanoStringGeoMxSet class #' diff --git a/man/NanoStringGeoMxSet-class.Rd b/man/NanoStringGeoMxSet-class.Rd index e69db89..97b010f 100644 --- a/man/NanoStringGeoMxSet-class.Rd +++ b/man/NanoStringGeoMxSet-class.Rd @@ -25,6 +25,10 @@ \alias{svarLabels,NanoStringGeoMxSet-method} \alias{dimLabels,NanoStringGeoMxSet-method} \alias{dimLabels<-,NanoStringGeoMxSet,character-method} +\alias{featureType} +\alias{featureType,NanoStringGeoMxSet-method} +\alias{featureType<-} +\alias{featureType<-,NanoStringGeoMxSet,character-method} \alias{signatures} \alias{signatures,NanoStringGeoMxSet-method} \alias{signatures<-} @@ -56,6 +60,7 @@ NanoStringGeoMxSet(assayData, dimLabels=c("TargetName", "SampleID"), signatures=SignatureSet(), design=NULL, + featureType="Probe", \ldots) } @@ -80,6 +85,8 @@ NanoStringGeoMxSet(assayData, signature definitions.} \item{design}{An optional one-sided formula representing the experimental design based on columns from \code{\link[Biobase]{phenoData}}} + \item{featureType}{A \code{character} string indicating if features are on + \code{"Probe"} or \code{"Target"} level} \item{\ldots}{Additional arguments for \code{\link{ExpressionSet}}.} } @@ -103,6 +110,14 @@ NanoStringGeoMxSet(assayData, \code{dimLabels(object) <- value}: replaces the \code{dimLabels} of the \code{object}. } + \item{}{ + \code{featureType(object)}: extracts the \code{featureType} + of the \code{object}. + } + \item{}{ + \code{featureType(object) <- value}: replaces the \code{featureType} + of the \code{object}. + } \item{}{ \code{signatures(object)}: extracts the \code{\link{SignatureSet}} of the \code{object}. @@ -172,6 +187,7 @@ dccSet <- readNanoStringGeoMxSet(dccFiles=dccFiles, # Accessing sample data and column names head(sData(dccSet)) svarLabels(dccSet) +featureType(dccSet) # Accessing number of samples and features dim(dccSet) diff --git a/man/normalize-NanoStringGeoMxSet-method.Rd b/man/normalize-NanoStringGeoMxSet-method.Rd index 3c13a0b..515432c 100644 --- a/man/normalize-NanoStringGeoMxSet-method.Rd +++ b/man/normalize-NanoStringGeoMxSet-method.Rd @@ -24,6 +24,8 @@ \item{fromElt}{name of the assayDataElement to normalize} \item{toElt}{name of the assayDataElement to store normalized values} + +\item{housekeepers}{optional vector of housekeeper target names} } \value{ a NanoStringGeoMxSet object with normalized counts and normalized factors @@ -33,7 +35,8 @@ normalize GeoMxSet using different normalization methods } \examples{ datadir <- system.file("extdata", "DSP_NGS_Example_Data", - package="GeomxTools") + package = "GeomxTools" +) demoData <- readRDS(file.path(datadir, "/demoData.rds")) norm_object <- normalize(demoData) } diff --git a/man/readNanoStringGeoMxSet.Rd b/man/readNanoStringGeoMxSet.Rd index 05bad16..a9e339a 100644 --- a/man/readNanoStringGeoMxSet.Rd +++ b/man/readNanoStringGeoMxSet.Rd @@ -12,8 +12,8 @@ \usage{ readNanoStringGeoMxSet(dccFiles, pkcFiles, phenoDataFile, phenoDataSheet, phenoDataDccColName = "Sample_ID", - phenoDataColPrefix = "", protocolDataColNames = c("slide name"), - experimentDataColNames = c("panel")) + phenoDataColPrefix = "", protocolDataColNames = NULL, + experimentDataColNames = NULL) } \arguments{ diff --git a/man/setQCFlags-NanoStringGeoMxSet-method.Rd b/man/setQCFlags-NanoStringGeoMxSet-method.Rd deleted file mode 100644 index 9ff42f5..0000000 --- a/man/setQCFlags-NanoStringGeoMxSet-method.Rd +++ /dev/null @@ -1,31 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/NanoStringGeoMxSet-qc.R -\name{setQCFlags,NanoStringGeoMxSet-method} -\alias{setQCFlags,NanoStringGeoMxSet-method} -\title{Add QC flags to feature or protocol data} -\usage{ -\S4method{setQCFlags}{NanoStringGeoMxSet}(object, qcCutoffs = DEFAULTS, ...) -} -\arguments{ -\item{object}{name of the object class to perform QC on -\enumerate{ - \item{NanoStringGeoMxSet, use the NanoStringGeoMxSet class} -}} - -\item{qcCutoffs}{list of cutoffs and thresholds to use for QC} - -\item{dataDim}{the dimension of the object to QC on -\enumerate{ - \item{sample, QC data on the AOI level} - \item{feature, QC data on probe or target level} -}} -} -\value{ -the object that QC was performed on -} -\description{ -Add QC flags to feature or protocol data -} -\examples{ - -} From c83d73651ccad5fa0b70ea1b8dbbcd97fe90cf61 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 25 Aug 2021 13:14:34 -0700 Subject: [PATCH 23/28] Man updates --- R/NanoStringGeoMxSet-qc.R | 12 +++---- man/checkQCFlags-NanoStringGeoMxSet-method.Rd | 7 +++-- man/checkQCFlags.Rd | 3 +- man/mixedModelDE.Rd | 2 ++ man/normalize-NanoStringGeoMxSet-method.Rd | 2 ++ man/setQCFlags-NanoStringGeoMxSet-method.Rd | 31 +++++++++++++++++++ man/shiftCountsOne.Rd | 3 ++ 7 files changed, 51 insertions(+), 9 deletions(-) create mode 100644 man/setQCFlags-NanoStringGeoMxSet-method.Rd diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index ea5cfb5..26d3951 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -515,24 +515,24 @@ appendFeatureFlags <- function(object, currFlags) { } ############Not used########################################################## -#' Add QC flags to feature or protocol data +#' Add QC flags to feature and protocol data simultaneously #' #' @param object name of the object class to perform QC on #' \enumerate{ #' \item{NanoStringGeoMxSet, use the NanoStringGeoMxSet class} #' } -#' @param dataDim the dimension of the object to QC on -#' \enumerate{ -#' \item{sample, QC data on the AOI level} -#' \item{feature, QC data on probe or target level} -#' } #' @param qcCutoffs list of cutoffs and thresholds to use for QC +#' @param \code{ldots} optional parameters to pass #' #' @return the object that QC was performed on #' #' @examples +#' datadir <- system.file("extdata", "DSP_NGS_Example_Data", +#' package="GeomxTools") +#' demoData <- readRDS(file.path(datadir, "/demoData.rds")) #' setQCFlags(object) #' +#' @export setMethod("setQCFlags", signature(object="NanoStringGeoMxSet"), function(object, qcCutoffs=DEFAULTS, ...) { diff --git a/man/checkQCFlags-NanoStringGeoMxSet-method.Rd b/man/checkQCFlags-NanoStringGeoMxSet-method.Rd index 1ca64f6..0a245ef 100644 --- a/man/checkQCFlags-NanoStringGeoMxSet-method.Rd +++ b/man/checkQCFlags-NanoStringGeoMxSet-method.Rd @@ -7,9 +7,11 @@ \S4method{checkQCFlags}{NanoStringGeoMxSet}(object, removeLowLocalOutliers = FALSE, ...) } \arguments{ +\item{object}{name of the NanoStringGeoMxSet object to check the QC Flags} + \item{removeLowLocalOutliers}{logical, if TRUE it sets outlier counts to zero, default is FALSE,} -\item{NanoStringGeoMxSet}{} +\item{\code{ldots}}{optional arguments} } \value{ NanoStringGeoMxSet @@ -19,7 +21,8 @@ checkQCFlags } \examples{ datadir <- system.file("extdata", "DSP_NGS_Example_Data", - package="GeomxTools") + package = "GeomxTools" +) demoData <- readRDS(file.path(datadir, "/demoData.rds")) QCobject <- checkQCFlags(demoData) } diff --git a/man/checkQCFlags.Rd b/man/checkQCFlags.Rd index 9b1e1c1..503bbb8 100644 --- a/man/checkQCFlags.Rd +++ b/man/checkQCFlags.Rd @@ -19,7 +19,8 @@ Check QC Flags in the GeoMxSet and removes the probe or sample from the object } \examples{ datadir <- system.file("extdata", "DSP_NGS_Example_Data", - package="GeomxTools") + package = "GeomxTools" +) demoData <- readRDS(file.path(datadir, "/demoData.rds")) QCobject <- checkQCFlags(demoData) } diff --git a/man/mixedModelDE.Rd b/man/mixedModelDE.Rd index f4b0a50..00d6da5 100644 --- a/man/mixedModelDE.Rd +++ b/man/mixedModelDE.Rd @@ -32,6 +32,8 @@ mixedModelDE( \item{multiCore}{= TRUE, set to TRUE to use multiCore, FALSE to run in cluster mode} \item{pAdjust}{= "BY" method for p-value adjustment} + +\item{pairwise}{boolean to calculate least-square means pairwise differences} } \value{ mixed model output list diff --git a/man/normalize-NanoStringGeoMxSet-method.Rd b/man/normalize-NanoStringGeoMxSet-method.Rd index 515432c..ef05e16 100644 --- a/man/normalize-NanoStringGeoMxSet-method.Rd +++ b/man/normalize-NanoStringGeoMxSet-method.Rd @@ -26,6 +26,8 @@ \item{toElt}{name of the assayDataElement to store normalized values} \item{housekeepers}{optional vector of housekeeper target names} + +\item{\code{ldots}}{optional arguments} } \value{ a NanoStringGeoMxSet object with normalized counts and normalized factors diff --git a/man/setQCFlags-NanoStringGeoMxSet-method.Rd b/man/setQCFlags-NanoStringGeoMxSet-method.Rd new file mode 100644 index 0000000..10b1838 --- /dev/null +++ b/man/setQCFlags-NanoStringGeoMxSet-method.Rd @@ -0,0 +1,31 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/NanoStringGeoMxSet-qc.R +\name{setQCFlags,NanoStringGeoMxSet-method} +\alias{setQCFlags,NanoStringGeoMxSet-method} +\title{Add QC flags to feature and protocol data simultaneously} +\usage{ +\S4method{setQCFlags}{NanoStringGeoMxSet}(object, qcCutoffs = DEFAULTS, ...) +} +\arguments{ +\item{object}{name of the object class to perform QC on +\enumerate{ + \item{NanoStringGeoMxSet, use the NanoStringGeoMxSet class} +}} + +\item{qcCutoffs}{list of cutoffs and thresholds to use for QC} + +\item{\code{ldots}}{optional parameters to pass} +} +\value{ +the object that QC was performed on +} +\description{ +Add QC flags to feature and protocol data simultaneously +} +\examples{ +datadir <- system.file("extdata", "DSP_NGS_Example_Data", + package="GeomxTools") +demoData <- readRDS(file.path(datadir, "/demoData.rds")) +setQCFlags(object) + +} diff --git a/man/shiftCountsOne.Rd b/man/shiftCountsOne.Rd index 65ee3b2..42b83eb 100644 --- a/man/shiftCountsOne.Rd +++ b/man/shiftCountsOne.Rd @@ -11,6 +11,9 @@ shiftCountsOne(object, elt = "exprs", useDALogic = FALSE) \item{elt}{expression matrix element in \code{assayDataElement} to shift all counts by} + +\item{useDALogic}{boolean to use the same logic in DA (impute 0s to 1s) +or set to FALSE to shift all counts by 1} } \value{ object of NanoStringGeoMxSet class From fd465b2a7acc94ee1084e6decfd456ed8afdb102 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 25 Aug 2021 15:49:38 -0700 Subject: [PATCH 24/28] Final 1.1.2 document updates --- DESCRIPTION | 10 ++++++---- NAMESPACE | 2 ++ R/NanoStringGeoMxSet-de.R | 11 +++++------ R/NanoStringGeoMxSet-normalize.R | 4 ++-- R/NanoStringGeoMxSet-qc.R | 11 +++++------ man/checkQCFlags-NanoStringGeoMxSet-method.Rd | 2 +- man/normalize-NanoStringGeoMxSet-method.Rd | 2 +- man/setQCFlags-NanoStringGeoMxSet-method.Rd | 4 ++-- 8 files changed, 24 insertions(+), 22 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 5cc9c8f..61a16fa 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -6,13 +6,15 @@ Version: 1.1.2 Encoding: UTF-8 Authors@R: c(person("Nicole", "Ortogero", email = "nortogero@nanostring.com", role = c("cre", "aut")), person("Zhi", "Yang", email = "zyang@nanostring.com", role = c("aut"))) -Depends: R (>= 3.6), NanoStringNCTools -Imports: Biobase, S4Vectors, rjson, readxl, EnvStats, reshape2, methods, - utils, stats, data.table, outliers, lmerTest, parallel, BiocGenerics +Depends: R (>= 3.6), Biobase, NanoStringNCTools, S4Vectors +Imports: BiocGenerics, rjson, readxl, EnvStats, reshape2, methods, + utils, stats, data.table, outliers, lmerTest, dplyr Suggests: rmarkdown, knitr, - testthat (>= 3.0.0) + testthat (>= 3.0.0), + parallel, + ggiraph License: Artistic-2.0 Collate: DccMetadata.R NanoStringGeoMxSet-class.R diff --git a/NAMESPACE b/NAMESPACE index 6d8c404..81d9808 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -4,6 +4,7 @@ import(S4Vectors) import(Biobase) import(NanoStringNCTools) importClassesFrom(S4Vectors, DataFrame) +importClassesFrom(NanoStringNCTools, SignatureSet) importFrom(rjson, fromJSON) importFrom(readxl, read_xlsx) importFrom(EnvStats, geoMean) @@ -18,6 +19,7 @@ importFrom(methods, callGeneric) importFrom(methods, callNextMethod) importFrom(methods, is) importFrom(methods, validObject) +importFrom(dplyr, bind_rows) importFrom(utils, write.table) importFrom(outliers, grubbs.test) importFrom(data.table, data.table, .SD) diff --git a/R/NanoStringGeoMxSet-de.R b/R/NanoStringGeoMxSet-de.R index 3791ffb..3469e33 100644 --- a/R/NanoStringGeoMxSet-de.R +++ b/R/NanoStringGeoMxSet-de.R @@ -62,7 +62,6 @@ mixedModelDE <- function(object, elt = "exprs", modelFormula = NULL, } } if (nCores > 1) { - require(parallel) deFunc <- function(i, groupVar, pDat, modelFormula, exprs, pairwise = TRUE) { dat <- data.frame(expr = exprs$exprs[i, ], pDat) lmOut <- suppressWarnings(lmerTest::lmer(modelFormula, dat)) @@ -79,12 +78,12 @@ mixedModelDE <- function(object, elt = "exprs", modelFormula = NULL, exprs <- new.env() exprs$exprs <- assayDataElement(object, elt = elt) if (multiCore) { - mixedOut <- mclapply(featureNames(object), deFunc, groupVar, pDat, formula(paste("expr", as.character(modelFormula)[2], sep = " ~ ")), exprs, mc.cores = nCores) + mixedOut <- parallel::mclapply(featureNames(object), deFunc, groupVar, pDat, formula(paste("expr", as.character(modelFormula)[2], sep = " ~ ")), exprs, mc.cores = nCores) } else { - cl <- makeCluster(getOption("cl.cores", nCores)) - mixedOut <- parLapply(cl, featureNames(object), deFunc, groupVar, pDat, formula(paste("expr", as.character(modelFormula)[2], sep = " ~ ")), exprs, pairwise) - suppressWarnings(stopCluster(cl)) + cl <- parallel::makeCluster(getOption("cl.cores", nCores)) + mixedOut <- parallel::parLapply(cl, featureNames(object), deFunc, groupVar, pDat, formula(paste("expr", as.character(modelFormula)[2], sep = " ~ ")), exprs, pairwise) + suppressWarnings(parallel::stopCluster(cl)) } mixedOut <- rbind(array(lapply(mixedOut, function(x) x[["anova"]])), array(lapply(mixedOut, function(x) x[["lsmeans"]]))) @@ -100,7 +99,7 @@ mixedModelDE <- function(object, elt = "exprs", modelFormula = NULL, } else { lsm <- lmerTest::ls_means(lmOut, which = groupVar, pairwise = TRUE) } - lmOut <- matrix(anova(lmOut)[groupVar, "Pr(>F)"], ncol = 1, dimnames = list(groupVar, "Pr(>F)")) + lmOut <- matrix(stats::anova(lmOut)[groupVar, "Pr(>F)"], ncol = 1, dimnames = list(groupVar, "Pr(>F)")) lsmOut <- matrix(cbind(lsm[,"Estimate"], lsm[,"Pr(>|t|)"]), ncol = 2, dimnames = list(gsub(groupVar, "", rownames(lsm)), c("Estimate", "Pr(>|t|)"))) return(list(anova = lmOut, lsmeans = lsmOut)) } diff --git a/R/NanoStringGeoMxSet-normalize.R b/R/NanoStringGeoMxSet-normalize.R index 23f40e1..59171f5 100644 --- a/R/NanoStringGeoMxSet-normalize.R +++ b/R/NanoStringGeoMxSet-normalize.R @@ -10,7 +10,7 @@ HOUSEKEEPERS <- c( #' @param fromElt name of the assayDataElement to normalize #' @param toElt name of the assayDataElement to store normalized values #' @param housekeepers optional vector of housekeeper target names -#' @param \code{ldots} optional arguments +#' @param ... optional arguments #' @return a NanoStringGeoMxSet object with normalized counts and normalized factors #' @examples #' datadir <- system.file("extdata", "DSP_NGS_Example_Data", @@ -185,7 +185,7 @@ setGeneric("checkQCFlags", #' checkQCFlags #' @param object name of the NanoStringGeoMxSet object to check the QC Flags #' @param removeLowLocalOutliers logical, if TRUE it sets outlier counts to zero, default is FALSE, -#' @param \code{ldots} optional arguments +#' @param ... optional arguments #' @return NanoStringGeoMxSet #' @export #' diff --git a/R/NanoStringGeoMxSet-qc.R b/R/NanoStringGeoMxSet-qc.R index 26d3951..8e38467 100644 --- a/R/NanoStringGeoMxSet-qc.R +++ b/R/NanoStringGeoMxSet-qc.R @@ -522,7 +522,7 @@ appendFeatureFlags <- function(object, currFlags) { #' \item{NanoStringGeoMxSet, use the NanoStringGeoMxSet class} #' } #' @param qcCutoffs list of cutoffs and thresholds to use for QC -#' @param \code{ldots} optional parameters to pass +#' @param ... optional parameters to pass #' #' @return the object that QC was performed on #' @@ -530,7 +530,7 @@ appendFeatureFlags <- function(object, currFlags) { #' datadir <- system.file("extdata", "DSP_NGS_Example_Data", #' package="GeomxTools") #' demoData <- readRDS(file.path(datadir, "/demoData.rds")) -#' setQCFlags(object) +#' setQCFlags(demoData) #' #' @export setMethod("setQCFlags", @@ -544,8 +544,6 @@ setMethod("setQCFlags", object <- setBioProbeQCFlags(object=object, qcCutoffs=qcCutoffs) } else if (featureType(object) == "Target") { # object <- setTargetFlags(object=object, qcCutoffs=qcCutoffs) - } else { - valid(object) } return(object) }) @@ -563,8 +561,9 @@ setTargetFlags <- function(object, qcCutoffs=DEFAULTS) { return(object) } -setHighCountFlags <- function(object, qcCutoffs=DEFAULTS) { - object +setHighCountFlags <- function(object, cutoff=DEFAULTS[["highCountCutoff"]]) { + cutoff <- checkCutoffs(cutoff) + return(object) } setLOQFlags <- function(object, cutoff=DEFAULTS[["loqCutoff"]]) { diff --git a/man/checkQCFlags-NanoStringGeoMxSet-method.Rd b/man/checkQCFlags-NanoStringGeoMxSet-method.Rd index 0a245ef..4142c76 100644 --- a/man/checkQCFlags-NanoStringGeoMxSet-method.Rd +++ b/man/checkQCFlags-NanoStringGeoMxSet-method.Rd @@ -11,7 +11,7 @@ \item{removeLowLocalOutliers}{logical, if TRUE it sets outlier counts to zero, default is FALSE,} -\item{\code{ldots}}{optional arguments} +\item{...}{optional arguments} } \value{ NanoStringGeoMxSet diff --git a/man/normalize-NanoStringGeoMxSet-method.Rd b/man/normalize-NanoStringGeoMxSet-method.Rd index ef05e16..5c60fda 100644 --- a/man/normalize-NanoStringGeoMxSet-method.Rd +++ b/man/normalize-NanoStringGeoMxSet-method.Rd @@ -27,7 +27,7 @@ \item{housekeepers}{optional vector of housekeeper target names} -\item{\code{ldots}}{optional arguments} +\item{...}{optional arguments} } \value{ a NanoStringGeoMxSet object with normalized counts and normalized factors diff --git a/man/setQCFlags-NanoStringGeoMxSet-method.Rd b/man/setQCFlags-NanoStringGeoMxSet-method.Rd index 10b1838..89a12a7 100644 --- a/man/setQCFlags-NanoStringGeoMxSet-method.Rd +++ b/man/setQCFlags-NanoStringGeoMxSet-method.Rd @@ -14,7 +14,7 @@ \item{qcCutoffs}{list of cutoffs and thresholds to use for QC} -\item{\code{ldots}}{optional parameters to pass} +\item{...}{optional parameters to pass} } \value{ the object that QC was performed on @@ -26,6 +26,6 @@ Add QC flags to feature and protocol data simultaneously datadir <- system.file("extdata", "DSP_NGS_Example_Data", package="GeomxTools") demoData <- readRDS(file.path(datadir, "/demoData.rds")) -setQCFlags(object) +setQCFlags(demoData) } From 1aeae88de60c08059aac57e285306e1671322371 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Wed, 25 Aug 2021 16:04:52 -0700 Subject: [PATCH 25/28] Add line --- R/utils.R | 2 +- tests/testthat/test_GeoMxSet_probeCollapse.R | 2 +- tests/testthat/test_de.R | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/R/utils.R b/R/utils.R index e5d1b98..3337ef3 100644 --- a/R/utils.R +++ b/R/utils.R @@ -126,4 +126,4 @@ collapseCounts <- function(object) { collapsedCounts <- probeCounts[, lapply(.SD, ngeoMean), by=c("TargetName", "Module")] return(collapsedCounts) -} \ No newline at end of file +} diff --git a/tests/testthat/test_GeoMxSet_probeCollapse.R b/tests/testthat/test_GeoMxSet_probeCollapse.R index a4b2871..b4f92e1 100644 --- a/tests/testthat/test_GeoMxSet_probeCollapse.R +++ b/tests/testthat/test_GeoMxSet_probeCollapse.R @@ -69,4 +69,4 @@ testthat::test_that("Other aggregation functions work", { return(all(aggdCounts == targSums)) }) expect_true(all(unlist(subList))) -}) \ No newline at end of file +}) diff --git a/tests/testthat/test_de.R b/tests/testthat/test_de.R index e27229d..239b688 100644 --- a/tests/testthat/test_de.R +++ b/tests/testthat/test_de.R @@ -127,4 +127,4 @@ if (Sys.info()['sysname'] != "Windows") { ) }) -} \ No newline at end of file +} From b80d975a4efb43ef24be1dec334f93fd4a07bfe9 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 30 Aug 2021 14:20:06 -0700 Subject: [PATCH 26/28] Update test note --- tests/testthat/test_GeoMxSet_replacers.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/testthat/test_GeoMxSet_replacers.R b/tests/testthat/test_GeoMxSet_replacers.R index 734e66a..5cca36c 100644 --- a/tests/testthat/test_GeoMxSet_replacers.R +++ b/tests/testthat/test_GeoMxSet_replacers.R @@ -42,7 +42,7 @@ testthat::test_that("test that dimLabels(testData) matches the value assigned", # req 2: test that design(testData) matches the value assigned:------ -testthat::test_that("test that dimLabels(testData) matches the value assigned", { +testthat::test_that("test that design(testData) matches the value assigned", { des <- design(testData) expect_null(des) design(testData) <- "x ~ y" @@ -51,7 +51,7 @@ testthat::test_that("test that dimLabels(testData) matches the value assigned", }) # req 3: test that featureType(testData) matches the value assigned:------ -testthat::test_that("test that dimLabels(testData) matches the value assigned", { +testthat::test_that("test that featureType(testData) matches the value assigned", { featType <- featureType(testData) featureType(testData) <- "Target" expect_true(all(featureType(testData) == "Target")) From 3fccba79f1529a15db8d4946410d61952996a258 Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Mon, 30 Aug 2021 14:43:15 -0700 Subject: [PATCH 27/28] Remove old comments --- tests/testthat/test_readDccFile.R | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/tests/testthat/test_readDccFile.R b/tests/testthat/test_readDccFile.R index 986df2d..6d20919 100644 --- a/tests/testthat/test_readDccFile.R +++ b/tests/testthat/test_readDccFile.R @@ -15,9 +15,8 @@ testthat::test_that("test that the names of DCC files have the right formats", { expect_true(all(names(dccFile$Header) == c("FileVersion", "SoftwareVersion", "Date"))) expect_true(all(names(dccFile$Scan_Attributes) == c("SampleID", "Plate_ID", "Well"))) # QuickBase: "ID", "Plate_ID", "Well" expect_true(all(names(dccFile$NGS_Processing_Attributes) == c("SeqSetId", "tamperedIni", "Raw", "Trimmed", "Stitched", "Aligned", "umiQ30", - "rtsQ30", "DeduplicatedReads"))) - # QuickBase: c("SeqSetId", "tamperedIni", "trimGaloreOpts", "flash2Opts", "umiExtractOpts", "bowtie2Opts", "umiDedupOpts", "Raw", "Trimmed", "Stitched", "Aligned", "umiQ30", "rtsQ30") - expect_true(all(names(dccFile$Code_Summary) == c("RTS_ID", "Count"))) # QuickBase: "RNAID", "Count" + "rtsQ30", "DeduplicatedReads"))) + expect_true(all(names(dccFile$Code_Summary) == c("RTS_ID", "Count"))) expect_true(all(rownames(dccFile$Code_Summary) == dccFile$Code_Summary$RTS_ID)) }) From 24e9b6f545a33045dfe7ab725c2adcdf6f92d40e Mon Sep 17 00:00:00 2001 From: NicoleEO Date: Thu, 2 Sep 2021 14:59:39 -0700 Subject: [PATCH 28/28] Update agg accessing feature data and removing probe specific qc cols --- R/NanoStringGeoMxSet-aggregate.R | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/R/NanoStringGeoMxSet-aggregate.R b/R/NanoStringGeoMxSet-aggregate.R index 7e4b3d1..718b68b 100644 --- a/R/NanoStringGeoMxSet-aggregate.R +++ b/R/NanoStringGeoMxSet-aggregate.R @@ -51,11 +51,12 @@ aggregateCounts <- function(object, FUN=ngeoMean) { targetCounts <- rbind(targetCounts, singleProbeCounts) targetCounts <- targetCounts[unique(fData(object)[["TargetName"]]), ] - targetFeats <- featureData(object)@data + targetFeats <- fData(object) targetFeats <- targetFeats[!duplicated(targetFeats[["TargetName"]]), ] rownames(targetFeats) <- targetFeats[, "TargetName"] - probeColumns <- c("RTS_ID", "QCFlags", "ProbeID") + probeColumns <- c("RTS_ID", "QCFlags", "ProbeID", + "ProbeRatio", "OutlierFrequency") targetFeats <- targetFeats[, !colnames(targetFeats) %in% probeColumns] targetFeats <-