Skip to content

Commit

Permalink
Add new LogNormalizeUsingAlternateAssay method
Browse files Browse the repository at this point in the history
  • Loading branch information
bbimber committed Dec 20, 2023
1 parent 0065bcc commit 3c9b513
Show file tree
Hide file tree
Showing 4 changed files with 84 additions and 1 deletion.
1 change: 1 addition & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -31,6 +31,7 @@ export(GetSeed)
export(GetXYDataFromPlot)
export(HighlightCellsOnSeuratPlot)
export(InspectSeurat)
export(LogNormalizeUsingAlternateAssay)
export(MakeEnrichmentDotPlot)
export(MergeSeuratObjs)
export(NanoString_Housekeeping_Normalization)
Expand Down
36 changes: 36 additions & 0 deletions R/Preprocessing.R
Original file line number Diff line number Diff line change
Expand Up @@ -293,4 +293,40 @@ PerformEmptyDrops <- function(seuratRawData, emptyDropNIters, fdrThreshold=0.001
}

stop(paste0('Unable to find matrix file in: ', dataDir, ' or ', dirWithFeatureMatrix))
}

#' @title LogNormalizeUsingAlternateAssay
#'
#' @param seuratObject The seurat object
#' @param assayToNormalize The name of the assay to normalize
#' @param assayForLibrarySize The name of the assay from which to derive library sizes. This will be added to the library size of assayToNormalize.
#' @param scale.factor A scale factor to be applied in normalization
#' @param maxLibrarySizeRatio This normalization relies on the assumption that the library size of the assay being normalized in negligible relative to the assayForLibrarySize. To verify this holds true, the method will error if librarySize(assayToNormalize)/librarySize(assayForLibrarySize) exceeds this value
#' @export
LogNormalizeUsingAlternateAssay <- function(seuratObj, assayToNormalize, assayForLibrarySize = 'RNA', scale.factor = 1e4, maxLibrarySizeRatio = 0.006) {
toNormalize <- Seurat::GetAssayData(seuratObj, assayToNormalize, slot = 'counts')
assayForLibrarySizeData <- Seurat::GetAssayData(seuratObj, assay = assayForLibrarySize, slot = 'counts')

if (any(colnames(toNormalize) != colnames(assayForLibrarySize))) {
stop(paste0('The assayToNormalize and assayForLibrarySize do not have the same cell names!'))
}

margin <- 2
ncells <- dim(x = toNormalize)[margin]

for (i in seq_len(length.out = ncells)) {
x <- toNormalize[, i]
librarySize <- sum(x) + sum(assayForLibrarySizeData[, i])

if ((sum(x) / librarySize) > maxLibrarySizeRatio) {
stop(paste0('The ratio of library sizes was above maxLibrarySizeRatio for cell: ', colnames(assayForLibrarySizeData)[i], '. was: ', (sum(x) / librarySize), ' (', sum(x), ' / ', librarySize, ')'))
}

xnorm <- log1p(x = x / librarySize * scale.factor)
toNormalize[, i] <- xnorm
}

seuratObj <- Seurat::SetAssayData(seuratObj, assay = assayToNormalize, slot = 'data', new.data = toNormalize)

return(seuratObj)
}
28 changes: 28 additions & 0 deletions man/LogNormalizeUsingAlternateAssay.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

20 changes: 19 additions & 1 deletion tests/testthat/test-seurat.R
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
library(DropletUtils);
library(DropletUtils)

context("scRNAseq")

Expand Down Expand Up @@ -124,4 +124,22 @@ test_that("Serat SCTransform works as expected", {

expect_equal(length(rownames(seuratObjSCT@assays$SCT@scale.data)), length(rownames(seuratObjSCT@assays$SCT@counts)))
expect_equal(ncol(seuratObjSCT), ncol(seuratObj))
})


test_that("LogNormalizeUsingAlternateAssay works as expected", {
seuratObj <- suppressWarnings(Seurat::UpdateSeuratObject(readRDS('../testdata/seuratOutput.rds')))

assayToAdd <- Seurat::GetAssayData(seuratObj, assay = 'RNA', layer = 'counts')
assayToAdd <- assayToAdd[1:10,]
assayToAdd <- assayToAdd / 50
rownames(assayToAdd) <- paste0('Feature', LETTERS[1:10])

seuratObj[['Norm']] <- Seurat::CreateAssayObject(assayToAdd)

seuratObj <- LogNormalizeUsingAlternateAssay(seuratObj, assayToNormalize = 'Norm', assayForLibrarySize = 'RNA')

nd <- Seurat::GetAssayData(seuratObj, assay = 'Norm', layer = 'data')
expect_equal(max(nd[,2]), 0.3087937)
expect_equal(max(nd[,100]), 0.7821186)
})

0 comments on commit 3c9b513

Please sign in to comment.