Question: Cufflinks Gene Transcripts
gravatar for rbj
8 months ago by
rbj30 wrote:

Hi I have been experimenting with Bowtie and Cufflinks and I now reach a point where it appears to work; but I think there is a subtle error. I don't find ENOUGH gene hits.
So (1) I got the same dataset professionally analysed by GATC and they returned 20k genes, but when I do it I get 2k. (2) When I inspected the tabular output from Cufflinks, all my 2k hits are in the first 4000 genes in the list. The Tabular file displays 33000 rows, and all hits in the first 4000 rows? That's implausible. Putting two and two together it seems like I am only really properly mapping the first bit of the files? Any ideas anyone?

I think it is a cufflinks thing, because even when I used the BAM from the GATC service, I still only got 10% of the hits they returned. I guess (a) Is there a way that cufflinks would trip up mid way through a run and NOT return any errors? (b) Are there, like, some essential options I should have selected in my cufflinks run. I have tried a couple of different options, but the manual to me is not really idiot proof enough! Thanks Everyone Richard

rna-seq genes cufflinks bowtie • 228 views
ADD COMMENTlink modified 8 months ago by Jennifer Hillman Jackson24k • written 8 months ago by rbj30
gravatar for Jennifer Hillman Jackson
8 months ago by
United States
Jennifer Hillman Jackson24k wrote:

Hi Richard,

1 - Are you using a reference annotation dataset with Cufflinks (GTF, GFF3)? Double check that the chromosome identifiers are an exact match to for the reference genome used for mapping.

2 - Also be sure that you have sorted the Bowtie BAM dataset and annotation GTF/GFF3 (if using one) prior to running Cufflinks.

3 - If using a Custom Reference genome, it may need to be reformatted to remove description content. This is somewhat related to item 1 above, as this particular tool suite does not automatically clip description content off of fasta title lines and leading to false "mismatches" for the identifier portion of the title line.

How to check for a reference genome mismatch, how to reformat custom genomes, and how to sort are included in FAQs here:

Other issues can be a factor but these three are the most common and are a good place to start troubleshooting.

Let us know how that works out. If working at http:/ or if you can reproduce the problem there, and cannot determine the problem, a bug report can be sent in for direct feedback. This is how:

Thanks! Jen, Galaxy team

ADD COMMENTlink written 8 months ago by Jennifer Hillman Jackson24k

Thank you very much Jen, I am still diagnosing; firstly (for everyone's info), this is rat tissue so the genome is (release 6) is installed in, but the annotation reference is one I prepared. So that is likely to be the error. I believe it is sorted, and I am pretty sure the bam is. I suspect some other error with this. So this is not a custom genome though. The chromosome label issue was dealt with by me with a script to convert "1" to "chr1", but I wonder if, in doing this I mucked something else up. I presume it is sorted, because it the locations look ascending, by eye. The suspected error is totally reproducible. When I just re-ran cufflinks on my bam WITHOUT annotation reference, I got a more realistic 20k+ hits ("genes"). but in the un-annotated cuff.1, cuff.2 format. So basically, it appears that something must be wrong in the midst's of my annotation file? So I am looking at this, but if anyone has a clue from their own experience of an annotation that partly works, but misses the latter half, please let me know!!! R.

ADD REPLYlink written 8 months ago by rbj30

...and, if anyone is following this saga in real time :-) I have now established that when I converted my "1"s to "chr1" etc.... (in Matlab) I only seem to have done the chromosome 1's not all the others. Correcting this now. I expect that was my fault. I will update when I have fixed my appending program/script.

ADD REPLYlink written 8 months ago by rbj30

So to conclude with this particular issue I've had. The error was indeed in the rat.GFF3. Because the file comes as chromosomes labelled simply by a number or letter (not chrx) as required, I wrote a script to pre-pend the "chr"s. For some reason this did the first few and not the whole lot, so Cufflinks appeared to work, but only delivered a subset of actual gene hits. I think I have corrected this now and am back in action. Thanks very much for your help.

ADD REPLYlink written 8 months ago by rbj30

Glad you worked this out! The help in this link has some help in case others run into this problem:

One last bit of advice now that you know what to look at - be sure to check the mito chromosome name - usually adding "chr" is not enough. Some sources have this as "MT" in the original and this will need to be "chrM" to work with UCSC-sourced genomes. Change this line-command or within Galaxy (whichever is easier for you).

Thanks! Jen

ADD REPLYlink written 8 months ago by Jennifer Hillman Jackson24k

Thank you Jen, OK had missed the MT chromosomes. I didn't find an easy safe way to change this (even following the link!), so wrote a Matlab script to change all the starting chromosome, hopefully appropriately. Seems to work now. Everything seems to have worked, but will cross check with the dataset analysed "properly". Kind Regards Richard

ADD REPLYlink written 8 months ago by rbj30
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 71 users visited in the last hour