Multiplex_Major_Patch

.github/workflows/check-bioc.yml

@@ -40,7 +40,7 @@ env:
   run_covr: 'false'
   run_pkgdown: 'false'
   has_RUnit: 'false'
-  cache-version: 'cache-v3'
+  cache-version: 'cache-v4'
   run_docker: 'false'
 
 jobs:
@@ -54,9 +54,9 @@ jobs:
       fail-fast: false
       matrix:
         config:
-          - { os: ubuntu-latest, r: '4.3', bioc: '3.18', cont: "bioconductor/bioconductor_docker:RELEASE_3_18", rspm: "https://packagemanager.rstudio.com/cran/__linux__/focal/latest" }
-          - { os: macOS-latest, r: '4.3', bioc: '3.18'}
-            ##- { os: windows-latest, r: '4.3', bioc: '3.18'}
+          - { os: ubuntu-latest, r: '4.4', bioc: '3.19', cont: "bioconductor/bioconductor_docker:RELEASE_3_19", rspm: "https://packagemanager.rstudio.com/cran/__linux__/focal/latest" }
+          - { os: macOS-latest, r: '4.4', bioc: '3.19'}
+          - { os: windows-latest, r: '4.4', bioc: '3.19'}
           ## Check https://github.com/r-lib/actions/tree/master/examples
           ## for examples using the http-user-agent
     env:
@@ -81,7 +81,7 @@ jobs:
       ## https://github.com/r-lib/actions/blob/master/examples/check-standard.yaml
       ## If they update their steps, we will also need to update ours.
       - name: Checkout Repository
-        uses: actions/checkout@v3
+        uses: actions/checkout@v4
 
       ## R is already included in the Bioconductor docker images
       - name: Setup R from r-lib
@@ -104,15 +104,15 @@ jobs:
 
       - name: Restore R package cache
         if: "!contains(github.event.head_commit.message, '/nocache') && runner.os != 'Linux'"
-        uses: actions/cache@v3
+        uses: actions/cache@v4
         with:
           path: ${{ env.R_LIBS_USER }}
           key: ${{ env.cache-version }}-${{ runner.os }}-biocversion-RELEASE-r-4.3-${{ hashFiles('.github/depends.Rds') }}
           restore-keys: ${{ env.cache-version }}-${{ runner.os }}-biocversion-RELEASE-r-4.3-
 
       - name: Cache R packages on Linux
         if: "!contains(github.event.head_commit.message, '/nocache') && runner.os == 'Linux' "
-        uses: actions/cache@v3
+        uses: actions/cache@v4
         with:
           path: /home/runner/work/_temp/Library
           key: ${{ env.cache-version }}-${{ runner.os }}-biocversion-devel-r-4.3-${{ hashFiles('.github/depends.Rds') }}
@@ -337,9 +337,9 @@ jobs:
 
       - name: Upload check results
         if: failure()
-        uses: actions/upload-artifact@v2
+        uses: actions/upload-artifact@v4
         with:
-          name: ${{ runner.os }}-biocversion-RELEASE-r-4.3-results
+          name: ${{ runner.os }}-biocversion-RELEASE-r-4.4-results
           path: check
 
       - uses: docker/build-push-action@v1

.github/workflows/lint.yaml

@@ -15,7 +15,7 @@ jobs:
     env:
       GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
     steps:
-      - uses: actions/checkout@v2
+      - uses: actions/checkout@v4
 
       - uses: r-lib/actions/setup-r@v2
         with:

.github/workflows/pr-commands.yaml

@@ -14,7 +14,7 @@ jobs:
     env:
       GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
     steps:
-      - uses: actions/checkout@v2
+      - uses: actions/checkout@v4
 
       - uses: r-lib/actions/pr-fetch@v2
         with:
@@ -49,7 +49,7 @@ jobs:
     env:
       GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
     steps:
-      - uses: actions/checkout@v2
+      - uses: actions/checkout@v4
 
       - uses: r-lib/actions/pr-fetch@v2
         with:

.github/workflows/test-coverage.yaml

@@ -16,7 +16,7 @@ jobs:
       GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
 
     steps:
-      - uses: actions/checkout@v3
+      - uses: actions/checkout@v4
 
       - uses: r-lib/actions/setup-r@v2
         with:

DESCRIPTION

@@ -86,7 +86,8 @@ Imports:
     Rsamtools,
     methods,
     Rcpp,
-    xgboost
+    xgboost,
+    Matrix
 VignetteBuilder: 
     knitr
 LazyData: true

NAMESPACE

@@ -7,7 +7,9 @@ export(readFromGTF)
 export(transcriptToGeneExpression)
 export(writeBambuOutput)
 export(writeToGTF)
+export(writeAnnotationsToGTF)
 export(trainBambu)
+export(setNDR)
 export(compareTranscripts)
 importFrom(stats,predict)
 importFrom(BiocGenerics,basename)
@@ -76,4 +78,5 @@ import(SummarizedExperiment)
 import(S4Vectors, except=c(rename, setequal, setdiff, intersect,cor))
 useDynLib(bambu, .registration = TRUE)
 import(xgboost)
+import(Matrix)
 import(BSgenome)

R/bambu-assignDist.R

@@ -0,0 +1,93 @@
+#' Create equivilence classes and assign to transcripts
+#' @inheritParams bambu
+#' @import data.table
+#' @noRd
+assignReadClasstoTranscripts <- function(readClassList, annotations, isoreParameters, 
+                                        verbose, demultiplexed, spatial, 
+                                        returnDistTable = FALSE, trackReads = TRUE) {
+    metadata(readClassList)$readClassDist <- calculateDistTable(readClassList, annotations, isoreParameters, verbose)
+    readClassList = splitReadClassFiles(readClassList)
+    readClassDt <- genEquiRCs(metadata(readClassList)$readClassDist, annotations, verbose) 
+    readClassDt$eqClass.match = match(readClassDt$eqClassById,metadata(readClassList)$eqClassById)
+    readClassDt <- simplifyNames(readClassDt)
+    readClassDt = readClassDt %>% group_by(eqClassId, gene_sid) %>% 
+        mutate(multi_align = length(unique(txid))>1) %>% ungroup() %>% mutate(aval = 1) %>%
+        data.table()
+
+    #return non-em counts
+    ColData = generateColData(colnames(metadata(readClassList)$countMatrix), clusters = NULL, demultiplexed, spatial)
+    quantData <- SummarizedExperiment(assays = SimpleList(
+        counts = generateUniqueCounts(readClassDt, metadata(readClassList)$countMatrix, annotations)),
+        rowRanges = annotations,
+        colData = ColData)
+    colnames(quantData) = ColData$id
+    if(sum(metadata(readClassList)$incompatibleCountMatrix)==0){
+        metadata(quantData)$incompatibleCounts = NULL
+    } else{
+        metadata(quantData)$incompatibleCounts = generateIncompatibleCounts(metadata(readClassList)$incompatibleCountMatrix, annotations)       
+    }
+    metadata(quantData)$nonuniqueCounts = generateNonUniqueCounts(readClassDt, metadata(readClassList)$countMatrix, annotations)
+    metadata(quantData)$readClassDt = readClassDt
+    metadata(quantData)$countMatrix = metadata(readClassList)$countMatrix
+    metadata(quantData)$incompatibleCountMatrix  = metadata(readClassList)$incompatibleCountMatrix 
+    metadata(quantData)$sampleNames = metadata(readClassList)$sampleNames 
+    if(returnDistTable){
+        metadata(quantData)$distTable = metadata(readClassList)$readClassDist
+    }
+    if (trackReads){
+        metadata(quantData)$readToTranscriptMap = 
+            generateReadToTranscriptMap(readClassList, metadata(readClassList)$readClassDist, 
+                                  annotations)
+    } 
+    return(quantData)     
+
+}
+
+
+generateUniqueCounts <- function(readClassDt, countMatrix, annotations){
+    x = readClassDt %>% filter(!multi_align & !is.na(eqClass.match))
+    uniqueCounts = countMatrix[x$eqClass.match,]
+    uniqueCounts.tx = sparse.model.matrix(~ factor(x$txid) - 1)
+    uniqueCounts = t(uniqueCounts.tx) %*% uniqueCounts
+    rownames(uniqueCounts) = names(annotations)[match(as.numeric(levels(factor(x$txid))),mcols(annotations)$txid)]
+    counts = sparseMatrix(length(annotations), ncol(uniqueCounts), x = 0)
+    rownames(counts) = names(annotations)
+    counts[rownames(uniqueCounts),] = uniqueCounts
+    return(counts)
+
+    counts.total = colSums(countMatrix) + colSums(incompatibleCountMatrix)
+    counts.total[counts.total==0] = 1
+    counts.CPM = counts/counts.total * 10^6
+
+}
+
+generateIncompatibleCounts <- function(incompatibleCountMatrix, annotations){
+    genes = levels(factor(unique(mcols(annotations)$GENEID)))
+    rownames(incompatibleCountMatrix) = genes[as.numeric(rownames(incompatibleCountMatrix))]
+    geneMat = sparseMatrix(length(genes), ncol(incompatibleCountMatrix), x = 0)
+    rownames(geneMat) = genes
+    geneMat[rownames(incompatibleCountMatrix),] = incompatibleCountMatrix
+    return(geneMat)
+}
+
+generateNonUniqueCounts <- function(readClassDt, countMatrix, annotations){
+    #fuse multi align RCs by gene
+    x = readClassDt %>% filter(multi_align & !is.na(eqClass.match))
+    x = x %>% distinct(eqClassId, .keep_all = TRUE)
+    nonuniqueCounts = countMatrix[x$eqClass.match,, drop = FALSE]
+    if(nrow(x)>1){
+    nonuniqueCounts.gene = sparse.model.matrix(~ factor(x$gene_sid) - 1)
+    nonuniqueCounts = t(nonuniqueCounts.gene) %*% nonuniqueCounts
+    }
+    #covert ids into gene ids
+    geneids = as.numeric(levels(factor(x$gene_sid)))
+    geneids = x$txid[match(geneids, x$gene_sid)]
+    geneids = mcols(annotations)$GENEID[as.numeric(geneids)]
+    rownames(nonuniqueCounts) = geneids
+    #create matrix for all annotated genes
+    genes = levels(factor(unique(mcols(annotations)$GENEID)))
+    geneMat = sparseMatrix(length(genes), ncol(nonuniqueCounts), x = 0)
+    rownames(geneMat) = genes
+    geneMat[rownames(nonuniqueCounts),] = nonuniqueCounts
+    return(geneMat)
+}

R/bambu-extendAnnotations-utilityCombine.R

@@ -19,6 +19,7 @@ isore.combineTranscriptCandidates <- function(readClassList,
         min.readCount, min.readFractionByGene, 
         min.txScore.multiExon, min.txScore.singleExon, verbose) %>% data.table()
     combinedSplicedTranscripts[,confidenceType := "highConfidenceJunctionReads"]
+    if (min.txScore.singleExon < 1) {return(combinedSplicedTranscripts)}
     combinedUnsplicedTranscripts <- 
         combineUnsplicedTranscriptModels(readClassList, bpParameters, 
         stranded, min.readCount, min.readFractionByGene,