count non-overlapping gene length from gff file

[tompoes]

Hi

I used htseq-count (using default settings) to count all mapped (single-end) reads in a .bam file.
I would like to convert these raw-read count data to RPKM values. Therefore I need to have the non-overlapping transcript (sum of exons) length from all genes in a gff3 file.

You can find a simple example below:

##gff-version
scaffold_1 test gene 47300 48200 . + . ID=A90;Name=A90;
scaffold_1 test mRNA 47300 48200 . + . ID=A90.1;Parent=A90;Name=A90.1;
scaffold_1 test exon 47300 48200 . + . ID=A90.1.1;Parent=A90.1;
scaffold_1 test CDS 47400 48000 . + 0 ID=CDS:A90.1.1;Parent=A90.1;Name=A90.1;
scaffold_1 test gene 48124 49903 . - . ID=A100;Name=A100;
scaffold_1 test mRNA 48124 49903 . - . ID=A100.1;Parent=A100;Name=A100.1;
scaffold_1 test exon 48124 49903 . - . ID=A100.1.1;Parent=A100.1;
scaffold_1 test three_prime_UTR 48124 48464 . - . ID=three_prime_UTR:A100.1.1;Parent=A100.1;Name=A100.1;
scaffold_1 test CDS 48465 49841 . - 0 ID=CDS:A100.1.1;Parent=A100.1;Name=A100.1;
scaffold_1 test five_prime_UTR 49842 49903 . - . ID=five_prime_UTR:A100.1.1;Parent=A100.1;Name=A100.1;

So non overlapping transcript length of A90 is 825 bp (instead of 901) and A100 is 1704 bp (instead of 1780). Are there any tools available that can calculate such unique transcript length or should a script be written?

Thank you for helping!

[Wallysb01]

If all you’re going for is an accurate FPKM, why not just use cufflinks, which assigns FPKMs based on a most likely explanation of the reads method rather than some non-overlapping business gene length combined with htseq counts?

I love htseq/DESeq2 for DE analysis, but for just getting FPKMs cufflinks works great.