R、awk、sed:合并 bin 并输出中心重叠,然后是中心 + 相邻重叠

R, awk, sed: merge bins and output central overlaps and then central + neighboring overlaps

我有两个文件,我想重叠一些范围,并根据完全匹配和部分匹配检索结果。举个例子就清楚了。

文件A:

chr1    200     400     E1
chr1    400     600     E2
chr1    600     800     E3
chr2    200     300     E4

文件 B:

chr1    100     250   TF1   G1
chr1    250     650   TF2   G2
chr1    450     850   TF3   G3

输出:

chr1    100 250 TF1 G1  chr1    200     400     E1
chr1    250 650 TF2 G2  chr1    200     400     E1
chr1    250 650 TF2 G2  chr1    400     600     E2
chr1    250 650 TF2 G2  chr1    600     800     E3
chr1    450 850 TF3 G3  chr1    400     600     E2
chr1    450 850 TF3 G3  chr1    600     800     E3

到这一步我可以做一些事情,但下一步我需要你的帮助。

在这里我想首先对这些行进行子集化

  1. 只有 1 个匹配项(例如输出文件的第 1 行,与重叠大小无关)
  2. 如果有两个匹配项(例如输出的第 5 行和第 6 行),则 'central row' 重叠最多(第 6 行重叠为 200,而第 5 行重叠为 150 )
  3. 如果有 3 个或超过 3 个匹配项(例如,输出的第 3 行完全重叠,但第 2 行和第 4 行是部分重叠的相邻行,分别为 150 和 50)那么我想 return 只有中心行,在这种情况下将是第 3 行。

稍后,我想检索第一个邻居,然后是第二个邻居,依此类推,因为在实际数据集中,文件 B 中的一个 bin 可能会与文件 A 中最多 5 或 7 个 bin 重叠。

所以,基本上我想要的是首先获得所有中心重叠,然后是中心+第一个邻居,然后是中心+第二个邻居,依此类推。

根据这个原理,我的第一个结果将是:

结果 1(中央重叠):

chr1    100 250 TF1 G1  chr1    200     400     E1
chr1    250 650 TF2 G2  chr1    400     600     E2
chr1    450 850 TF3 G3  chr1    600     800     E3

结果 2(中央 + 第一个邻居):

chr1    100 250 TF1 G1  chr1    200     400     E1
chr1    250 650 TF2 G2  chr1    200     400     E1
chr1    250 650 TF2 G2  chr1    400     600     E2
chr1    250 650 TF2 G2  chr1    600     800     E3
chr1    450 850 TF3 G3  chr1    400     600     E2
chr1    450 850 TF3 G3  chr1    600     800     E3

如果可能的话,我想单独检索相邻的行而不是中间的行。

任何帮助将不胜感激。谢谢。

这不是完整的解决方案,因为我无法理解对我的时间预算的额外要求,但也许这可以帮助您入门。

假设文件按第一个键排序...

join fileB fileA | 
awk '{diff=(<?:)-(>?:)} diff>0{print [=10=], diff}' | 
sort -k1,1 -k9nr | 
awk '!a[,,]++'

chr1 250 650 TF2 G2 400 600 E2 200
chr1 450 850 TF3 G3 600 800 E3 200
chr1 100 250 TF1 G1 200 400 E1 50

最后一列显示重叠量,也许对后续步骤也有用。

更新

稍微修改最后一个awk你也可以获得第二个和第三个邻居

$ join fileB fileA | ...| awk '!(a[,,]++-1)'
chr1 250 650 TF2 G2 200 400 E1 150
chr1 450 850 TF3 G3 400 600 E2 150

$ join fileB fileA | ... | awk '!(a[,,]++-2)'
chr1 250 650 TF2 G2 600 800 E3 50

在你的输出中你 chr1 250 650 列出了 3 次,也许这是一个错字,或者我完全误解了你在这里想做什么...

或者,您可以在记录上标记顺序并根据顺序进行过滤。

$ join fileB fileA | ... | awk '{print a[,,]++, [=12=]}' | sort -k1n

0 chr1 100 250 TF1 G1 200 400 E1 50
0 chr1 250 650 TF2 G2 400 600 E2 200
0 chr1 450 850 TF3 G3 600 800 E3 200
1 chr1 250 650 TF2 G2 200 400 E1 150
1 chr1 450 850 TF3 G3 400 600 E2 150
2 chr1 250 650 TF2 G2 600 800 E3 50

此处第一列表示邻居编号,其中0是中心。

将所有内容拉到一起,您可以将所需的字段提取到单独的文件中

join fileB fileA                            | 
awk '    {diff=(<?:)-(>?:)} 
  diff>0 {print [=13=],diff}'                   | 
sort -k1,1 -k9nr                            | 
awk '{print a[,,]++, [=13=]}'             | 
sort -k1n                                   | 
awk '{file=(==0)?"central":"neighbor"; 
      print ,,,,,,, > file}'

创建这些文件。

==> central <==
chr1 100 250 TF1 G1 200 400 E1
chr1 250 650 TF2 G2 400 600 E2
chr1 450 850 TF3 G3 600 800 E3

==> neighbor1 <==
chr1 250 650 TF2 G2 200 400 E1
chr1 450 850 TF3 G3 400 600 E2

==> neighbor2 <==
chr1 250 650 TF2 G2 600 800 E3

请注意,所有这些都可以合并到一个 awk 脚本中,但我认为以这种形式更容易理解(并在需要时进行更新)。