# HG changeset patch # User petr-novak # Date 1656419602 0 # Node ID 9de392f2fc022d969299cdd248f9bc9a24c12bbd # Parent c33d6583e54865cd8feb7bf4a69f1ab45ab024b7 "planemo upload commit d6433b48c9bae079edb06364147f19500501c986" diff -r c33d6583e548 -r 9de392f2fc02 R/ltr_utils.R --- a/R/ltr_utils.R Fri Jun 24 14:19:48 2022 +0000 +++ b/R/ltr_utils.R Tue Jun 28 12:33:22 2022 +0000 @@ -138,17 +138,20 @@ gr <- GRanges(seqnames = sapply(gx, function(x)x$seqnames[1]), IRanges(start = S - offset, end = E + offset)) } -get_ranges_left <- function(gx, offset = OFFSET, offset2 = 300) { +get_ranges_left <- function(gx, offset = OFFSET, offset2 = 10) { + ## offset2 - how many nt cen LTR extend to closes protein domain + ## this is necassary as some detected proteins domains does not have correct bopundaries + ## if LTR retrotransposons insters to other protein domain. S <- sapply(gx, function(x)min(x$start)) - max_offset <- S - sapply(gx, function(x)min(x$upstream_domain)) + max_offset <- S - sapply(gx, function(x)min(x$upstream_domain)) + offset2 offset_adjusted <- ifelse(max_offset < offset, max_offset, offset) gr <- GRanges(seqnames = sapply(gx, function(x)x$seqnames[1]), IRanges(start = S - offset_adjusted, end = S + offset2)) return(gr) } -get_ranges_right <- function(gx, offset = OFFSET, offset2 = 300) { +get_ranges_right <- function(gx, offset = OFFSET, offset2 = 10) { E <- sapply(gx, function(x)max(x$end)) - max_offset <- sapply(gx, function(x)max(x$downstream_domain)) - E + max_offset <- sapply(gx, function(x)max(x$downstream_domain)) - E + offset2 offset_adjusted <- ifelse(max_offset < offset, max_offset, offset) gr <- GRanges(seqnames = sapply(gx, function(x)x$seqnames[1]), IRanges(start = E - offset2, end = E + offset_adjusted)) return(gr) @@ -189,7 +192,7 @@ tg_P <- firstTG(cons) ca_P <- lastCA(cons) e_dist <- bl$length[i] - ca_P - max_dist = 50 # was 25 + max_dist <- 50 # was 25 no_match <- any(tg_P == 0, ca_P == 0) if (!no_match & tg_P < max_dist & @@ -684,6 +687,17 @@ return(Rank) } +dante_filtering <- function(dante_gff, min_similarity=0.4, + min_identity=0.2, Relative_Length=0.6, + min_relat_interuptions=8) { + include <- as.numeric(dante_gff$Similarity) >= min_similarity & + as.numeric(dante_gff$Identity) >= min_identity & + as.numeric(dante_gff$Relat_Length) >= Relative_Length & + as.numeric(dante_gff$Relat_Interruptions) <= min_relat_interuptions + + include[is.na(include)] <- FALSE + return(dante_gff[include]) +} get_te_statistics <- function(gr, RT){ Ranks <- c("D", "DL", "DLT", "DLP", "DLTP") diff -r c33d6583e548 -r 9de392f2fc02 README.md --- a/README.md Fri Jun 24 14:19:48 2022 +0000 +++ b/README.md Tue Jun 28 12:33:22 2022 +0000 @@ -1,15 +1,33 @@ -# dante_ltr +# DANTE_LTR + +Tool for identifying complete LTR retrotransposons based on analysis of protein domains identified with the [DANTE tool](https://github.com/kavonrtep/dante). Both DANTE and DANTE_LTR are available on [Galaxy server](ttps://repeatexplorer-elixir.cerit-sc.cz/). -Tool for identification of complete LTR retrotransposons based on analysis of protein -domains identified by DANTE tool. +## Principle of DANTE _LTR +Complete retrotransposons are identified as clusters of protein domains recognized by the DANTE tool. The domains in the clusters must be assigned to a single retrotransposon lineage by DANTE. In addition, the orientation and order of the protein domains, as well as the distances between them, must conform to the characteristics of elements from REXXdb database [Neumann et al. (2019)](https://mobilednajournal.biomedcentral.com/articles/10.1186/s13100-018-0144-1). +In the next step, the 5' and 3' regions of the putative retrotransposon are examined for the presence of 5' and 3' long terminal repeats. If 5'- and 3'-long terminal repeats are detected, detection of target site duplication (TSD) and primer binding site (PSB) is performed. The detected LTR retrotranspsons are classified into 5 categories: +- Elements with protein domains, 5'LTR, 3'LTR, TSD and PBS - rank **DLTP**. +- Elements with protein domains, 5'LTR, 3'LTR, and PBS (TSD was not found) Rank **DLP** +- Elements with protein domains, 5' LTR, 3'LTR, TSD (PBS was not found) - rank **DTL** +- Elements with protein domains, 5'LTR and 3'LTR (PBS and TDS were not found) - rank **DL** +- Elements as clusters of protein domains with the same classification, no LTRs - rank **D**. + +![dante_ltr_workflow.png](dante_ltr_workflow.png) + ## Installation: ```shell conda create -n dante_ltr -c bioconda -c conda-forge -c petrnovak dante_ltr ``` + +## Input data +One input is a reference sequence in fasta fromat. The second input is an annotation of the reference genome using the tool DANTE in GFF3 format. For better results, use the unfiltered full output of the DANTE pipeline. + + ## Usage +### Detection of complete LTR retrotransposons + ```shell Usage: ./extract_putative_ltr.R COMMAND [OPTIONS] @@ -27,25 +45,63 @@ -c NUMBER, --cpu=NUMBER Number of cpu to use [default 5] + -M NUMBER, --max_missing_domains=NUMBER + Maximum number of missing domains is retrotransposon [default 0] + + -L NUMBER, --min_relative_length=NUMBER + Minimum relative length of protein domain to be considered for retrostransposon detection [default 0.6] -h, --help Show this help message and exit + ``` -## Example +#### Example: + ```shell mkdir -p tmp ./extract_putative_ltr.R -g test_data/sample_DANTE.gff3 -s test_data/sample_genome.fasta -o tmp/ltr_annotation ``` -## Output files - - -### Output of script `extract_putative_ltr.R`: - +#### Files in the output of `extract_putative_ltr.R`: - `prefix.gff3` - annotation of all identified elements +- `prefix_D.fasta` - partial elements with protein **d**omains - `prefix_DL.fasta` - elements with protein **d**omains and **L**TR - `prefix_DLTP.fasta` - elements with **d**omains, **L**TR, **T**SD and **P**BS - `prefix_DLP.fasta` - elements with **d**omains, **L**TR and **P**BS - `prefix_DLT.fasta` - elements with **d**omains, **L**TR, **T**SD - `prefix_statistics.csv` - number of elements in individual categories + + + +### Validation of LTR retrotransposons detected un previous step: + +```shell +./clean_ltr.R --help +Usage: ./clean_ltr.R COMMAND [OPTIONS] + + +Options: + -g GFF3, --gff3=GFF3 + gff3 with LTR Transposable elements + + -s REFERENCE_SEQUENCE, --reference_sequence=REFERENCE_SEQUENCE + reference sequence as fasta + + -o OUTPUT, --output=OUTPUT + output file prefix + + -c NUMBER, --cpu=NUMBER + Number of cpu to use [default 5] + + -h, --help + Show this help message and exit +``` + +This script check for potentially chimeric elements and removes them from GFF3 file. + +#### Example +```shell +./clean_ltr.R -g test_data/sample_DANTE_LTR_annotation.gff3 -s test_data/sample_genome.fasta -o tmp/ltr_annotation_clean +``` + diff -r c33d6583e548 -r 9de392f2fc02 clean_dante_ltr.xml --- a/clean_dante_ltr.xml Fri Jun 24 14:19:48 2022 +0000 +++ b/clean_dante_ltr.xml Tue Jun 28 12:33:22 2022 +0000 @@ -1,4 +1,4 @@ - + r-optparse diff -r c33d6583e548 -r 9de392f2fc02 dante_ltr_search.xml --- a/dante_ltr_search.xml Fri Jun 24 14:19:48 2022 +0000 +++ b/dante_ltr_search.xml Tue Jun 28 12:33:22 2022 +0000 @@ -1,4 +1,4 @@ - + blast r-optparse diff -r c33d6583e548 -r 9de392f2fc02 dante_ltr_workflow.odg Binary file dante_ltr_workflow.odg has changed diff -r c33d6583e548 -r 9de392f2fc02 dante_ltr_workflow.png Binary file dante_ltr_workflow.png has changed diff -r c33d6583e548 -r 9de392f2fc02 extract_putative_ltr.R --- a/extract_putative_ltr.R Fri Jun 24 14:19:48 2022 +0000 +++ b/extract_putative_ltr.R Tue Jun 28 12:33:22 2022 +0000 @@ -17,9 +17,11 @@ help = "Number of cpu to use [default %default]", metavar = "number"), make_option(c("-M", "--max_missing_domains"), type = "integer", default = 0, help = "Maximum number of missing domains is retrotransposon [default %default]", + metavar = "number"), + make_option(c("-L", "--min_relative_length"), type = "numeric", default = 0.6, + help = "Minimum relative length of protein domain to be considered for retrostransposon detection [default %default]", metavar = "number") - ) description <- paste(strwrap("")) @@ -71,6 +73,7 @@ s <- readDNAStringSet("/mnt/raid/454_data/cuscuta/Ceuropea_assembly_v4/0_final_asm_hifiasm+longstitch/asm.bp.p.ctg_scaffolds.short_names.fa") lineage_info <- read.table("/mnt/raid/users/petr/workspace/ltr_finder_test/lineage_domain_order.csv", sep = "\t", header = TRUE, as.is = TRUE) + g <- rtracklayer::import("/mnt/raid/users/petr/workspace/dante_ltr/test_data/sample_DANTE_unfiltered.gff3") g <- rtracklayer::import("/mnt/raid/users/petr/workspace/ltr_finder_test/test_data/DANTE_filtered_part.gff3") s <- readDNAStringSet("/mnt/raid/users/petr/workspace/ltr_finder_test/test_data/Rbp_part.fa") @@ -122,6 +125,8 @@ paste(strsplit(d, " ")[[1]], ":", l, sep = "") }, d = lineage_domain, l = names(lineage_domain))) +# filter g gff3 +g <- dante_filtering(g, Relative_Length = opt$min_relative_length) # default seqlengths(g) <- seqlengths(s)[names(seqlengths(g))] g <- add_coordinates_of_closest_neighbor(g) @@ -156,7 +161,7 @@ ) g$Ndomains_in_cluster <- count_occurences_for_each_element(g$Cluster) g$Parent <- paste0("TE_partial_", g$Cluster) -g$Rank="D" +g$Rank <- "D" # keep only partial TE with more than one domain TE_partial_with_more_than_one_domain <- TE_partial[TE_partial$Ndomains > 1]