Number of transcripts in the gff3 and position factor output don't match

[gsarfo-boateng]

Dear @taliaferrojm,

Thank you for LABRAT! Super easy to install and use.

I am working with Drosophila melanogaster genome which has over 33k transcripts in its 3UTR. However there are about ~3k transcripts in the number of position file produced by running calculatepsi mode.

Is there a reason for this discrepancy?

See below my script

script:
maksTFfasta step:

python3 LABRAT-master/LABRAT_dm6annotation.py --mode makeTFfasta --gff LABRAT-master/annot/Drosophila_melanogaster.BDGP6.88.chr.gff3 --genomefasta LABRAT-master/annot/Drosophila_melanogaster.BDGP6.genome.fa.gz --lasttwoexons --librarytype 3pseq

SAlmon quant step:

python3 LABRAT-master/LABRAT_dm6annotation.py --mode runSalmon --txfasta TFseqs.fasta --reads1 $READS1_STR --samplename $SAMPLES_STR --threads $THREADS

psi:
python3 LABRAT-master/LABRAT_dm6annotation1.py --mode calculatepsi --salmondir salmondir/ --sampconds sampconds --conditionA $CONDITION_A --conditionB $CONDITION_B --gff LABRAT-master/annot/Drosophila_melanogaster.BDGP6.88.chr.gff3 --librarytype 3pseq

[taliaferrojm]

Hello,

I think you are running things correctly.

This discrepancy is I believe expected. For LABRAT, what is more important than the number of transcripts is the number of distinct 3' ends. So LABRAT is not necessarily calculating differences in transcript usage. Rather, it is calculating differences in polyadenylation site usage. You may have two distinct transcripts, A and B, that differ in an upstream feature (say, an alternative 5' UTR or an alternative coding exon), but if they share the same 3' end, then they share the same polyadenylation site and LABRAT will therefore treat them as a single entity.

Further, if all transcripts for a gene share the same 3' end, then no alternative polyadenylation is possible, and that gene will be excluded from the file you are talking about and ignored by LABRAT.

Finally, there are a number of filters to exclude transcripts whose existence is supported by only weak evidence, and by default only protein coding transcripts are interrogated. You can see these filters (and modify them if you want) in the genefilters and transcriptfilters functions.

Happy to talk about this more if you have more questions.

[gsarfo-boateng]

Thank you for the prompt reply, @taliaferrojm,

I am just wondering if it's typical for LABRAT to annotate some transcripts as having 3′UTRs in numberofposfactors.txt, even if they are not annotated with 3′UTRs in the provided GFF3 file?

Please see the comparison below from one of my output file:

Total unique LABRAT transcripts in numberofposfactors.txt: 7,208

Total transcripts in GFF3: 30,908

Transcripts with 3′UTR annotations in GFF3: 29,849

=== COMPARISON SUMMARY ===

LABRAT transcripts found in GFF3: 4,870
(Of the 7,208 LABRAT transcripts, 4,870 are present in the GFF3 file)

LABRAT transcripts not found in GFF3: 2,338

LABRAT transcripts with 3′UTR annotations: 4,853
(Out of the 4,870 that are found in GFF3, 4,853 have 3′UTRs)

LABRAT transcripts without 3′UTR annotations: 2,355

[taliaferrojm]

Yes, LABRAT is not paying attention at all to 3'UTR annotations in the gff3, mainly because many (most?) gff3 annotations do not have 3'UTR features explicitly annotated. It is looking at the 3' end of 'transcript' features, without regard to whether or not that 'transcript' feature has a '3'UTR' subfeature or not.

[gsarfo-boateng]

Ok. Thank you for the clarification.

[taliaferrojm]

Of course. Let me know if you have more questions about running LABRAT or interpreting results.