drop genomic-unstranded option

rnabioco · Aug 28, 2023 · 5883d51 · 5883d51
1 parent 5442866
commit 5883d51
Show file tree

Hide file tree

Showing 11 changed files with 39 additions and 50 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: raer
 Type: Package
 Title: RNA editing tools in R
-Version: 0.99.8
+Version: 0.99.9
 Authors@R: c(
     person("Kent", "Riemondy", , "[email protected]", role = c("aut", "cre"),
            comment = c(ORCID = "0000-0003-0750-1273")),

diff --git a/NEWS.md b/NEWS.md
@@ -1,3 +1,7 @@
+# raer 0.99.9
+
+* The `genomic-unstranded` option for the `library-type` argument in `FilterParam()` has been renamed to `unstranded`, and the `unstranded` option has been removed. 
+
 # raer 0.99.8
 
 * Function arguments involving a fasta file have been renamed to all be `fasta`

diff --git a/R/calc_AEI.R b/R/calc_AEI.R
@@ -19,8 +19,8 @@
 #'   calculating the AEI, typically ALU repeats.
 #' @param txdb A [TxDb] object, if supplied, will be used to subset the alu_ranges
 #'   to those found overlapping genes. Alternatively a [GRanges] object with gene
-#'   coordinates.  If the `library_type` specified by
-#'   `FilterParam`, is `unstranded` or `genomic-unstranded`, then the [TxDb] will
+#'   coordinates.  If the `library_type`, specified by
+#'   `FilterParam`, is `unstranded` then the [TxDb] will
 #'   be used to correct the strandness relative to the reference and is a required
 #'   parameter.
 #' @param snp_db either a [SNPlocs], [GPos], or [GRanges] object. If supplied,
@@ -92,11 +92,11 @@ calc_AEI <- function(bamfiles,
 
     genes_gr <- NULL
     alu_bed_fn <- NULL
-    if (param@library_type %in% c(0, 3)) {
+    if (param@library_type == 0) {
         if (is.null(txdb)) {
             cli::cli_abort("txdb required for processing unstranded data")
         }
-        param@library_type <- 0L
+
         if (is(txdb, "TxDb")) {
             genes_gr <- suppressWarnings(GenomicFeatures::genes(txdb))
         } else {
@@ -364,14 +364,14 @@ get_overlapping_snps <- function(gr,
 #'
 #' @description Gene annotations are used to infer the likely strand of editing
 #'   sites. This function will operate on unstranded datasets which have been
-#'   processed using "genomic-unstranded" library type which reports variants
+#'   processed using "unstranded" library type which reports variants
 #'   with respect to the + strand for all sites. The strand of the editing site
 #'   will be assigned the strand of overlapping features in the `genes_gr`
 #'   object. Sites with no-overlap, or overlapping features with conflicting
 #'   strands (+ and -) will be removed.
 #'
 #' @param rse RangedSummarizedExperiment object containing editing sites processed with
-#'   "genomic-unstranded" setting
+#'   "unstranded" setting
 #' @param genes_gr GRanges object containing reference features to annotate the
 #'   strand of the editing sites.
 #'
@@ -383,7 +383,7 @@ get_overlapping_snps <- function(gr,
 #'
 #' bamfn <- raer_example("SRR5564269_Aligned.sortedByCoord.out.md.bam")
 #' fafn <- raer_example("human.fasta")
-#' fp <- FilterParam(library_type = "genomic-unstranded")
+#' fp <- FilterParam(library_type = "unstranded")
 #' rse <- pileup_sites(bamfn, fafn, param = fp)
 #'
 #' genes <- GRanges(c(

diff --git a/R/pileup.R b/R/pileup.R
@@ -628,21 +628,18 @@ setMethod(show, "FilterParam", function(object) {
 })
 
 
-.encode_libtype <- function(library_type = c("genomic-unstranded",
+.encode_libtype <- function(library_type = c("unstranded",
                                              "fr-first-strand",
-                                             "fr-second-strand",
-                                             "unstranded"),
+                                             "fr-second-strand"),
                             n_files) {
     # encode libtype as integer
-    # 0 = genomic-unstranded  all reads on + strand
+    # 0 = unstranded  all reads on + strand
     # 1 = fr-first-strand     strand based on R1/antisense, R2/sense
     # 2 = fr-second-strand    strand based on R1/sense, R2/antisense
-    # 3 = unstranded          strand based on alignment
     lib_values <- c(
-        "genomic-unstranded",
+        "unstranded",
         "fr-first-strand",
-        "fr-second-strand",
-        "unstranded"
+        "fr-second-strand"
     )
     lib_code <- match(library_type, lib_values)
     if (any(is.na(lib_code))) {
@@ -712,9 +709,8 @@ setMethod(show, "FilterParam", function(object) {
 #' @param min_mapq minimum required MAPQ score. Values for each input BAM file
 #' can be provided as a vector.
 #' @param library_type read orientation, one of `fr-first-strand`,
-#' `fr-second-strand`, `unstranded`, and `genomic-unstranded`. Unstranded library
-#' type will be reported based on read alignment. genomic-unstranded will report
-#' all variants w.r.t the + strand. Values for each
+#' `fr-second-strand`, and `unstranded`. Unstranded library
+#' type will be reported with variants w.r.t the + strand. Values for each
 #' input BAM file can be provided as a vector.
 #' @param only_keep_variants if TRUE, then only variant sites will be reported
 #' (FALSE by default). Values for each input BAM file can be provided as a vector.

diff --git a/man/calc_AEI.Rd b/man/calc_AEI.Rd
diff --git a/man/correct_strand.Rd b/man/correct_strand.Rd
diff --git a/man/pileup_sites.Rd b/man/pileup_sites.Rd
diff --git a/src/plp_utils.c b/src/plp_utils.c
@@ -365,10 +365,6 @@ int invert_read_orientation(bam1_t* b, int libtype) {
         invert = 1;
       }
     }
-  } else if (libtype == 3) {
-    if (is_neg) {
-      invert = 1;
-    }
   } else {
     invert = -1;
   }

diff --git a/tests/testthat/test_aei.R b/tests/testthat/test_aei.R
@@ -56,7 +56,7 @@ test_that("calc_aei basic options work", {
     ssr3_aei <- pc[pc$chrom == "SPCS3" & startsWith(pc$allele, "A_"), ]
     expect_true(all(ssr3_aei[c("alt", "ref")] == 0))
 
-    fp <- FilterParam(library_type = "genomic-unstranded")
+    fp <- FilterParam(library_type = "unstranded")
     expect_error(calc_AEI(bams, fafn, mock_alu_ranges, param = fp))
 
     us_aei <- calc_AEI(bams, fafn, mock_alu_ranges, mock_genes, param = fp)

diff --git a/tests/testthat/test_pileup_sites.R b/tests/testthat/test_pileup_sites.R
@@ -55,7 +55,7 @@ test_that("n-bam pileup works", {
     res <- pileup_sites(c(bamfn, bamfn), fafn, sites,
         param = FilterParam(library_type = c(
             "fr-first-strand",
-            "genomic-unstranded"
+            "unstranded"
         ))
     )
     expect_equal(dim(res), c(214, 2))
@@ -137,15 +137,10 @@ test_that("library types are respected", {
     expect_true(all(strand(res[seqnames(res) == "DHFR"]) == "+"))
     expect_true(all(strand(res[seqnames(res) == "SPCS3"]) == "-"))
 
-    res <- pileup_sites(bamfn, fafn,
-        param = FilterParam(library_type = "genomic-unstranded")
-    )
-    expect_true(all(strand(res) == "+"))
-
     res <- pileup_sites(bamfn, fafn,
         param = FilterParam(library_type = "unstranded")
     )
-    expect_true(all(strand(res) %in% c("+", "-")))
+    expect_true(all(strand(res) == "+"))
 
     expect_error(pileup_sites(bamfn, fafn, sites,
         param = FilterParam(library_type = "unknown-string")
@@ -194,9 +189,8 @@ check_nRef_calc <- function(input, nts_in = nts) {
 
 test_that("pileup check nRef and nAlt", {
     strds <- c(
-        "fr-first-strand", "fr-second-strand",
-        "unstranded", "genomic-unstranded"
-    )
+        "fr-first-strand", "fr-second-strand","unstranded"
+        )
 
     for (strd in strds) {
         res <- pileup_sites(bamfn, fafn,
@@ -476,7 +470,7 @@ test_that("poor quality reads get excluded with read_bqual", {
     res_no_filter <- pileup_sites(bamfn, fafn,
         param = FilterParam(
             min_base_quality = 0L,
-            library_type = c("genomic-unstranded")
+            library_type = c("unstranded")
         )
     )
 
@@ -496,7 +490,7 @@ test_that("poor quality reads get excluded with read_bqual", {
         param = FilterParam(
             read_bqual = bqual_cutoff,
             min_base_quality = 0L,
-            library_type = c("genomic-unstranded")
+            library_type = c("unstranded")
         )
     )
 
@@ -517,7 +511,7 @@ test_that("poor quality reads get excluded with read_bqual", {
 test_that("sites within short splice overhangs can be excluded", {
     # 5 spliced reads, 2 unspliced
     fp <- FilterParam(
-        library_type = "genomic-unstranded",
+        library_type = "unstranded",
         min_splice_overhang = 0,
         min_base_quality = 0
     )
@@ -527,7 +521,7 @@ test_that("sites within short splice overhangs can be excluded", {
 
     # excludes 1 read 5' of splice site (SRR5564269.2688337, cigar 51M120N100M)
     fp <- FilterParam(
-        library_type = "genomic-unstranded",
+        library_type = "unstranded",
         min_splice_overhang = 52,
         min_base_quality = 0
     )
@@ -536,7 +530,7 @@ test_that("sites within short splice overhangs can be excluded", {
     expect_equal(as.vector(assay(res, "nRef")), 6)
 
     fp <- FilterParam(
-        library_type = "genomic-unstranded",
+        library_type = "unstranded",
         min_splice_overhang = 0,
         min_base_quality = 0
     )
@@ -546,7 +540,7 @@ test_that("sites within short splice overhangs can be excluded", {
 
     # excludes 1 read 3' of splice site (SRR5564269.224850, cigar 73M120N78M)
     fp <- FilterParam(
-        library_type = "genomic-unstranded",
+        library_type = "unstranded",
         min_splice_overhang = 79,
         min_base_quality = 0
     )
@@ -556,7 +550,7 @@ test_that("sites within short splice overhangs can be excluded", {
 
     # excludes all 5 spliced reads
     fp <- FilterParam(
-        library_type = "genomic-unstranded",
+        library_type = "unstranded",
         min_splice_overhang = 101,
         min_base_quality = 0
     )

diff --git a/vignettes/novel-sites.Rmd b/vignettes/novel-sites.Rmd
@@ -74,7 +74,7 @@ fp <- FilterParam(
     min_depth = 1,
     min_base_quality = 30,
     min_mapq = c(255, 30),
-    library_type = c("fr-first-strand", "genomic-unstranded"),
+    library_type = c("fr-first-strand", "unstranded"),
     trim_5p = 5,
     trim_3p = 5,
     indel_dist = 4,