diff --git a/R/ComponentScores.R b/R/ComponentScores.R index 4cff378..6954385 100644 --- a/R/ComponentScores.R +++ b/R/ComponentScores.R @@ -30,6 +30,30 @@ ScoreUsingSavedComponent <- function(seuratObj, componentOrName, fieldName, assa names(geneWeights) <- savedComponent$feature ad <- Seurat::GetAssayData(seuratObj, assay = assayName, layer = layer) + toDrop <- names(geneWeights)[!names(geneWeights) %in% rownames(ad)] + + # Perform the check first for scale.data, since we can add missing genes: + if (length(toDrop) > 0 && layer == 'scale.data') { + dataLayer <- Seurat::GetAssayData(seuratObj, assay = assayName, layer = 'data') + toScale <- toDrop[toDrop %in% rownames(dataLayer)] + if (length(toScale) > 0) { + print(paste0('The following ', length(toScale), ' genes are missing from scale.data and will be added now: ', paste0(toScale, collapse = ','))) + scaled <- Seurat::ScaleData( + dataLayer, + features = toScale, + do.scale = TRUE, + do.center = TRUE + ) + + if (nrow(scaled) != length(toScale)) { + stop('The number of rows after scaling did not match the expected number') + } + + ad <- rbind(ad, scaled) + seuratObj <- suppressWarnings(Seurat::SetAssayData(seuratObj, assay = assayName, layer = layer, new.data = ad)) + } + } + toDrop <- names(geneWeights)[!names(geneWeights) %in% rownames(ad)] if (length(toDrop) > 0) { print(paste0('The following ', length(toDrop), ' genes were in the component but not the assay, skipping: ', paste0(toDrop, collapse = ','))) @@ -47,8 +71,6 @@ ScoreUsingSavedComponent <- function(seuratObj, componentOrName, fieldName, assa if (length(names(seuratObj@reductions)) > 0) { suppressMessages(print(FeaturePlot(seuratObj, features = fieldName, order = T) & ggplot2::scale_colour_gradientn(colours = c("navy", "dodgerblue", "gold", "red")))) - } else { - print('No reductions present, cannot plot') } graphics::hist(seuratObj@meta.data[[fieldName]], breaks = 300, main = fieldName) diff --git a/tests/testthat/test-classification.R b/tests/testthat/test-classification.R index 9b3d6af..64944ea 100644 --- a/tests/testthat/test-classification.R +++ b/tests/testthat/test-classification.R @@ -5,7 +5,7 @@ testthat::context("Classification") source('installSeuratData.R') -SimulateSeuratData <- function() { +SimulateSeuratDataForTcrPrediction <- function() { set.seed(GetSeed()) #union features from all activation component files @@ -52,7 +52,9 @@ SimulateSeuratData <- function() { )) seuratObj <- Seurat::NormalizeData(seuratObj, verbose = FALSE) - seuratObj <- Seurat::ScaleData(seuratObj, features = features, verbose = FALSE) + + # NOTE: only scale a subset of genes, so we can test whether PredictTcrActivation automatically adds others: + seuratObj <- Seurat::ScaleData(seuratObj, features = features[1:10], verbose = FALSE) return(seuratObj) } @@ -141,8 +143,15 @@ test_that("Cell type classification works", { test_that("Predict TCell Activation works ", { - SimulateSeuratObj <- SimulateSeuratData() + SimulateSeuratObj <- SimulateSeuratDataForTcrPrediction() + testthat::expect_equal(2995, nrow(Seurat::GetAssayData(SimulateSeuratObj, assay = 'RNA', layer = 'data'))) + testthat::expect_equal(10, nrow(Seurat::GetAssayData(SimulateSeuratObj, assay = 'RNA', layer = 'scale.data'))) SimulateSeuratObj <- PredictTcellActivation(SimulateSeuratObj) + + # This tests whether features were added to scale.data: + testthat::expect_equal(2995, nrow(Seurat::GetAssayData(SimulateSeuratObj, assay = 'RNA', layer = 'data'))) + testthat::expect_equal(2874, nrow(Seurat::GetAssayData(SimulateSeuratObj, assay = 'RNA', layer = 'scale.data'))) + SimulateSeuratObj <- PredictTcellActivation(SimulateSeuratObj, modelName = 'CD4') SimulateSeuratObj <- PredictTcellActivation(SimulateSeuratObj, modelName = 'CD8') @@ -190,7 +199,7 @@ test_that("Predict TCell Activation works ", { test_that("CombineTcellActivationClasses combines classes and writes versioned columns", { - seuratObj <- SimulateSeuratData() + seuratObj <- SimulateSeuratDataForTcrPrediction() seuratObj <- PredictTcellActivation(seuratObj, combineClasses = FALSE) classMapping <- list(