如何自定義袖釦輸出欄

Question

我正在用mm10分析鼠標的RNAseq數據。我從tophat2/bowtie開始。然後我跑袖釦生成基因和異構體的FPKM。我正在使用一個gtf文件，它在第二列以及與基因名稱相鄰的行內具有基因生物型（即是否爲假基因，蛋白質編碼，snRNA，lincRna等）。行我GTF的例子是：

 
1 unprocessed_pseudogene exon 3054233 3054733 . + . exon_id "ENSMUSE00000848981"; exon_number "1"; gene_biotype "pseudogene"; gene_id "ENSMUSG00000090025"; gene_name "Gm16088"; gene_source "havana"; tag "mRNA_start_NF"; transcript_id "ENSMUST00000160944"; transcript_name "Gm16088-001"; transcript_source "havana"; tss_id "TSS82763"; 
1 unprocessed_pseudogene transcript 3054233 3054733 . + . gene_biotype "pseudogene"; gene_id "ENSMUSG00000090025"; gene_name "Gm16088"; gene_source "havana"; tag "mRNA_start_NF"; transcript_id "ENSMUST00000160944"; transcript_name "Gm16088-001"; transcript_source "havana"; tss_id "TSS82763"; 
1 snRNA exon 3102016 3102125 . + . exon_id "ENSMUSE00000522066"; exon_number "1"; gene_biotype "snRNA"; gene_id "ENSMUSG00000064842"; gene_name "Gm26206"; gene_source "ensembl"; transcript_id "ENSMUST00000082908"; transcript_name "Gm26206-201"; transcript_source "ensembl"; tss_id "TSS81070"; 
1 snRNA transcript 3102016 3102125 . + . gene_biotype "snRNA"; gene_id "ENSMUSG00000064842"; gene_name "Gm26206"; gene_source "ensembl"; transcript_id "ENSMUST00000082908"; transcript_name "Gm26206-201"; transcript_source "ensembl"; tss_id "TSS81070"; 

我的袖釦基因和亞型跟蹤輸出文件看起來是這樣的：

 
tracking_id  class_code  nearest_ref_id gene_id gene_short_name tss_id locus length coverage  FPKM FPKM_conf_lo FPKM_conf_hi FPKM_status 
ENSMUSG00000090025  -  -  ENSMUSG00000090025  Gm16088 TSS82763  1:3054232-3054733  -  -  0  0  0  OK 
ENSMUSG00000064842  -  -  ENSMUSG00000064842  Gm26206 TSS81070  1:3102015-3102125  -  -  0  0  0  OK 
ENSMUSG00000025900  -  -  ENSMUSG00000025900  Rp1  TSS11475  1:4343506-4360314  -  -  0  0  0  OK 
ENSMUSG00000088333  -  -  ENSMUSG00000088333  Gm22848 TSS18078  1:3783875-3783933  -  -  0  0  0  OK 
ENSMUSG00000025902  -  -  ENSMUSG00000025902  Sox17 TSS56047,TSS74369  1:4490927-4496413  -  -  0.611985  0.394887  0.829082  OK 
ENSMUSG00000051951  -  -  ENSMUSG00000051951  Xkr4 TSS1201,TSS70682,TSS88403  1:3205900-3671498  -  -  0  0  0  OK 

正如你所看到的，它缺乏指示GTF的第二列基因產物的類型。反正有袖釦自動將其納入其輸出文件？似乎沒有一個簡單的命令，除非我錯過了它。請告知 -

Answer 1

我嘗試了幾種方法讓它直接運行袖釦但沒有任何工作（如果有人有解決方案，請告訴我）。不過，我找到了解決辦法。我做了以下事情： 1-運行所有袖釦命令後，將數據通過FPKM列合併到另一個文件中（我有一個爲此編寫的腳本 - 如果有人願意，請告訴我）。 2-下載從UCSC瀏覽器表函數表 3-添加使用以下到合併輸出的命令該列上面的輸出的

Comb_Iso = pd.read_csv('Combined.Isoforms.txt', sep='\t') # reads in combined cufflinks output 
fnc_library = pd.read_csv(args.library, sep='\t', usecols=[0,1]) #only takes first 2 columns #reads in the path to the function library 
Comb_Iso_fnc = Comb_Iso.merge(fnc_library, on='TranscriptID', how='left') #merges them 
ncol_CIF = len(Comb_Iso_fnc.columns) 
df_iso=Comb_Iso_fnc.iloc[:,[0, 1, ncol_CIF-1] + list(range(2,ncol_CIF-1))] 
df_iso.to_csv(args.output_combined.strip()+'.Isoforms.txt', index=False, sep="\t") #generates the output. 

實施例目前是：

 
TranscriptID GeneID Function GeneName TSS-ID Locus-ID ExampleHeader1 ExampleHeader2 ExampleHeader3 
ENSMUST00000115488 ENSMUSG00000088000 protein_coding Gm25493 TSS84820 1:4723276-4723379 9 9 9 
ENSMUST00000027042 ENSMUSG00000096126 protein_coding Gm22307 TSS39629 1:4529016-4529123 8 8 8 
ENSMUST00000140302 ENSMUSG00000025902 Test2 Sox17 TSS56047,TSS74369 1:4490927-4496413 5 5 5 
ENSMUST00000115484 ENSMUSG00000051951 protein_coding Xkr4 TSS1201,TSS70682,TSS88403 1:3205900-3671498 6 6 6 
ENSMUST00000135046 ENSMUSG00000089699 protein_coding Gm1992 TSS23899 1:3466586-3513553 7 -7 -7 
ENSMUST00000132064 ENSMUSG00000025900 linc Rp1 TSS11475 1:4343506-4360314 3 -3 -3 
ENSMUST00000072177 ENSMUSG00000090025  Gm16088 TSS82763 1:3054232-3054733 1 1 1 
ENSMUST00000082125 ENSMUSG00000064842 protein_coding Gm26206 TSS81070 1:3102015-3102125 2 2 2 

可以看到第3欄中的功能。 很高興能分享我的任何腳本來合併袖釦（儘管這些可能不是最好的，仍然在優化工作），如果你想要的話。