@@ -112,7 +112,7 @@ filterTranscriptsByAnnotation <- function(rowDataCombined, annotationGrangesList
112112 } else if (is.null(NDR )) {
113113 NDR <- 0.5
114114 }
115- filterSet <- (rowDataCombined $ NDR < = NDR | rowDataCombined $ readClassType == " equal:compatible" )
115+ filterSet <- (( ! is.na( rowDataCombined $ NDR ) & rowDataCombined $ NDR < = NDR ) | rowDataCombined $ readClassType == " equal:compatible" )
116116 lowConfidenceTranscripts <- combindRowDataWithRanges(
117117 rowDataCombined [! filterSet ,],
118118 exonRangesCombined [! filterSet ])
@@ -224,7 +224,7 @@ calculateNDROnTranscripts <- function(combinedTranscripts, useTxScore = FALSE){
224224 } else {
225225 combinedTranscripts $ NDR <- calculateNDR(combinedTranscripts $ maxTxScore , equal )
226226 }
227- combinedTranscripts $ NDR [combinedTranscripts $ maxTxScore == - 1 ] <- 1
227+ combinedTranscripts $ NDR [combinedTranscripts $ maxTxScore == - 1 ] <- NA
228228 return (combinedTranscripts )
229229}
230230
@@ -827,15 +827,14 @@ addGeneIdsToReadClassTable <- function(readClassTable, distTable,
827827# ' @description This function train a model for use on other data
828828# ' @param extendedAnnotations A GRangesList object produced from bambu(quant = FALSE) or rowRanges(se)
829829# ' @param NDR The maximum NDR for novel transcripts to be in extendedAnnotations (0-1). If not provided a recommended NDR is calculated.
830- # ' @param includeRef A boolean which if TRUE will also filter out reference annotations based on their NDR
831830# ' @param prefix A string which determines which transcripts are considered novel by bambu and will be filtered (by default = 'Bambu')
832831# ' @param baselineFDR a value between 0-1. Bambu uses this FDR on the trained model to recommend an equivilent NDR threshold to be used for the sample. By default, a baseline FDR of 0.1 is used. This does not impact the analysis if an NDR is set.
833832# ' @param defaultModels a bambu trained model object that bambu will use when fitReadClassModel==FALSE or the data is not suitable for training, defaults to the pretrained model in the bambu package
834833# ' Output - returns a similiar GRangesList object with entries swapped into or out of metadata(extendedAnnotations)$lowConfidenceTranscripts
835834# ' @details
836835# ' @return extendedAnnotations with a new NDR threshold
837836# ' @export
838- setNDR <- function (extendedAnnotations , NDR = NULL , includeRef = FALSE , prefix = ' Bambu' , baselineFDR = 0.1 , defaultModels2 = defaultModels ){
837+ setNDR <- function (extendedAnnotations , NDR = NULL , prefix = ' Bambu' , baselineFDR = 0.1 , defaultModels2 = defaultModels ){
839838 # Check to see if the annotations/gtf are dervived from Bambu
840839 if (is.null(mcols(extendedAnnotations )$ NDR )){
841840 warning(" Annotations were not extended by Bambu (or the wrong prefix was provided). NDR can not be set" )
@@ -852,17 +851,10 @@ setNDR <- function(extendedAnnotations, NDR = NULL, includeRef = FALSE, prefix =
852851 message(" Recommending a novel discovery rate (NDR) of: " , NDR )
853852 }
854853
855- # If reference annotations should be filtered too (note that reference annotations with no read support arn't filtered)
856- if (includeRef ){
857- toRemove <- (! is.na(mcols(extendedAnnotations )$ NDR ) & mcols(extendedAnnotations )$ NDR > NDR )
858- toAdd <- ! is.na(mcols(metadata(extendedAnnotations )$ lowConfidenceTranscripts )$ NDR ) &
859- mcols(metadata(extendedAnnotations )$ lowConfidenceTranscripts )$ NDR < = NDR
860- } else {
861- toRemove <- (mcols(extendedAnnotations )$ NDR > NDR &
862- grepl(prefix , mcols(extendedAnnotations )$ TXNAME ))
863- toAdd <- (mcols(metadata(extendedAnnotations )$ lowConfidenceTranscripts )$ NDR < = NDR &
864- grepl(prefix , mcols(metadata(extendedAnnotations )$ lowConfidenceTranscripts )$ TXNAME ))
865- }
854+ toRemove <- (mcols(extendedAnnotations )$ NDR > NDR &
855+ grepl(prefix , mcols(extendedAnnotations )$ TXNAME ))
856+ toAdd <- (mcols(metadata(extendedAnnotations )$ lowConfidenceTranscripts )$ NDR < = NDR &
857+ grepl(prefix , mcols(metadata(extendedAnnotations )$ lowConfidenceTranscripts )$ TXNAME ))
866858
867859 temp <- c(metadata(extendedAnnotations )$ lowConfidenceTranscripts [! toAdd ], extendedAnnotations [toRemove ])
868860 extendedAnnotations <- c(extendedAnnotations [! toRemove ], metadata(extendedAnnotations )$ lowConfidenceTranscripts [toAdd ])
@@ -880,7 +872,7 @@ setNDR <- function(extendedAnnotations, NDR = NULL, includeRef = FALSE, prefix =
880872
881873# ' Extend annotations by clusters
882874# ' @noRd
883- isore.extendAnnotations.clusters <- function (readClassList , annotations , clusters , NDR , isoreParameters , stranded , bpParameters , fusionMode , verbose = FALSE ){
875+ isore.extendAnnotations.clusters <- function (readClassList , annotations , clusters , NDR , discoveryParameters , stranded , bpParameters , fusionMode , verbose = FALSE ){
884876 message(" --- Start extending annotations for clusters ---" )
885877 # if clustering is a csv, create a list with the barcodes for each cluster
886878 # csv must have two cols with heading barcode, cluster
@@ -912,21 +904,21 @@ isore.extendAnnotations.clusters <- function(readClassList, annotations, cluster
912904 rowData(rcf.filt )$ startSD <- 0
913905 rowData(rcf.filt )$ endSD <- 0
914906 rowData(rcf.filt )$ readCount.posStrand <- 0
915- thresholdIndex <- which(rowData(rcf.filt )$ readCount > = isoreParameters $ min.readCount )
916- model <- trainBambu(rcf.filt , verbose = verbose , min.readCount = isoreParameters $ min.readCount )
907+ thresholdIndex <- which(rowData(rcf.filt )$ readCount > = discoveryParameters $ min.readCount )
908+ model <- trainBambu(rcf.filt , verbose = verbose , min.readCount = discoveryParameters $ min.readCount )
917909 txScore <- getTranscriptScore(rowData(rcf.filt )[thresholdIndex ,], model ,
918910 defaultModels )
919911 rowData(rcf.filt )$ txScore <- rep(NA ,nrow(rcf.filt ))
920912 rowData(rcf.filt )$ txScore [thresholdIndex ] <- txScore
921913 # txScores = cbind(txScores, rowData(rcf.filt)$txScore)
922914 rcfs.clusters [[names(clusters )[i ]]] <- rcf.filt
923915 annotations.clusters [[names(clusters )[i ]]] <- bambu.extendAnnotations(list (rcf.filt ), annotations , NDR ,
924- isoreParameters , stranded , bpParameters , fusionMode , verbose )
916+ discoveryParameters , stranded , bpParameters , fusionMode , verbose )
925917 }
926918 if (length(rcfs.clusters )> 0 ){
927919 print(" --- Merging all individual clusters ---" )
928920 annotations.clusters [[" merged" ]] <- bambu.extendAnnotations(rcfs.clusters , annotations , NDR ,
929- isoreParameters , stranded , bpParameters , fusionMode , verbose )
921+ discoveryParameters , stranded , bpParameters , fusionMode , verbose )
930922 }
931923
932924 return (annotations.clusters )
0 commit comments