This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How to check RNAseq support for annotated genes?

Hello All, I have a set of annotated genes in gff3 format and corresponding RNA-seq data. What is the recommended approach and are there specific tools and parameters to determine the percentage of genes supported by the RNA-seq data?"

Regards, B

gene rna-seq annotation
• 2,502 views • link
updated by Michael Dondrup Michael Dondrup • written by BioinfoBee • BioinfoBee

2 answers

To check for support of gene models by RNAseq, there are several proven workflows. I recommend using them all:

  • Standard RNAseq analysis based on alignments: Reference-based alignment by Hisat or STAR -> featureCounts/HTseq/RSEM -> check the coverage of the gene models of interest. There is no universally accepted cut-off for the number of reads or coverage per gene.
  • If you are interested in a computationally fast approach that also quantifies transcript isoforms (relevant if more than a single isoform is annotated in your GFF) use either Salmon or Kallisto
  • Use reference-guided transcriptome assembly (e.g. in Trinity) and Blast the transcripts. Check for overlaps of transcripts with gene-models
  • Visually inspect genes of interest for anomalies like shadowing, overlapping genes, duplicated sequences, inhomogeneous coverage, etc.
• 0 views • link
written by Michael Dondrup Michael Dondrup

Michael Thanks. I tried it using STAR to map the transcript to the annotation gtf file, and was planning to use featureCounts/gene Counts afterwards to check the read coverage for each genes. For STAR, I used the following command, but it gave me an empty outputs with no transcripts mapped. It is due to the length of the transcripts (>200bp) being used?.

STAR --genomeDir index_directory/ \ --runThreadN 33 \ --runMode alignReads \ --readFilesCommand zcat \ --readFilesIn ${INPUT}/transcriptome.fasta \ --outFileNamePrefix hybrid_star \ --outSAMtype BAM SortedByCoordinate \ --quantMode GeneCounts

Kindly suggest!,

Regards, B

• 0 views • link
written by BioinfoBee • BioinfoBee

I tried it using STAR to map the transcript to the annotation gtf file

That is not what STAR is for. STAR is for aligning short reads in fastq format to a reference genome. If you have assembled transcripts already (from where?) use e.g. GMAP, BBMap, or maybe exonerate.

• 0 views • link
written by Michael Dondrup Michael Dondrup

Michael Thanks again!, Just curious, if its possible to align transcripts to the genic regions only using the gene annotation gtf/gff file. I couldn't find any option in GMAP, unlike STAR which can use GTF annotation. Kindly suggest!

• 0 views • link
written by BioinfoBee • BioinfoBee

I am not sure you can restrict the regions of alignment in STAR. The GTF file is to provide information on splice sites and possibly read counting. During read counting, only genic/exonic regions will be counted.

• 0 views • link
written by Michael Dondrup Michael Dondrup

You could also do transcriptome assembly and then map the transcripts on the annotation via MAKER that will add an AED score for concordance between evidence (here transcript) and the annotation. AED=0 is perfect concordance, AED of 1 is complete absence of support.

• 0 views • link
written by Juke34 Juke34

Juke34 Thanks for the suggestion. I wonder if bedtool coverage or intersect can be used to check the transcript support or coverage for each gene?

• 0 views • link
written by BioinfoBee • BioinfoBee

Log in to answer this question.

More posts like this