2017-06-06 48 views
0

當基於三個單獨列合併文件時,我遇到一個小問題。首先,我的代碼和文件結構以及下面有關我的問題的更多詳細信息。這裏是我的設置至今:通過匹配三列合併兩個文件+計算輸出的均值

#making directories 
subprocess.call(""" mkdir %s/junctions """%(temp_dir), shell=True) 
subprocess.call(""" mkdir %s/out """%(temp_dir), shell=True) 

#opening a file with paths to other files 
with open(sys.argv[2],"r") as j: 
    #creating the output file 
    subprocess.call(""" touch %s/junctions/catjunc.txt """%(temp_dir), shell=True) 

    for line in j: 
     #reformatting the input file (not important for this question) 
     command = """awk 'BEGIN{OFS="\t"}{print $1, $2-20-1, $3+20, "JUNCBJ"NR"%s", $7, ($4 == 1)? "+":"-",$2-20-1, $3+20, "255,0,0", 2, "20,20", "0,300", $7, $8 ,$5 , $6}' %s > %s"""%(line[:-1].split(".")[0].split("/")[-1],line[:-1],temp_dir + "/junctions/junc.bed") 
     subprocess.call(command, shell=True) 

     # So here i basically concatenate the files. However, I also want to reduce them. 
     subprocess.call(""" cat %s >> %s/junctions/catjunc.txt """%(temp_dir + "/junctions/junc.bed", temp_dir), shell=True) 

的文件:

chrom start stop ID                        count1 count2 
1  14809 14989 JUNCBJ1adipose_HS110_50bp_SJ 0  -  14809 14989 255,0,0 2  20,20 0,300 0  59  2  1 
1  14809 15815 JUNCBJ2adipose_HS110_50bp_SJ 0  -  14809 15815 255,0,0 2  20,20 0,300 0  2  2  1 
1  15018 15815 JUNCBJ4adipose_HS110_50bp_SJ 0  -  15018 15815 255,0,0 2  20,20 0,300 0  76  2  1 
1  15927 16626 JUNCBJ5adipose_HS110_50bp_SJ 0  -  15927 16626 255,0,0 2  20,20 0,300 0  4  2  1 
1  16745 16873 JUNCBJ6adipose_HS110_50bp_SJ 0  -  16745 16873 255,0,0 2  20,20 0,300 0  2  2  1 

我所有的文件看起來是這樣的。第1,2和3列是染色體和起始和終止座標。在第12和13欄中有一些重要的內容。第4列是行ID。 現在我的問題。在連接兩個文件時,我想繼續使用所有的獨特行。但是,如果有兩行,其中三個第一列全部相同,那麼我只想包括其中兩個計數列(12和13)的計數值被均值替換的兩個行中的一個。我也想追加兩個ID(現在不太重要)。

實施例:

file1的:

1  14809 14989 JUNCBJ1adipose_HS110_50bp_SJ 0  -  14809 14989 255,0,0 2  20,20 0,300 10  59  2  1 
1  14809 15815 JUNCBJ2adipose_HS110_50bp_SJ 0  -  14809 15815 255,0,0 2  20,20 0,300 0  2  2  1 

file2的:

1  14809 14989 JUNCBG2adipose_HS110_50bp_SJ 0  -  14809 14989 255,0,0 2  20,20 0,300 20  41  2  1 

這裏,第一行出現在這兩個文件用計數10和59以及20和41在輸出文件的行會出現一次,計數爲15和50.第二行是簡單的複製。

輸出:

                                     replaced means 
1  14809 14989 JUNCBJ1adipose_HS110_50bp_SJ,JUNCBG2adipose_HS110_50bp_SJ 0  -  14809 14989 255,0,0 2  20,20 0,300 15  50  2  1 
1  14809 15815 JUNCBJ2adipose_HS110_50bp_SJ 0  -  14809 15815 255,0,0 2  20,20 0,300 0  2  2  1 

如果有我錯過任何細節,我很高興地對帖子進行編輯。請不要判斷我的python-bash混合,並且

乾杯。

回答

0

使用awk。文件中的字段分隔符是空格,因此值10和59位於字段13和14(awk中的$13$14)。如果字段分隔符是製表符,而是在第一個{(未經測試但受過教育的猜測,其詞語是演示效果的邀請)之前添加BEGIN{FS=OFS="\t"},並用13(以該順序;)替換全部13與12和14。

$ awk ' 
{ 
    a13[$1,$2,$3]+=$13    # sum field 13 
    a14[$1,$2,$3]+=$14    # sum field 14 
    $13="XX"       # set a replaceable key to $13 
    $14="YY"       # see above 
    a0[$1,$2,$3]=$0     # now store the record 
    an[$1,$2,$3]++     # count instances of $1,$2,$3 
} 
END {         # after all the records 
    for(i in a0) {     # loop all stored records 
     sub(/XX/,a13[i]/an[i],a0[i]) # replace XXs with means 
     sub(/YY/,a14[i]/an[i],a0[i]) # and YY 
     print a0[i] }     # output 
}' file1 file2 
1 14809 15815 JUNCBJ2adipose_HS110_50bp_SJ 0 - 14809 15815 255,0,0 2 20,20 0,300 0 2 2 1 
1 14809 14989 JUNCBG2adipose_HS110_50bp_SJ 0 - 14809 14989 255,0,0 2 20,20 0,300 15 50 2 1 

當然,如果關鍵字XXYY出現在數據您遇到的問題,但相應地選擇它們。

相關問題