1

我正在寻找一种有效的解决方案来过滤两个数据集。基本上，我只想保留一个文件中没有缺失的行作为它们的键列，并且在两个文件中都没有值“0/0”。

输入数据（对于那些感兴趣的人，我为这个问题简化的基因组 VCF 文件）具有以下特征：

  • 第 1 列和第 2 列一起是按数字排序的唯一标识符
  • 第 3 列以值 0/0、0/1 或 1/1 开头

理想情况下，该脚本将执行以下操作：

  1. 遍历 sample1.dat 中的每一行并在 sample2.dat 中查找相同的标识符
  2. 如果在 sample2.dat 中找不到来自 sample1.dat 的标识符，则什么也不做
  3. 如果两行都包含“0/0”，则什么也不做
  4. 如果一行或两行不包含“0/0”，则将这两行写入各自的输出。

输入

样本1.dat

1      1001  0/0:8:8:99:PASS
1      1002  0/0:8:8:99:PASS
1      1003  0/1:5,3:8:99:PASS,PASS
2      1234  0/0:8:8:99:PASS # not present in sample2
2      2345  1/1:8:8:99:PASS

样本2.dat

1      1001   0/0:8:8:99:PASS
1      1002   0/1:5,3:8:99:PASS,PASS
1      1003   0/0:8:8:99:PASS
2      2345   1/1:8:8:99:PASS
2      3456   0/1:8:8:99:PASS # not present in sample1

输出

sample1_out.dat

1      1002  0/0:8:8:99:PASS
1      1003  0/1:5,3:8:99:PASS,PASS
2      2345  1/1:8:8:99:PASS

sample2_out.dat

1      1002   0/1:5,3:8:99:PASS,PASS
1      1003   0/0:8:8:99:PASS
2      2345   1/1:8:8:99:PASS

在这种情况下，不会打印 1-1001，因为它们都具有值“0/0”，并且不会打印 2-1234 和 2-3456，因为它们在两个文件中都不存在。

一些注意事项：

  • 这些文件大约 260GB，但我可以轻松地将它们拆分为多个最大 18GB 的​​文件（我基本上将它们拆分为染色体）
  • 我机器上的可用内存约为 128GB
  • 第 1 列和第 2 列一起已经按数字顺序排序

任何帮助是极大的赞赏！

4

1 回答 1

1

awk救援！可能您需要先拆分文件，对于每个块，您都可以这样做

$ function f { awk -v OFS='\t' '{print $1"~"$2,$0}' $1; }; 
  join <(f file1) <(f file2) | 
  awk -v OFS='\t' '$4!~/0\/0/ || $7!~/0\/0/ 
                     {print $2,$3,$4 > "file1.out"; 
                      print $5,$6,$7 > "file2.out"}'

解释让我们join来做匹配对应记录的工作，但是需要通过合并前两个字段来创建一个合成键。输出包含我们需要的所有信息，将结果的相应字段和输出部分中的“0/0”逻辑应用到相应的输出文件。

$ head file{1,2}.out                          

==> file1.out <==                                                                                                     
1       1002    0/0:8:8:99:PASS
1       1003    0/1:5,3:8:99:PASS,PASS
2       2345    1/1:8:8:99:PASS

==> file2.out <==
1       1002    0/1:5,3:8:99:PASS,PASS
1       1003    0/0:8:8:99:PASS
2       2345    1/1:8:8:99:PASS

您可能会更好地参数化文件名，无论是输入还是输出。

于 2016-07-01T12:37:54.110 回答