天天看點

生物資訊學習——cufflinks一. 簡介二. 安裝三. Cufflinks的使用四. Cuffmerge的使用五. Cuffcompare的使用六. Cuffdiff的使用八 Cuffquant九  Cuffnorm

原文位址: http://blog.sina.com.cn/s/blog_751bd9440102v72b.html

一. 簡介

Cufflinks下主要包含cufflinks,cuffmerge,cuffcompare和cuffdiff等幾支主要的程式。主要用于基因表達量的計算和差異表達基因的尋找。

二. 安裝

Cufflinks下載下傳網頁。

1. 為了安裝Cufflinks,必須有Boost C++libraries。下載下傳Boost并安裝。預設安裝在/usr/local。

$ tar jxvf boost_1_53_0.tar.bz2
$ cd boost_1_53_0
$ ./bootstrap.sh
$ sudo ./b2 install
      
2.安裝SAM tools。
下載下傳SAM tools。
$ tar jxvf samtools-0.1.18.tar.bz2
$ cd samtools-0.1.18
$ make
$ sudo su 
# mkdir /usr/local/include/bam
# cp libbam.a /usr/local/lib
# cp *.h /usr/local/include/bam/
# cp samtools /usr/bin/
      
3. 安裝 Eigen libraries。
下載下傳Eigen
$ tar jxvf 3.1.2.tar.bz2
$ cd eigen-eigen-5097c01bcdc4
$ sudo cp -r Eigen/ /usr/local/include/
      
4. 安裝Cufflinks。
$ tar zxvf cufflinks-2.0.2.tar.gz
$ cd cufflinks-2.0.2
$ ./configure --prefix=/path/to/cufflinks/install --with-boost=/usr/local/ --with-eigen=/usr/local/include//Eigen/
$ make
$ make install
      
5. 可以直接下載下傳Linux x86_64 binary。不需要上述繁瑣步驟,解壓後的程式直接可用。(推薦)

三. Cufflinks的使用

1. Cufflinks簡介

Cufflinks程式主要根據Tophat的比對結果,依托或不依托于參考基因組的GTF注釋檔案,計算出(各個gene的)isoform的FPKM值,并給出trascripts.gtf注釋結果(組裝出轉錄組)。

注意:

1.fragment的長度的估測,若為pair-end測序,則cufflinks自己會有一套算法,算出結果。若為single-end測序,則cufflinks預設的是高斯分布,或者你自己提供相關的參數設定。

2. cufflinks計算multi-mapped reads,一般a read map到10個位置,則每個位置記為10%。aread mapping to 10 positions will count as 10% of a read at eachposition.

3. 一般不推薦用cufflinks拼接細菌的轉錄組,推薦 Glimmer。但是,若有注釋檔案,可以用cufflinks和cuffdiff來檢測基因的表達和差異性。

4. cufflinks/cuffdiff不能計算出exon或splicing event的FPKM

5.cuffdiff處理時間序列data:采用參數-t

6.當你使用cufflinks時,在最後出現了99%,然後一直不動。因為cuffdiff需要更多的CPU來處理一些比對很多reads的loci。而這些位點一般要等其他位點全部解決了後,才由cuffdiff來處理。可以用參數-M來提供相關的檔案,過濾掉rRNA或者線粒體RNA。

7. 當使用cufflinks或cuffdiff出現了“crash with a ‘bad_alloc'error”,cuffdiff和cufflinks運作了很長時間才結束————這表明計算機拼接一個高表達的基因或定量分析一個高表達的基因,運作的記憶體使用玩盡了!解決方法:修改選項“-max-bundle-frags”,可以先嘗試500000,若錯誤依舊在,可以繼續下調!

8. cuffdiff報道的結果裡面所有的基因和轉錄本的FPKM=0,這表明GTF中的染色體名字和BAM裡的名字不比對。

9. cuffdiff和cufflinks的缺點:存在一定的假基因和轉錄本(原因:測序深度,測序品質,測序樣本的測序次數,以及注釋的錯誤)

10. large foldchange表達量不代表資料的明顯性(這些基因的isform多或這些基因測序測到的少,整體較低的表達)。cuffdiff中明顯表達倍數改變的基因,存在不确定性。

11. 通過cufflinks産生的結果中transcript.gtf檔案中cuff辨別的轉錄本就是新的轉錄本。相應的,其他子產品輸出中CUFF辨別代表着新的轉錄本。

12. 若出現了如下錯誤:

You  are  using  Cufflinks  v2.2.1,  which  is  the  most  recent  release.

open:  No  such  file  or  directory

File  30  doesn't  appear  to  be  a  valid  BAM  file,  trying  SAM...

Error:  cannot  open  alignment  file  30  for  reading

這表明,你的參數有問題。例如“--min-intron-length”,你設定為了:“-min-intron-length”

2. 使用方法

$ cufflinks [options]* 

一個常用的例子:
$ cufflinks -p 8 -G transcript.gtf --library-type fr-unstranded -o cufflinks_output tophat_out/accepted_hits.bam
      

3. 普通參數

   -o | --output-dir   default: ./
    設定輸出的檔案夾名稱 
      
 
-p | --num-threads  default: 1
    用于比對reads的CPU線程數 
      
 
-G | --GTF 
    提供一個GFF檔案,以此來計算isoform的表達。此時,将不會組裝新的transcripts,
程式會忽略和reference transcript不相容的比對結果 
      
 
-g | --GTF-guide 
    提供GFF檔案,以此來指導轉錄子組裝(RABT assembly)。此時,輸出結果會包含reference transcripts和novel genes and isforms。 
      
 
-M | --mask-file 
    提供GFF檔案。Cufflinks将忽略比對到該GTF檔案的transcripts中的reads。該
檔案中常常是rRNA的注釋,也可以包含線立體和其它希望忽略的transcripts的注釋。将這些不需要的RNA去除後,對計算mRNA的表達量是有利的。 
      
 
-b | --frag-bias-correct 
    提供一個fasta檔案來指導Cufflinks運作新的bias detection and correction algorithm。這樣能明顯提高轉錄子豐度計算的精确性。 
      
 
-u | --multi-read-correct
    讓Cufflinks來做initial estimation步驟,進而更精确衡量比對到genome多個位點的reads。 
      
 
--library-type  default:fr-unstranded
    處理的reads具有鍊特異性。比對結果中将會有個XS标簽。一般Illumina資料的lib
rary-type為 fr-unstranded。

--library-norm-method    具體參考官網,三種方式:classic-fpkm  預設的方式。geometric  針對DESeq。quartile  計算時,fragments和總的map的count取75%

      

4. 豐度評估參數

-m | --frag-len-mean default: 200
插入片段的平均長度。不過現在Cufflinks能learns插入片段的平均長度,是以不推薦自主
設定此值。 
      
 
-s | --frag-len-std-dev default: 80
插入片段長度的标準差。不過現在Cufflinks能learns插入片段的平均長度,是以不推薦自
主設定此值。 
      
 
-N | --upper-quartile-form
使用75%分為數的值來代替總的值(比對到單一位點的fragments的數值),作normalize。這樣有利于在低豐度基因和轉錄子中尋找差異基因。 
      
 
--total-hits-norm default: TRUE
Cufflinks在計算FPKM時,算入所有的fragments和比對上的reads。和下一個參數
對立。預設激活該參數。 
      
 
--compatible-hits-norm 
Cufflinks在計算FPKM時,隻針對和reference transcripts相容的fragments以及比對上的reads。該參數預設不激活,隻能在有 --GTF 參數下有效,并且作 RABT
或 ab initio 的時候無效。
      

5. 組裝常用參數

-L | --label  default: CUFF
    Cufflink以GTF格式來報告轉錄子片段(transfrags),該參數是GTF檔案的字首 
      
-a/--junc-alpha <0.0-1.0>    剪接比對過濾中假陽性的二項檢驗中的 alpha value。預設為 0.001
      
--min-frags-per-transfrag   default: 10
    組裝出的transfrags被支援的RNA-seq的fragments數少于該值則不被報道。 
      
--min-intron-length   default: 50
    最小的intron大小。 
      
--trim-3-dropoff-frac    最低百分比的拼接的轉錄本的3‘端的平均覆寫程度。預設0.1  The fraction of average coverage below which to trim the 3' end of an assembled          transcript.  The default is 0.1.
      
--max-multiread-fraction <0.0-1.0>   若一個轉錄本Transfrags的reads能比對到基因組的多個位置,其中該轉錄本的reads有超過該百分比是multireads,則不會報告這個轉錄本。預設為75%   The fraction a transfrag's supporting reads that may be multiply mapped to the genome. A transcript composed of more than this fraction will not be reported by the assembler.  Default: 0.75 (75% multireads or more is suppressed). 
      
--overlap-radius   default: 50
    Transfrags之間的距離少于該值,則将其連到一起。
      
--intron-overhang-tolerance    當決定一個拼接的轉錄本(這個轉錄本可能不是新的轉錄本)和一個參考轉錄本是否合并時,參考轉錄本的外顯子允許延伸的長度。預設50bp   The number of bp allowed to enter the intron of a reference transcript when determining if an     assembled transcript should be merged with it (ie, the assembled transcript is not novel).      The default is 50 bp.
      
其他參數(無關緊要)
      
-v/--verbose   顯示版本資訊等等
      

6. Cufflinks輸出結果

cufflinks的輸入檔案是sam或bam格式。并且sam或bam格式的檔案必須排好序。(The SAM file supplied to Cufflinks must be sorted by 
          reference position.)Tophat的輸出結果sam或bam已經排好了序。針對其他的未排序的sam或bam檔案采用如下排序方式:

sort -k 3,3 -k 4,4n hits.sam > hits.sam.sorted



1. 
transcripts.gtf

該檔案包含Cufflinks的組裝結果isoforms。前7列為标準的GTF格式,最後一列為attributes。其每一列的意義:
      
列數   列的名稱  例子         描述


1     序列名    chrX        染色體或contig名
; 2     來源      Cufflinks   産生該檔案的程式名
; 3     類型      exon        記錄的類型,一般是transcript或exon
; 4     起始      1           1-base的值
; 5     結束      1000        結束位置
; 6     得分      1000        
; 7     鍊        +          Cufflinks猜測isoform來自參考序列的那一條鍊,
一般是'+','-'或'.';  
8     frame    .           Cufflinks不去預測起始或終止密碼子框的位置
; 9     attributes  ...      詳見下
      
每一個GTF記錄包含如下attributes:
      
Attribute      例子       描述


gene_id          CUFF.1      Cufflinks的gene id
;  transcript_id    CUFF.1.1  Cufflinks的轉錄子 id  
; FPKM           101.267   isoform水準上的豐度, Fragments Per Kilobase
 of exon model per Million mapped fragments
; frac           0.7647    保留着的一項,忽略即可,以後可能會取消這個; 
conf_lo        0.07      isoform豐度的95%置信區間的下邊界,即 下邊界值 =
 FPKM * ( 1.0 - conf_lo )
;  conf_hi        0.1102    isoform豐度的95%置信區間的上邊界,即 上邊界值 =
 FPKM * ( 1.0 + conf_hi )
; cov            100.765   計算整個transcript上read的覆寫度; 
full_read_support   yes  當使用 RABT assembly 時,該選項報告所有的intr
ons和exons是否完全被reads所覆寫
      
2. ispforms.fpkm_tracking

isoforms(可以了解為gene的各個外顯子)的fpkm計算結果

3. 
genes.fpkm_tracking

gene的fpkm計算結果
      

四. Cuffmerge的使用

1. Cuffmerge簡介

Cuffmerge将各個Cufflinks生成的transcripts.gtf檔案融合稱為一個更加全面的transcripts注釋結果檔案merged.gtf。以利于用Cuffdiff來分析基因差異表達。
      

2. 使用方法

$ cuffmerge [options]* 
輸入檔案為一個文本檔案,是包含着GTF檔案路徑的list。常用例子:
$ cuffmerge -o ./merged_asm -p 8 assembly_list.txt
      

3. 使用參數

-h | --help


-o  default: ./merged_asm
将結果輸出至該檔案夾。

 -g | --ref-gtf
将該reference GTF一起融合到最終結果中。



-p | --num-threads  defautl: 1
使用的CPU線程數


-s | --ref-sequence /
該參數指向基因組DNA序列。如果是一個檔案夾,則每個contig則是一個fasta檔案;如果是
一個fasta檔案,則所有的contigs都需要在裡面。Cuffmerge将使用該ref-sequence來
幫助對transfrags分類,并排除repeats。比如transcripts包含一些小寫堿基的将歸類
到repeats.
      

4. Cuffmerge輸出結果

輸出的結果檔案預設為 /merged.gtf
      

五. Cuffcompare的使用

1. Cuffcompare簡介

Cuffcompare使用Cufflinks的GTF結果,對GTF結果進行比較。和reference gtf比較尋找novel轉錄本等。
      

2. Cuffcompare的使用方法

$ cuffcompare [options]*  [cuff2.gtf] ... [cuffN.gtf]

使用例子:
$ cuffcompare -o cuffcmp cuff1.gtf cuff2.gtf
      

3. 使用參數

-h                -V  顯示程序    

-C  
預設,表示"contained" transcripts 也會寫入.combined.gtf中。
-o  default: cuffcmp
輸出檔案的字首


-r 
參考的GFF檔案。用來評估輸入的gtf檔案中gene models的精确性。每一個輸入的gtf的isoforms将和該參考檔案進行比較,并被标注為 overlapping, matching 或 novel。

 -R
當有了 -r 參數時,指定該參數時,将忽略參考GFF檔案中的一些transcripts。這些transcripts不和任何輸入的GTF檔案overlapped。



-s   該參數指向基因組DNA序列。如果是一個檔案夾,則每個contig則是一個fasta檔案;如果是
一個fasta檔案,則所有的contigs都需要在裡面。小寫字母的堿基用來将相應的transcripts作為repeats處理。
      

4. 輸出結果

在目前目錄下輸出3個檔案:

.stats, 報告與參考注釋比較時,各種與準确性相關的資料。其中,Sn和Sp展示的是specificity and sensitivity values。 fSn and fSp 列展示的 "fuzzy" variants of these same accuracy calculations。允許存在變動。(-o 沒有設定,預設為cuffcmp為檔案字首)

.combined.gtf    報告每個樣本的所有的 transfrags 的資訊。若一個transfrag在多個樣本中,它隻報道一次。

 .tracking      比對到樣本間的轉錄本。this file matches transcripts up between samples.  Each row contains 
                a transcript structure that is present in one or more input GTF files. 
                Because the transcripts will generally have different IDs (unless you 
                assembled your RNA-Seq reads against a reference transcriptome), 
                cuffcompare examines the structure of each the transcripts, 
                matching transcripts that agree on the coordinates and order of all of
                their introns, as well as strand.  Matching transcripts are allowed to 
                differ on the length of the first and last exons, since these lengths
                will naturally vary from sample to sample due to the random nature of 
                sequencing.

例子;



TCONS_00000045 XLOC_000023 Tcea|uc007afj.1    j            q1:exp.115|exp.115.0|100|3.061355|0.350242|0.350207      q2:60hr.292|60hr.292.0|100|4.094084|0.000000|0.000000
      
In this example, a transcript present in the two input files,called exp.115.0 in the first and60hr.292.0 inthe second, doesn't match any reference transcript exactly, butshares exons withuc007afj.1, an isoform of the gene Tcea,as indicated by the classcodej. The first three columns are as follows:
其中,1 Cufflinks transfrag id  TCONS_00000045  内部的transfrag id;2  Cufflinks locus id  XLOC_000023  内部的locus id; 3  Reference gene id   Tcea   參考的注釋的gene的id或者“-”表示沒有比對到參考的轉錄本; 4  Reference transcript id  uc007afj.1  參考的注釋的轉錄本的id或者“-”表示沒有比對到參考的轉錄本 ; 5 Class code  c  轉錄本和參考轉錄本之間的比對類型。第五列之後如下:

      
qJ: | | | | | | |
在輸入的GTF的同目錄下輸出.refmap 和 
.tmap 檔案。

.refmap  具體内容如下:

1  Reference gene name   參考注釋的gtf中的基因名字 2 Reference transcript id 參考的轉錄本id  3  Class code 表示cufflinks拼接的轉錄本和參考轉錄本間的比對情況:c 表示部分比對;= 表示全部比對

4  Cufflinks matches  比對到參考轉錄本的cufflinks拼接的轉錄本的id



.tmap  具體内容如下:

1  Reference gene name   參考注釋的gtf中的基因名字 2 Reference transcript id 參考的轉錄本id  3  Class code 表示cufflinks拼接的轉錄本和參考轉錄本間的比對情況:c 表示部分比對;= 表示全部比對

4 Cufflinks gene id  ; 5 Cufflinks transcript id;  6 Fraction of major isofor m (FMI) ; 7  FPKM ; 8 FPKM_conf_lo; 9  FPKM_conf_hi  ; 10 Coverage ; 11 Length; 12  Major isoform ID



class cord :

      
Priority Code Description
1 = Complete match of intron chain
2 c Contained
3 j Potentially novel isoform (fragment): at least onesplice junction is shared with a reference transcript
4 e Single exon transfrag overlapping a reference exonand at least 10 bp of a reference intron, indicating a possiblepre-mRNA fragment.
5 i A transfrag falling entirely within a referenceintron
6 o Generic exonic overlap with a referencetranscript
7 p Possible polymerase run-on fragment (within2Kbases of a reference transcript)
8 r Repeat. Currently determined by looking at thesoft-masked reference sequence and applied to transcripts where atleast 50% of the bases are lower case
9 u Unknown, intergenic transcript
10 x Exonic overlap with reference on the oppositestrand
11 s An intron of the transfrag overlaps a referenceintron on the opposite strand (likely due to read mappingerrors)
12 . (.tracking file only, indicates multipleclassifications)

六. Cuffdiff的使用

1. Cuffdiff簡介

用于尋找轉錄子表達的顯著性差異。
      

2. Cuffdiff使用方法

cuffdiff主要是發現轉錄本表達,剪接,啟動子使用的明顯變化。

cuffdiff [options]* ...[sampleN.sam_replicate1.sam[,...,sample2_replicateM.sam]]

$ cuffdiff [options]*   ...[sampleN_1.sam[,...,sampleN_M.sam]]
其中transcripts.gtf是由cufflinks,cuffcompare,cuffmerge所生成的檔案,或是由其它程式生成的。一個樣本有多個replicate,用逗号隔開。sample多于一個時,cuffdiff将比較samples間的基因表達的差異性。

一個常用例子:
$ cuffdiff --lables lable1,lable2 -p 8 --time-series --multi-read-correct --library-type fr-unstranded --poisson-dispersion transcripts.gtf sample1.sam sample2.sam

cuffdiff接受bam/sam或cuffquant的CXB檔案,同時也可以接受bam與sam的混合檔案,不能接受bam/sam和CXB的混合檔案。

      

3. 使用參數

-h | --help


-o | --output-dir  default: ./
輸出的檔案夾目錄。


-L | --lables   default: q1,q2,...qN
給每個sample一個樣品名或者一個環境條件一個lable



-p | --num-threads  default: 1
使用的CPU線程數



-T | --time-series
讓Cuffdiff來按樣品順序來比對樣品,而不是對所有的samples都進行兩兩比對。即第二個
SAM和第一個SAM比;第三個SAM和第二個SAM比;第四個SAM和第三個SAM比...



-N | --upper-quartile-form
使用75%分為數的值來代替總的值(比對到單一位點的fragments的數值),作normalize。
這樣有利于在低豐度基因和轉錄子中尋找差異基因。



--total-hits-norm 
Cufflinks在計算FPKM時,算入所有的fragments和比對上的reads。和下一個參數對立。
預設不激活該參數。

 --compatible-hits-norm
Cufflinks在計算FPKM時,隻針對和reference transcripts相容的fragments以及
比對上的reads。該參數預設激活,使用該參數可以降低核糖體rna的reads對基因表達的幹擾。

 -b | --frag-bias-correct(一般是genome.fa)
提供一個fasta檔案來指導Cufflinks運作新的bias detection and correction 
algorithm。這樣能明顯提高轉錄子豐度計算的精确性。

 -u | --multi-read-correct
讓Cufflinks來做initial estimation步驟,進而更精确衡量比對到genome多個位點
的reads。



-c | --min-alignment-count   default: 10
如果比對到某一個位點的fragments數目少于該值,則不做該位點的顯著性分析。認為該位點的表達量沒有顯著性差異。



-M | --mask-file 
提供GFF檔案。Cufflinks将忽略比對到該GTF檔案的transcripts中的reads。該檔案中常常是rRNA的注釋,也可以包含線立體和其它希望忽略的transcripts的注釋。将這些不需要的RNA去除後,對計算mRNA的表達量是有利的。


-FDR  default: 0.05
允許的false discovery rate.


--library-type default:fr-unstranded
處理的reads具有鍊特異性。比對結果中将會有個XS标簽。一般Illumina資料的library-
type為 fr-unstranded。


--dispersion-method   


其他進階參數:

-m | --frag-len-mean default: 200
插入片段的平均長度。不過現在Cufflinks能learns插入片段的平均長度,是以不推薦自主
設定此值。



-s | --frag-len-std-dev default: 80
插入片段長度的标準差。不過現在Cufflinks能learns插入片段的平均長度,是以不推薦自
主設定此值。


-v/--verbose   顯示版本資訊等等

 -q/--quiet     除了警告和錯誤外,其他資訊将不會print


--no-update-check   關系cufflinks自動更新的能力


-F/--min-isoform-fraction <0.0-1.0>   建議不要更改,主要的isorform豐度若低于這個分數,可變的isoform将四舍五入為0.預設為1e-5

--max-bundle-frags   一個skipped locus/loci在skipped前可以擁有的最大的fragment片段。預設為1000000  

--max-frag-count-draws (預設為100)和--max-frag-assign-draws (預設為50)
--min-reps-for-js-test      一個針對不同調控的基因做test的最小的複制次數。Cuffdiff won't test genes for differential regulation unless the 
conditions in question have at least this many replicates.  Default: 3. 

--no-effective-length-correction   Cuffdiff will not employ its "effective" length normalization to transcript FPKM. Cufflinks将不會使用它的“effective” 長度标準化去計算轉錄的FPKM

--no-length-correction    cufflinks将根本不會使用轉錄本的長度去标準化fragment的數目。當fragment的數目和the features being quantified的size是獨立的,可以使用(例如for small RNA libraries, where no fragmentation takes place, or 3 prime end sequencing, where sampled RNA fragments are all essentially the same length).小心使用

--max-mle-iterations       極大似然法的疊代次數,預設5000
--poisson-dispersion
Use the Poisson fragment dispersion model instead of learning one 
in each condition.



      

4. Cuffdiff輸出

1. FPKM tracking files   cuffdiff計算每個樣本中的轉錄本,初始轉錄本和基因的FPKM。其中,基因和初始轉錄本的FPKM的計算是在每個轉錄本group和基因group中的轉錄本的FPKM的求和。

      
isoforms.fpkm_tracking Transcript FPKMs
genes.fpkm_tracking Gene FPKMs. Tracks the summed FPKM of transcriptssharing each gene_id
cds.fpkm_tracking Coding sequence FPKMs. Tracks the summed FPKM oftranscripts sharing eachp_id, independent oftss_id
tss_groups.fpkm_tracking Primary transcript FPKMs. Tracks the summed FPKMof transcripts sharing eachtss_id
2. Count tracking files    評估每個樣本中來自每個 transcript, primary transcript, 
                and gene的fragment數目。其中primary transcript, 
                and gene的fragment數目是每個primary transcript group或gene group中trancript的數目之和。

      
isoforms.count_tracking Transcript counts
genes.count_tracking Gene counts. Tracks the summed counts oftranscripts sharing eachgene_id
cds.count_tracking Coding sequence counts. Tracks the summed countsof transcripts sharing eachp_id, independent oftss_id
tss_groups.count_tracking Primary transcript counts. Tracks the summedcounts of transcripts sharing eachtss_id
 3. Read group tracking 
files   計算在每個repulate中每個transcript, primary transcript和gene的表達量和frage數目

      
isoforms.read_group_tracking Transcript read group tracking
genes.read_group_tracking Gene read group tracking. Tracks the summedexpression and counts of transcripts sharing eachgene_idin each replicate
cds.read_group_tracking Coding sequence FPKMs. Tracks the summedexpression and counts of transcripts sharing eachp_id,independent of tss_id in each replicate
tss_groups.read_group_tracking Primary transcript FPKMs. Tracks the summedexpression and counts of transcripts sharing eachtss_idin each replicate
4. Differential expression test    對于splicing transcript,
                primary transcripts, genes, and coding sequences.樣本之間的表達差異檢驗。對于每一對樣本x和y,都會有以下四個檔案:

      
isoform_exp.diff Transcript differential FPKM.
gene_exp.diff Gene differential FPKM. Tests difference sin thesummed FPKM of transcripts sharing eachgene_id
tss_group_exp.diff Primary transcript differential FPKM. Testsdifferences in the summed FPKM of transcripts sharing eachtss_id
cds_exp.diff Coding sequence differential FPKM. Testsdifferences in the summed FPKM of transcripts sharing eachp_id independent oftss_id
每個檔案的樣式如下:

      
Column number Column name Example Description
1 Tested id XLOC_000001 A unique identifier describing the transcipt,gene, primary transcript, or CDS being tested
2 gene Lypla1 The gene_name(s) or gene_id(s)being tested
3 locus chr1:4797771-4835363 Genomic coordinates for easy browsing to the genesor transcripts being tested.
4 sample 1 Liver Label (or number if no labels provided) of thefirst sample being tested
5 sample 2 Brain Label (or number if no labels provided) of thesecond sample being tested
6 Test status NOTEST Can be one of OK (test successful), NOTEST (notenough alignments for testing), LOWDATA (too complex or shallowlysequenced), HIDATA (too many fragments in locus), or FAIL, when anill-conditioned covariance matrix or other numerical exceptionprevents testing.
7 FPKMx 8.01089 FPKM of the gene in sample x
8 FPKMy 8.551545 FPKM of the gene in sample y
9 log2(FPKMy/FPKMx) 0.06531 The (base 2) log of the fold changey/x
10 test stat 0.860902 The value of the test statistic used to computesignificance of the observed change in FPKM
11 p value 0.389292 The uncorrected p-value of thetest statistic
12 q value 0.985216 The FDR-adjusted p-value of thetest statistic
13 significant no Can be either "yes" or "no", depending on whetherp is greater then the FDRafter Benjamini-Hochbergcorrection for multiple-testing
5. Differential splicing tests – 
splicing.diff     對于每個primary transcript,鑒定的不同的isoform的差異性。隻有2個或2個以上的isoforms的primary transcript存在

      
Column number Column name Example Description
1 Tested id TSS10015 A unique identifier describing the primarytranscript being tested.
2 gene name Rtkn The gene_name or gene_id thatthe primary transcript being tested belongs to
3 locus chr6:83087311-83102572 Genomic coordinates for easy browsing to the genesor transcripts being tested.
4 sample 1 Liver Label (or number if no labels provided) of thefirst sample being tested
5 sample 2 Brain Label (or number if no labels provided) of thesecond sample being tested
6 Test status OK Can be one of OK (test successful), NOTEST (notenough alignments for testing), LOWDATA (too complex or shallowlysequenced), HIDATA (too many fragments in locus), or FAIL, when anill-conditioned covariance matrix or other numerical exceptionprevents testing.
7 Reserved
8 Reserved
9 √JS(x,y) 0.22115 The splice overloading of the primary transcript,as measured by the square root of the Jensen-Shannon divergencecomputed on the relative abundances of the splice variants
10 test stat 0.22115 The value of the test statistic used to computesignificance of the observed overloading, equal to √JS(x,y)
11 p value 0.000174982 The uncorrected p-value of thetest statistic.
12 q value 0.985216 The FDR-adjusted p-value of thetest statistic
13 significant yes Can be either "yes" or "no", depending on whetherp is greater then the FDRafter Benjamini-Hochbergcorrection for multiple-testing
6. Differential coding output – cds.diff    對于每個基因,它的cds的鑒定。樣本間的輸出cds的差異性。隻有2個或2個以上的cds(multi-protein genes)列舉在檔案中。

      
Column number Column name Example Description
1 Tested id XLOC_000002-[chr1:5073200-5152501] A unique identifier describing the gene beingtested.
2 gene name Atp6v1h The gene_name or gene_id
3 locus chr1:5073200-5152501 Genomic coordinates for easy browsing to the genesor transcripts being tested.
4 sample 1 Liver Label (or number if no labels provided) of thefirst sample being tested
5 sample 2 Brain Label (or number if no labels provided) of thesecond sample being tested
6 Test status OK Can be one of OK (test successful), NOTEST (notenough alignments for testing), LOWDATA (too complex or shallowlysequenced), HIDATA (too many fragments in locus), or FAIL, when anill-conditioned covariance matrix or other numerical exceptionprevents testing.
7 Reserved
8 Reserved
9 √JS(x,y) 0.0686517 The CDS overloading of the gene, as measured bythe square root of the Jensen-Shannon divergence computed on therelative abundances of the coding sequences
10 test stat 0.0686517 The value of the test statistic used to computesignificance of the observed overloading, equal to √JS(x,y)
11 p value 0.00546783 The uncorrected p-value of thetest statistic
12 q value 0.985216 The FDR-adjusted p-value of thetest statistic
13 significant yes Can be either "yes" or "no", depending on whetherp is greater then the FDRafter Benjamini-Hochbergcorrection for multiple-testing
7. Differential 
promoter use – promoters.diff  樣本間啟動子使用的差異性。隻有表達2個或2個以上isoform的基因列舉在這裡。

8. Read group info – read_groups.info   每個repulate,在進行定量分析時,cuffdiff的關鍵屬性會列出。

      
Column number Column name Example Description
1 file mCherry_rep_A/accepted_hits.bam BAM or SAM file containing the data for the readgroup
2 condition mCherry Condition to which the read group belongs
3 replicate_num Replicate number of the read group
4 total_mass 4.72517e+06 Total number of fragments for the read group
5 norm_mass 4.72517e+06 Fragment normalization constant used duringcalculation of FPKMs.
6 internal_scale 1.23916 Internal scaling factor, used to transformreplicates of a single condition onto the "internal" common countscale.
7 external_scale 0.96 External scaling factor, used to transform countsfrom different conditions onto an internal common count scale.
9. Run 
info – run.info   運作的資訊。



其中:輸出檔案FPKM Tracking file的格式如下:

      

1 tracking_id TCONS_00000001 内部唯一object的id(識别基因,轉錄本,CDS,初始轉錄本)Aunique identifier describing the object (gene, transcript, CDS,primary transcript)

2 class_code = 内部定義的類别的id,“-”表明不是轉錄本。Theclass_codeattribute for the object, or "-" if not a transcript, or ifclass_codeisn't present

3 nearest_ref_id NM_008866.1 最接近的參考轉錄本Thereference transcript to which the class code refers, ifany

4 gene_id NM_008866 基因id Thegene_id(s)associated with the object

5 gene_short_name Lypla1 基因名字 Thegene_short_name(s)associated with the object

6 tss_id TSS1 初始轉錄本id,或者“-”表示沒有初始轉錄本。Thetss_idassociated with the object, or "-" if not a transcript/primarytranscript, or iftss_idisn'tpresent

7 locus chr1:4797771-4835363 基因組上的位置Genomiccoordinates for easy browsing to the object

8 length 2447 轉錄本的長度Thenumber of base pairs in the transcript, or '-' if not atranscript/primary transcript

9 coverage 43.4279 read覆寫深度的估測值 Estimate for the absolute depth of read coverage across theobject

10 q0_FPKM 8.01089 樣本0中object的FPKM  FPKMofthe object in sample 0

11 q0_FPKM_lo 7.03583 object在樣本0中FPKM的95%置信區間的下界thelower bound of the 95% confidence interval on the FPKM of theobject in sample 0

12 q0_FPKM_hi 8.98595 object在樣本0中FPKM的95%置信區間的上界theupper bound of the 95% confidence interval on the FPKM of theobject in sample 0

13 q0_status OK object在樣本0中的量化狀态,0K表示成功,LOWDATA:太複雜或測序深度不夠;HIDATA:在一個基因座上太多fragments,FAIL:失敗的協方差矩陣或其他數值阻止了去卷積Quantificationstatus for the object in sample 0. Can be one of OK(deconvolutionsuccessful), LOWDATA (too complex or shallowly sequenced), HIDATA(too many fragments in locus), or FAIL, when an ill-conditionedcovariance matrix or other numerical exception preventsdeconvolution.

Count tracking files 格式如下:

1 tracking_id TCONS_00000001  A unique identifierdescribing the object (gene, transcript, CDS, primarytranscript)

2 q0_count 201.334 Estimated (externally scaled) number of fragments generated by theobject in sample 0

3 q0_count_variance  5988.24  Estimated variance inthe number of fragments generated by the object in sample0

4 q0_count_uncertainty_var  170.21  Estimated variance inthe number of fragments generated by the object in sample 0 due tofragment assignment uncertainty.

5 q0_count_dispersion_var  4905.63  Estimated variance inthe number of fragments generated by the object in sample 0 due tocross-replicate variability.

6 q0_status OK Quantification status for the object in sample 0. Can be one of OK(deconvolutionsuccessful), LOWDATA (too complex or shallowly sequenced), HIDATA(too many fragments in locus), or FAIL, when an ill-conditionedcovariance matrix or other numerical exception preventsdeconvolution.

七. cufflinks使用中遇到的問題

      
1. 
使用cuffdiff時候,在最新版本下,無重複的RNA-seq樣作比較,結果中沒有差異表達基因?

在v2.0.1及之後的版本中cuffdiff貌似不支援無重複的RNA-seq資料了。使用之前的版本即可。

      

八 Cuffquant

cuffquant是cuffquant能夠對單個 BAM檔案的基因轉錄本表達水準進行定量分析。生成的是CXB檔案abundances.cxb,,可以作為cuffdiff的輸入,這會加快cuffdiff的運作速度。也可以作為Cuffnorm的輸入。

具體使用:Usage: cuffquant [options]*

它的參數:(和前面參數的含義是一樣的)

-h/--help;-o/--output-dir;-p/--num-threads ;-M/--mask-file;-b/--frag-bias-correct;-u/--multi-read-correct;--library-type;-m/--frag-len-mean;-s/--frag-len-std-dev ;--max-mle-iterations;--max-bundle-frags;--no-effective-length-correction;--no-length-correction;-v/--verbose;-q/--quiet;--no-update-check;

九  Cuffnorm

cuffnorm能夠用 cuffquant的輸出檔案作為輸入檔案,對基因和轉錄組,簡單計算标準化過的表達水準。當你想要的是一系列可比較的基因、轉錄組、CDS 組和 TSS組的表達值時,可是使用 cuffnorm。例如,當你僅僅想對單個基因的表達值做個熱圖或者點圖時。

cuffnorm [options]* ...[sampleN.sam_replicate1.sam[,...,sample2_replicateM.sam]]

具體參數:它的參數和前面的類似,可以看前面的相關參數。

-h/--help  ;-o/--output-dir;-L/--labels ;-p/--num-threads ;

--total-hits-norm(預設不激活);--compatible-hits-norm(預設激活);--library-type;--library-norm-method;--output-format;-v/--verbose;-q/--quiet;--no-update-check;

cuffnorm的輸出檔案是實驗中的each gene, transcript, TSS group, andCDSgroup的标準化的表達水準。不做表達差異的分析。cuffnorm的輸出檔案預設是“simple-table”的檔案。這些檔案和cuffdiff輸出的檔案格式不同。若你想要cuffdiff格式的檔案,你需要輸入指令:--output-format cuffdiff

cuffnorm 報道FPKM values and normalized,estimates for the number of fragments that originate from eachgene, transcript, TSS group, and CDSgroup.這些結果已經做了标準化處理。對于某些下遊軟體需要原始檔案,是不作為其輸入的。

可以建立一個檔案,例如sample_sheet.txt作為cuffdiff或cuffnorm的輸入(存入sam檔案的path)。檔案格式如下:

sample_id      group_label


C1_R1.sam       C1


C1_R2.sam       C1


C2_R1.sam       C2


C2_R2.sam       C2



      

輸出結果檔案如下:

FPKM tracking files:估測的基因的表達水準

Count tracking files:估測的基因的fragment count values

Read group tracking files:報道per-replicate expression and countdata.

對于每個genes, transcripts, TSS groups, and CDSgroups,cuffnorm會報道兩種檔案形式: *.fpkm_table files and *.count_tablefiles。

繼續閱讀