0

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

文件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 行),则重叠最多的“中心行”(将是第 6 行,因为重叠为 200,而第 5 行的重叠为 150)
  3. 如果有 3 个或超过 3 个匹配项(例如,输出的第 3 行完全重叠,但第 2 行和第 4 行是具有部分重叠的相邻行,分别为 150 和 50),那么我只想返回中心行,这将在这种情况下是第 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

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

任何帮助都感激不尽。谢谢你。

4

1 回答 1

1

这不是完整的解决方案,因为我无法理解对时间预算的额外要求,但也许这会让你开始。

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

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

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,$2,$3]++-1)'
chr1 250 650 TF2 G2 200 400 E1 150
chr1 450 850 TF3 G3 400 600 E2 150

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

在您的输出中,您chr1 250 650列出了三遍,也许这是一个错字,或者我完全误解了您在这里尝试做的事情......

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

$ join fileB fileA | ... | awk '{print a[$1,$2,$3]++, $0}' | 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=($3<$7?$3:$7)-($2>$6?$2:$6)} 
  diff>0 {print $0,diff}'                   | 
sort -k1,1 -k9nr                            | 
awk '{print a[$1,$2,$3]++, $0}'             | 
sort -k1n                                   | 
awk '{file=($1==0)?"central":"neighbor"$1; 
      print $2,$3,$4,$5,$6,$7,$8,$9 > 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脚本中,但我认为以这种形式更容易理解(并在需要时更新)。

于 2017-05-12T14:56:14.350 回答