Seqanswers Leaderboard Ad

Collapse

Announcement

Collapse
No announcement yet.
X
 
  • Filter
  • Time
  • Show
Clear All
new posts

  • Solution found: For Cuffdiff/links 2.0.2 Make cufflinks FPKM match Cuffdiff FPKM

    This is in reference to cuffdiff v2.0.2 (3522)

    Hi, I have been banging my head on why cuffdiff and cufflinks have different FPKM calculations. I finally stumbled upon the settings to make cufflinks match the FPKM values created with default settings of cuffdiff:

    cuffdiff (defaults) mygenes.gtf mysample1.bam mysample1.bam
    cufflinks --compatible-hits-norm --max-frag-multihits 1 -G mygenes.gtf mysample1.bam

    I think there is a mistake in the cuffdiff options list, where --max-frag-multihits is listed as unlimited in the default state, but it is probably actually set to 1.

    also, I add a pseudocount RPKM of 0.01 to both, otherwise cufflinks makes very small RPKMs, unlike cuffdiff, where it sets it to zero.



    At least now it is possible for me to reproduce the FPKMs that cuffdiff makes on default settings. I hope the authors fix this.

    Note the plot is on log2 scale for the RPKMs, but the match is exact when on the linear scale as well (with exception of a few - those genes with extremely small FPKMs ~ e-200 as calculated by cufflinks).

    btw, I carefully compared all the options and saw that the options for --compatible-hits-norm and --total-hits-norm are flipped between the two programs:

    CUFFLINKS: --compatible-hits-norm count hits compatible with reference RNAs only [ default: FALSE ]
    CUFFDIFF: --compatible-hits-norm count hits compatible with reference RNAs only [ default: TRUE ]

    CUFFLINKS: --total-hits-norm count all hits for normalization [ default: TRUE ]
    CUFFDIFF: --total-hits-norm count all hits for normalization [ default: FALSE ]

    Interestingly, there is another normalization only offered by Cuffdiff and set to true:

    CUFFDIFF only: --geometric-norm use geometric mean normalization [ default: TRUE ]


    Both Cufflinks and Cuffdiff say the default for max-frag-multihits is unlimited, but I have found setting it to 1 in Cufflinks brings the two into agreement:

    CUFFLINKS --max-frag-multihits Maximum number of alignments allowed per fragment [ default: unlim ]
    CUFFDIFF --max-frag-multihits Maximum number of alignments allowed per fragment [ default: unlim ]
    Last edited by NGSfan; 12-27-2012, 10:18 AM.

  • #2
    Thanks for sharing

    Comment


    • #3
      Alternatively, instead of pseudocounts to deal with some of the very small RPKMs generated by Cufflinks, you can also replace any RPKM < 0.0001 with 0 and this will give you the same values as reported by Cuffdiff.

      I don't know why on some genes, the RPKMs get so small (1e-300 for example) in Cufflinks, whereas in Cuffdiff they are 0 RPKM.

      Has anybody compared Cufflinks to eXpress? I notice that there are subsets of genes where cufflinks reports 1000 RPKM while eXpress reports 0 and vice versa... many look to be small RNAs.. just wondering what others are seeing out there. This happens with or without effective length correction.

      I will try RSEM next.

      I am a bit concerned though with aligning only to transcripts, wouldn't this result in some reads being forcibly aligned to a transcript reference? A lot of transcription goes on that is not annotated (as we know from the ENCODE project) and this could result in skewing RPKMs for some genes. Perhaps the genomic reference approach is better? With accurate aligners like STAR the argument for aligning solely to the transcriptome seems weaker.
      Last edited by NGSfan; 12-28-2012, 02:27 PM.

      Comment


      • #4
        Very Useful I will try it now. Just wondering why you make multi-hit only 1, this seems odd because the aligner will stop for looking of other possibility in mapping when it encounter the first hit.

        Comment


        • #5
          Originally posted by fuad193 View Post
          Very Useful I will try it now. Just wondering why you make multi-hit only 1, this seems odd because the aligner will stop for looking of other possibility in mapping when it encounter the first hit.

          If my interpretation or understanding of the --max-frag-multihits 1 option is correct, this tells cuffdiff to limit itself to using reads that align uniquely to one locus, and to ignore reads that map to more than 1 locus. In the case of tophat/STAR, reads that are uniquely mapped get a score of 255, so it should be counting those only.

          I have not played around with changing this value. I would assume that unique reads are better to use for differential gene expression.

          Comment

          Latest Articles

          Collapse

          • seqadmin
            Non-Coding RNA Research and Technologies
            by seqadmin




            Non-coding RNAs (ncRNAs) do not code for proteins but play important roles in numerous cellular processes including gene silencing, developmental pathways, and more. There are numerous types including microRNA (miRNA), long ncRNA (lncRNA), circular RNA (circRNA), and more. In this article, we discuss innovative ncRNA research and explore recent technological advancements that improve the study of ncRNAs.

            Nobel Prize for MicroRNA Discovery
            This week,...
            10-07-2024, 08:07 AM
          • seqadmin
            Recent Developments in Metagenomics
            by seqadmin





            Metagenomics has improved the way researchers study microorganisms across diverse environments. Historically, studying microorganisms relied on culturing them in the lab, a method that limits the investigation of many species since most are unculturable1. Metagenomics overcomes these issues by allowing the study of microorganisms regardless of their ability to be cultured or the environments they inhabit. Over time, the field has evolved, especially with the advent...
            09-23-2024, 06:35 AM

          ad_right_rmr

          Collapse

          News

          Collapse

          Topics Statistics Last Post
          Started by seqadmin, 10-02-2024, 04:51 AM
          0 responses
          104 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 10-01-2024, 07:10 AM
          0 responses
          112 views
          0 likes
          Last Post seqadmin  
          Started by seqadmin, 09-30-2024, 08:33 AM
          1 response
          115 views
          0 likes
          Last Post EmiTom
          by EmiTom
           
          Started by seqadmin, 09-26-2024, 12:57 PM
          0 responses
          22 views
          0 likes
          Last Post seqadmin  
          Working...
          X