comparison extract_putative_ltr.R @ 8:9de392f2fc02 draft

"planemo upload commit d6433b48c9bae079edb06364147f19500501c986"

7:c33d6583e548

              help = "output file path and prefix", default = NULL),
  make_option(c("-c", "--cpu"), type = "integer", default = 5,
              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")


)
description <- paste(strwrap(""))

epilogue <- ""

if (FALSE) {
  g <- rtracklayer::import("/mnt/raid/454_data/cuscuta/Ceuropea_assembly_v4/0_final_asm_hifiasm+longstitch/repeat_annotation/DANTE_on_CEUR_filtered_short_names.gff3")
  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/ltr_finder_test/test_data/DANTE_filtered_part.gff3")
  s <- readDNAStringSet("/mnt/raid/users/petr/workspace/ltr_finder_test/test_data/Rbp_part.fa")

  g <- rtracklayer::import("/mnt/raid/users/petr/workspace/dante_ltr/test_data
  /DANTE_Vfaba_chr5.gff3")

names(lineage_ltr_mean_length) <- ln
lineage_domains_sequence <- unlist(mapply(function(d,l) {
  paste(strsplit(d, " ")[[1]], ":", l, sep = "")
}, d = lineage_domain, l = names(lineage_domain)))


seqlengths(g) <- seqlengths(s)[names(seqlengths(g))]
g <- add_coordinates_of_closest_neighbor(g)

# add info about domain order:

                       IRanges(start = sapply(gcl_alt, function(x) min(x$start)),
                               end = sapply(gcl_alt, function(x) max(x$end)))
)
g$Ndomains_in_cluster <- count_occurences_for_each_element(g$Cluster)
g$Parent <- paste0("TE_partial_", g$Cluster)
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]
g_with_more_than_one_domain <- g[as.vector(g$Ndomains_in_cluster > 1)]

8:9de392f2fc02

              help = "output file path and prefix", default = NULL),
  make_option(c("-c", "--cpu"), type = "integer", default = 5,
              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(""))

epilogue <- ""

if (FALSE) {
  g <- rtracklayer::import("/mnt/raid/454_data/cuscuta/Ceuropea_assembly_v4/0_final_asm_hifiasm+longstitch/repeat_annotation/DANTE_on_CEUR_filtered_short_names.gff3")
  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")

  g <- rtracklayer::import("/mnt/raid/users/petr/workspace/dante_ltr/test_data
  /DANTE_Vfaba_chr5.gff3")

names(lineage_ltr_mean_length) <- ln
lineage_domains_sequence <- unlist(mapply(function(d,l) {
  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)

# add info about domain order:

                       IRanges(start = sapply(gcl_alt, function(x) min(x$start)),
                               end = sapply(gcl_alt, function(x) max(x$end)))
)
g$Ndomains_in_cluster <- count_occurences_for_each_element(g$Cluster)
g$Parent <- paste0("TE_partial_", g$Cluster)
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]
g_with_more_than_one_domain <- g[as.vector(g$Ndomains_in_cluster > 1)]