Merge branch 'fix_liftover' into 'dev'

Fix liftover

Closes #138

See merge request tron/bnt_neoants/splice2neo!132

R/liftover_junc_id.R

@@ -7,12 +7,13 @@
 #' @param junc_df a data.frame with at least one column `junc_id` containing junction IDs
 #' @param chain_file path to a chain file for the UCSC liftOver tool. See also \code{\link[rtracklayer]{liftOver}}
 #' @return a data.frame like the input `junc_df` with the following additional columns:
-#'   - `junc_id_lifted_is_unique` a logical vector indicating if the liftOver was successful and unique (1-to-1 correspondence).
+#'   - `liftover_successful` a logical vector indicating if the liftOver was successful and lifted junction positions build valid genomic intervals and not single positions.
+#'   - `liftover_unique` a logical vector indicating if the liftOver was unique (1-to-1 correspondence).
 #'   - `junc_id_lifted_collapsed` a character vector with the lifted junction IDs.
 #'   Multiple IDs are separated by `|`.
 #'   NA represent junc_ids that could not be lifted.
 #'   - `junc_id_lifted` a character vector with a unique lifted junction IDs.
-#'   Potentially multiple lifted IDs are combined by the minmal start and maximal
+#'   Potentially multiple lifted IDs are combined by the minimal start and maximal
 #'   end coordinate. NA represent junc_ids that could not be lifted.
 #'
 #' @examples
@@ -40,10 +41,13 @@ liftover_junc_id <- function(junc_df, chain_file){
       # add lifted junctions ranges as a list
       junc_id_lifted_lst = as.list(junc_id_lifted_grl),
 
-      # check if liftOver was successful and unique
-      junc_id_lifted_is_unique = purrr::map_lgl(junc_id_lifted_lst, ~ length(.x) == 1),
+      # check that lifted junctions are valid intervals (not single positions)
+      liftover_successful = purrr::map_lgl(junc_id_lifted_lst, ~ length(.x) > 0 & all(BiocGenerics::width(.x) >= 2)),
 
-      # covert back to junc_id and collapse multiple IDs with `|`, and replace empty junc_id with NA
+      # check if liftOver was unique
+      liftover_unique = purrr::map_lgl(junc_id_lifted_lst, ~ length(.x) == 1),
+
+      # covert back to junc_id, collapse multiple IDs with `|`, and replace empty junc_id with NA
       junc_id_lifted_collapsed = junc_id_lifted_lst %>%
         purrr::map(as.character) %>%
         purrr::map_chr(stringr::str_c, collapse = "|") %>%

inst/extdata/test.chain

@@ -0,0 +1,8 @@
+chain 170805908 chr6 170805979 + 0 170805979 chr6 170805968 + 0 170805968 6
+30704795	30	0
+4528475	41	0
+17760859	0	11
+48805465	0	18
+29911311	0	31
+39095003
+

tests/testthat/test-liftover_junc_id.R

@@ -8,7 +8,8 @@ test_that("liftover_junc_id works on toy example data", {
   junc_df_lifted <- liftover_junc_id(junc_df, chain_file)
 
   expect_equal(nrow(junc_df_lifted), nrow(junc_df))
-  expect_true(all(c("junc_id_lifted_is_unique", "junc_id_lifted") %in% names(junc_df_lifted)))
+  expect_true(all(c("liftover_successful", "liftover_unique", "junc_id_lifted") %in% names(junc_df_lifted)))
+
 
 })
 
@@ -27,8 +28,27 @@ test_that("liftover_junc_id works with non-unique mappings", {
   junc_df_lifted <- liftover_junc_id(junc_df, chain_file)
 
   expect_equal(nrow(junc_df_lifted), nrow(junc_df))
-  expect_false(all(junc_df_lifted$junc_id_lifted_is_unique))
-  expect_true(all(junc_df_lifted$junc_id_lifted_is_unique[1:3]))
+  expect_false(all(junc_df_lifted$liftover_unique))
+  expect_false(all(junc_df_lifted$liftover_successful))
+  expect_true(all(junc_df_lifted$liftover_unique[1:3]))
   expect_true(any(is.na(junc_df_lifted$junc_id_lifted)))
 
 })
+
+
+test_that("liftover_junc_id for deleted junction position", {
+
+  chain_file = system.file(package="splice2neo", "extdata", "test.chain")
+  stopifnot(chain_file != "")
+
+  junc_df <- dplyr::tibble(
+    junc_id = "chr6:35233301-35233342"
+  )
+
+  junc_df_lifted <- liftover_junc_id(junc_df, chain_file)
+
+  expect_equal(nrow(junc_df_lifted), nrow(junc_df))
+  expect_false(all(junc_df_lifted$liftover_successful))
+  expect_true(all(junc_df_lifted$liftover_unique))
+
+})