Mercurial > repos > mvdbeek > r_goseq_1_22_0

comparison get_length_and_gc_content.r @ 10:7f8d888e3355 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
planemo upload for repository https://github.com/ARTbio/tools-artbio/tree/master/tools/goseq_1_22_0 commit fdd0811efc61c31f88ff17096fbe8ee8cfacd766-dirty
author mvdbeek
date Thu, 25 Feb 2016 06:18:47 -0500
parents 15ce6435ab83
children 02e88556ce1d
comparison
equal deleted inserted replaced
9:d1eed994287a 10:7f8d888e3355
15 parser <- OptionParser(usage = "%prog [options] file", option_list=option_list) 15 parser <- OptionParser(usage = "%prog [options] file", option_list=option_list)
16 args = parse_args(parser) 16 args = parse_args(parser)
17 17
18 GTFfile = args$gtf 18 GTFfile = args$gtf
19 FASTAfile = args$fasta 19 FASTAfile = args$fasta
20 output = args$output 20 output_file = args$output
21 21
22 #Load the annotation and reduce it 22 #Load the annotation and reduce it
23 GTF <- import.gff(GTFfile, format="gtf", genome=NA, feature.type="exon") 23 GTF <- import.gff(GTFfile, format="gtf", genome=NA, feature.type="exon")
24 grl <- reduce(split(GTF, elementMetadata(GTF)$gene_id)) 24 grl <- reduce(split(GTF, elementMetadata(GTF)$gene_id))
25 reducedGTF <- unlist(grl, use.names=T) 25 reducedGTF <- unlist(grl, use.names=T)
40 c(width, nGCs/width) 40 c(width, nGCs/width)
41 } 41 }
42 output <- t(sapply(split(reducedGTF, elementMetadata(reducedGTF)$gene_id), calc_GC_length)) 42 output <- t(sapply(split(reducedGTF, elementMetadata(reducedGTF)$gene_id), calc_GC_length))
43 colnames(output) <- c("Length", "GC") 43 colnames(output) <- c("Length", "GC")
44 44
45 write.table(output, file="GC_lengths.tsv", sep="\t") 45 write.table(output, file=output_file, sep="\t")