I agree with Artur, you can simply map your RNAseq data to the reference genome. Use a splice-aware aligner such as Tophat2 and then load the output BAM file in a genome browser (eg. IGV). This will show you where split reads are mapping, which indicate intron sequences. The coverage over the intronic region indicates intron retention. This doesn't even require an annotation, although loading an annotation into the genome browser would make it much easier to see what's going on with read mapping patterns. Tophat also provides an output file with the positions of all introns, I believe it's in BED format. You can use these positions to extract intron sequences if that is your goal.

There are many other splice-aware aligners you could try, see this reference:

http://www.nature.com/nmeth/journal/v10/n12/full/nmeth.2722.html

Dear Muthusamy,

What do you need to do with your RNA seq data? Do you really need to remove those putative intronic segments from your RNA seq sequences? I ask this because you say that retention by alternative splicing might be happening, and thus you probably need that information in order to properly predict the product of your transcripts. Also, when you use BLAST to compare with your CDS, you can adjust the parameters of BLAST in order to allow longer gaps. That is normally done by adjusting the "gap_penalty" parameter (the lower the gap_penalty, the longer the gaps that are allowed). This might help you to get matches even in presence of intronic sequences.

Best regards.

Maybe he just doesn't trust very much the annotation of the current available genomes, so he is not sure of where the CDS start and end? this would not be a surprise to me -- bad/poor annotation of genomes is something rather common. So he might want to map the translated reads to a protein database instead, rather than mapping them to a genome that he isn't sure about its annotation quality? I don't know, just guessing here, that's why I started asking.. but I can think on a couple of reasons about why I wouldn't just map the reads to a genome sequence.

Best regards.

I do totally agree with you this time, Artur.

I hope you would clarify the above questions. However if you have your introns annotated and you have coordinates for your intronic regions then you could get your RNAseq read counts for each gene in the whole genome using BedTools  (http://bedtools.readthedocs.io/en/latest/index.html). This way you could see how prominant intron retention is in your species under the condition you performed the RNAseq.

Cheers,

Shajahan

I agree with Artur, you can simply map your RNAseq data to the reference genome. Use a splice-aware aligner such as Tophat2 and then load the output BAM file in a genome browser (eg. IGV). This will show you where split reads are mapping, which indicate intron sequences. The coverage over the intronic region indicates intron retention. This doesn't even require an annotation, although loading an annotation into the genome browser would make it much easier to see what's going on with read mapping patterns. Tophat also provides an output file with the positions of all introns, I believe it's in BED format. You can use these positions to extract intron sequences if that is your goal.

There are many other splice-aware aligners you could try, see this reference:

http://www.nature.com/nmeth/journal/v10/n12/full/nmeth.2722.html

 Dear Salvador Ramirez-Flandes  and All

First of all apologies for late response. I was out of station. Thank you so much for thought provoking questions and suggestions.

Sorry for not being clear in my question. I must try simple things as Salvador Ramirez-Flandes said like increasing the gap penalty score...

The work is to prove that drought modulated splicing pattern in banana transcriptome. All sequence reads were assembled and annotated . I am just trying to analyse and compare splicing pattern b/w two banana varieties.

As far as Reference banana genome annotation is concerned, still it needs manual curation. Moreover for my banana cultivar still not 100% of whole genome is available.

What i understand from my study is there are intron sequences present in my transcript reads. Now i want see what kind of introns, sequences features or they following uniform pattern or different across cultivars?, impact of intron retention on gene products etc.

So i would like to have each intron sequences present in transcriptome as prerequisite for my experiment

Hope i am clear in my question this time.

Once again thank you so much for your valuable time and suggestions.

Dear Muthusamy,

How did you assemble your sequence reads? paired end assembly? After this now I think you need to assemble your resulting reads somehow, and this can be done with the help of a reference genome (you've received a couple of interesting suggestions about this here), or with a de novo assembler. The latter alternative might be a good option for you because you are saying that you still don't have access to the 100% of the genome. A popular de novo assembler for RNA-seq data is Trinity (https://trinityrnaseq.github.io/). 

After the assembly (either by a reference genome or de novo, or even maybe both), you will have (hopefully) full-length transcriptomic sequences, that now you can process further, either for predicting putative protein products (alignment against a protein sequence database) or for determining alternative splicing sites, etc. In all these analyses proper adjustment of the parameters for the aligner (blast or any other) is important in my opinion, unless you use an already configured pipeline for all this (I always prefer to do each step by myself though, so I can't recommend any such pipeline if exists).

PS: There is also a recent article that compare de novo transcriptome assembly software if you are interested: Comprehensive evaluation of de novo transcriptome assembly programs and their effects on differential gene expression analysis (2017) (https://www.ncbi.nlm.nih.gov/pubmed/27694201)

Best regards.

As Salvador suggests, if the genome assembly is not very good or complete then you could try a de novo transcriptome assembly. Trinity is probably the standard for transcriptome assembly, but Bridger is also a good choice. 

I have tried doing exactly what is being suggested here and it was not easy or very successful. Theoretically, if there are reads covering the introns then the assembler should be able to assemble multiple isoforms, some with introns present, others with introns removed. Then you can simply select the longest assembled isoform and map reads back to that to look at splicing. It probably depends a lot on the quality and depth of your RNAseq data set, but when I tried this only a very small proportion of introns were found in assembled transcripts despite the fact that I know I had reads that covered those introns (I tested it on a data set where I did have a reference genome). Anyways, it is worth a try and you should probably try several assemblies with different parameters and/or programs to see which gives you the most transcripts and the longest (containing the most introns).

Let me know if you have more questions about this.

In the case of incomplete genomes, where are going to be mapped the reads associated to the missing parts? OK maybe I'm too used to prokaryotic genomes, that are so extremely dynamical, but perhaps here if the reference genome is only a close relative then the unmapped reads will be negligible, right?

Indeed de novo assembling is more complex and computer resources demanding than a "mapper", but servers with 512GB of main memory (a typical requirement to run almost any de novo assembler) today are not so expensive. Moreover, Trinity can be run on publicly available computer resources:

https://github.com/trinityrnaseq/trinityrnaseq/wiki/Accessing-Trinity-on-Publicly-Available-Compute-Resources

(the most interesting is the Galaxy initiative: https://usegalaxy.org -- a lot of other resource-demanding tools publicly available too)

But yes.. I agree, you can probably run 100 mappings in the time of a single de novo assembly process.

Best regards.