Relatively large proportion of "LOWDATA", "FAIL" of FPKM_status running cufflink

[ruben6um]

Hello,

I am relatively new to analyses of RNA-seq. I am right now analyzing human blood data from 22 biological samples using the tophat and cufflinks pipeline. The cufflinks command I used to analysed the *.bam files (generated from tophat using hg19 reference) is

cufflinks -p 16 -o S01 -G Homo_sapiensPlusChr.GRCh37.63.gtf -b hg19.fa -u -N --compatible-hits-norm ./Sample1_accepted_hits.bam

Both the genes.fpkm_tracking and the isoform.fpkm_tracking resulting output files generated seem to have a relatively large proportion of "LOWDATA" and "FAIL" calls for the FPKM_status attribute.

This proportion of these calls seems similar (~30%) across the multiple samples and also the genes getting these calls seem almost the same again across the multiple samples.

I am not sure if I am doing something wrong - or if this is the expected behavior of the algorithms. I am hoping that I am (or the algorithms) are doing something incorrect.

We have matching microarray data generated from these samples. Some of the highly expressing genes from the microarray data get this "FAIL" status even though the FPKM values seem relatively high.

Any help would be appreciated.

Thanks,
-Reuben

[Jon_Keats]

I've encountered the same issue recently with the newest version. I'm actually testing if it is a package distribution issue right now as old cufflinks version produced values. In my case I came across a dianosis/relapse pair were TP53 went from FPKM 50 to 0 but really identical but for some reason the status is Fail in the relapse sample. In my case this occurred using the pre-compiled linux binary but in the test I did on my laptop with 10 million reads it worked on my Mac laptop with a version compiled from source. So now testing the entire set of reads on our Mac workstation to compare the output from our linux HPC resource

[Jon_Keats]

In my hands the "FAIL" result is completely inconsistent. For TP53, with the full dataset (37.34 million read pairs) you get a "FAIL" but with the first 10 million pairs FPKM of 42 with an "OK". The odd thing is there is a clear correlation with the number of reads as there were 5546 genes listed as "FAIL" using 10 million reads and 8059 with the full dataset.... To make matters worse the number of "FAIL" genes moves around as you change the normalization parameters. Right now I'm not too impressed and I think I'll go back to calculating FPKM by hand and using other packages. I can understand there being computation issues with fragment handling but at the gene level it seems pretty straight forward.

[oliviera]

Dear Jon,
Which package would you suggest as an alternative? I am not convince I can use DESeq or EdgeR with paired end data. What do you think?

Olivier