Introduces setRcAssignmentParameters() with min.exonDistance, min.primarySecondaryDist, and min.primarySecondaryDistStartEnd2, so the assignment step can be tuned independently from opt.discovery. calculateDistTable() and assignReadClasstoTranscripts() now use rcAssignmentParameters instead of isoreParameters.

Author: SuiYue-2308

R/bambu-assignDist.R

@@ -2,11 +2,11 @@
 #' @inheritParams bambu
 #' @import data.table
 #' @noRd
-assignReadClasstoTranscripts <- function(readClassList, annotations, isoreParameters, 
-                                        verbose, sampleMetadata, demultiplexed,
+assignReadClasstoTranscripts <- function(readClassList, annotations,
+                                        rcAssignmentParameters, verbose, sampleMetadata, demultiplexed,
                                         returnDistTable = FALSE, trackReads = TRUE) {
     if (is.character(readClassList)) readClassList <- readRDS(file = readClassList)
-    metadata(readClassList)$readClassDist <- calculateDistTable(readClassList, annotations, isoreParameters, verbose, returnDistTable)
+    metadata(readClassList)$readClassDist <- calculateDistTable(readClassList, annotations, rcAssignmentParameters, verbose, returnDistTable)
     readClassList <- splitReadClassFiles(readClassList)
     readClassDt <- genEquiRCs(metadata(readClassList)$readClassDist, annotations, verbose) 
     readClassDt$eqClass.match = match(readClassDt$eqClassById,metadata(readClassList)$eqClassById)

R/bambu.R

@@ -62,6 +62,17 @@
 #'     \item{prefix}{specifying prefix for new gene Ids (genePrefix.number),
 #'     defaults to "Bambu"}
 #' }
+#' @param opt.rcAssignment A list of controlling parameters for the read class
+#' to transcript assignment process:
+#' \describe{
+#'     \item{min.exonDistance}{specifying minimum distance to known transcript
+#'     to be considered a valid match, defaults to 35bp}
+#'     \item{min.primarySecondaryDist}{specifying the minimum distance
+#'     threshold between primary and secondary assignments, defaults to 5bp}
+#'     \item{min.primarySecondaryDistStartEnd2}{specifying the minimum
+#'     distance threshold for start/end positions used for read assignment,
+#'     defaults to 5bp}
+#' }
 #' @param opt.em A list of controlling parameters for quantification
 #' algorithm estimation process:
 #' \describe{
@@ -141,8 +152,8 @@
 #'     genome = fa.file,  discovery = TRUE, quant = TRUE)
 #' @export
 bambu <- function(reads, annotations = NULL, genome = NULL, NDR = NULL,
-    mode = NULL, opt.discovery = NULL, opt.em = NULL, rcOutDir = NULL, discovery = TRUE, 
-    assignDist = TRUE, quant = TRUE, stranded = FALSE,  ncore = 1, yieldSize = NULL,  
+    mode = NULL, opt.discovery = NULL, opt.rcAssignment = NULL, opt.em = NULL, rcOutDir = NULL, discovery = TRUE,
+    assignDist = TRUE, quant = TRUE, stranded = FALSE,  ncore = 1, yieldSize = NULL,
     trackReads = FALSE, returnDistTable = FALSE, lowMemory = FALSE, sampleData = NULL,
     fusionMode = FALSE, verbose = FALSE, demultiplexed = FALSE, quantData = NULL,
     sampleNames = NULL, cleanReads = FALSE, dedupUMI = FALSE, barcodesToFilter = NULL, clusters = NULL,
@@ -183,7 +194,7 @@ bambu <- function(reads, annotations = NULL, genome = NULL, NDR = NULL,
     isoreParameters <- setIsoreParameters(isoreParameters = opt.discovery)
     #below line is to be compatible with earlier version of running bambu
     if(!is.null(isoreParameters$max.txNDR)) NDR = isoreParameters$max.txNDR
-    
+    rcAssignmentParameters <- setRcAssignmentParameters(rcAssignmentParameters = opt.rcAssignment)
     emParameters <- setEmParameters(emParameters = opt.em)
     bpParameters <- setBiocParallelParameters(reads, ncore, verbose, demultiplexed)
 	xgb.set.config(nthread = 1)
@@ -242,9 +253,9 @@ bambu <- function(reads, annotations = NULL, genome = NULL, NDR = NULL,
             quantData <- bplapply(seq_along(readClassList), function(i){
               assignReadClasstoTranscripts(
                 readClassList = readClassList[[i]],
-                annotations = annotations, 
-                isoreParameters = isoreParameters, 
-                verbose = verbose, 
+                annotations = annotations,
+                rcAssignmentParameters = rcAssignmentParameters,
+                verbose = verbose,
                 # for bulk data, there is one sampleData (keep sampleData[1]), for single-cell, there is one per sample
                 sampleMetadata = if(length(sampleData) == 1) sampleData[1] else sampleData[i],
                 demultiplexed = demultiplexed, 

R/bambu_utilityFunctions.R

@@ -41,6 +41,19 @@ setIsoreParameters <- function(isoreParameters){
 }
 
 
+#' setRcAssignmentParameters
+#' @noRd
+setRcAssignmentParameters <- function(rcAssignmentParameters){
+    rcAssignmentParameters.default <- list(
+        min.exonDistance = 35,
+        min.primarySecondaryDist = 5,
+        min.primarySecondaryDistStartEnd2 = 5)
+    rcAssignmentParameters <-
+        updateParameters(rcAssignmentParameters, rcAssignmentParameters.default)
+    return(rcAssignmentParameters)
+}
+
+
 #' setEmParameters
 #' @noRd
 setEmParameters <- function(emParameters){
@@ -242,11 +255,11 @@ handleWarnings <- function(readClassList, verbose){
 }
 
 #' Calculate the dist table used for Bambu Quantification
-calculateDistTable <- function(readClassList, annotations, isoreParameters, verbose, returnDistTable){
+calculateDistTable <- function(readClassList, annotations, rcAssignmentParameters, verbose, returnDistTable){
     readClassDist <- isore.estimateDistanceToAnnotations(readClassList, annotations,
-                                                            min.exonDistance = isoreParameters[["min.exonDistance"]],
-                                                            min.primarySecondaryDist = isoreParameters[['min.primarySecondaryDist']],
-                                                            min.primarySecondaryDistStartEnd = isoreParameters[['min.primarySecondaryDistStartEnd2']],
+                                                            min.exonDistance = rcAssignmentParameters[["min.exonDistance"]],
+                                                            min.primarySecondaryDist = rcAssignmentParameters[['min.primarySecondaryDist']],
+                                                            min.primarySecondaryDistStartEnd = rcAssignmentParameters[['min.primarySecondaryDistStartEnd2']],
                                                             verbose = verbose)
         metadata(readClassDist)$distTable <- modifyIncompatibleAssignment(metadata(readClassDist)$distTable)
         if(returnDistTable) metadata(readClassDist)$distTableOld <- metadata(readClassDist)$distTable