Mercurial > repos > petr-novak > dante_ltr

comparison R/ltr_utils.R @ 2:f131886ea194 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload commit 891bfe9acf7349c2b887aff6d7e52a7f4ebf3b3a"
author petr-novak
date Tue, 12 Apr 2022 12:55:32 +0000
parents 7b0bbe7477c4
children 6ae4a341d1f3
comparison
equal deleted inserted replaced
1:c1498f679b50 2:f131886ea194
90 }else { 90 }else {
91 return(max(start(x))) 91 return(max(start(x)))
92 } 92 }
93 } 93 }
94 94
95
96
95 trim2TGAC <- function(bl) { 97 trim2TGAC <- function(bl) {
96 for (i in 1:nrow(bl)) { 98 for (i in 1:nrow(bl)) {
97 tg_L <- firstTG(bl$qseq[i]) 99 cons <- consensusString(c(bl$qseq[i], bl$sseq[i]), ambiguityMap="?")
98 tg_R <- firstTG(bl$sseq[i]) 100 tg_P <- firstTG(cons)
99 ca_L <- lastCA(bl$qseq[i]) 101 ca_P <- lastCA(cons)
100 ca_R <- lastCA(bl$sseq[i]) 102 e_dist <- bl$length[i] - ca_P
101 e_dist <- bl$length[i] - ca_R 103 max_dist = 50 # was 25
102 no_match <- any(tg_L == 0, tg_R == 0, ca_L == 0, ca_R == 0) 104 no_match <- any(tg_P == 0, ca_P == 0)
103 if (!no_match & 105 if (!no_match &
104 tg_L == tg_R & 106 tg_P < max_dist &
105 ca_L == ca_R & 107 e_dist < max_dist) {
106 tg_L < 8 &
107 e_dist < 8) {
108 ## trim alignment 108 ## trim alignment
109 bl[i,] <- trim_blast_table(bl[i,], tg_L, e_dist - 1) 109 bl[i,] <- trim_blast_table(bl[i,], tg_P, e_dist - 1)
110 } 110 }
111
111 } 112 }
112 return(bl) 113 return(bl)
113 } 114 }
114 115
115 trim_blast_table <- function(b, T1, T2) { 116 trim_blast_table <- function(b, T1, T2) {
562 } 563 }
563 te <- append(te, pbs_exact_gr) 564 te <- append(te, pbs_exact_gr)
564 unlink(tmp) 565 unlink(tmp)
565 unlink(tmp_out) 566 unlink(tmp_out)
566 return(te) 567 return(te)
568 }
569
570
571 get_te_rank <- function (gr){
572 DL_id <- gr$ID[gr$type == "transposable_element" &
573 gr$TSD == "not_found" &
574 is.na(gr$trna_id)]
575 DLT_id <- gr$ID[gr$type == "transposable_element" &
576 gr$TSD != "not_found" &
577 is.na(gr$trna_id)]
578 DLTP_id <- gr$ID[gr$type == "transposable_element" &
579 gr$TSD != "not_found" &
580 !is.na(gr$trna_id)]
581 DLP_id <- gr$ID[gr$type == "transposable_element" &
582 gr$TSD == "not_found" &
583 !is.na(gr$trna_id)]
584 Rank <- character(length(gr))
585 ID <- unlist(ifelse(is.na(gr$ID), gr$Parent, gr$ID))
586
587 Rank[ID %in% DL_id] <- "DL"
588 Rank[ID %in% DLT_id] <- "DLT"
589 Rank[ID %in% DLP_id] <- "DLP"
590 Rank[ID %in% DLTP_id] <- "DLTP"
591 return(Rank)
567 } 592 }
568 593
569 get_te_statistics <- function(gr, RT) { 594 get_te_statistics <- function(gr, RT) {
570 DOMAINS_LTR <- gr[gr$type == "transposable_element" & 595 DOMAINS_LTR <- gr[gr$type == "transposable_element" &
571 gr$TSD == "not_found" & 596 gr$TSD == "not_found" &