2017-06-06 32 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出现在数据您遇到的问题,但相应地选择它们。

相关问题