按行列号对文件进行子集筛选

Question

我们希望从一个文本文件中按行和列进行子集筛选，其中行号和列号从另一个文件中读取。需要排除表头（第一行）和行名（第一列）。

inputFile.txt 是一个制表符分隔的文本文件。

header  62  9   3   54  6   1
25  1   2   3   4   5   6
96  1   1   1   1   0   1
72  3   3   3   3   3   3
18  0   1   0   1   1   0
82  1   0   0   0   0   1
77  1   0   1   0   1   1
15  7   7   7   7   7   7
82  0   0   1   1   1   0
37  0   1   0   0   1   0
18  0   1   0   0   1   0
53  0   0   1   0   0   0
57  1   1   1   1   1   1

subsetCols.txt是一个没有空格的逗号分隔一行数字文件。在实际数据中，我们有500K列，并需要对其进行子集筛选，只保留约10K。

1,4,6

subsetRows.txt 是一个逗号分隔的一行数字文件，实际数据中有 20K 行，需要取其中约 ~300 行。

1,3,7

目前使用 cut 和 awk 循环的解决方案（相关帖子：使用 awk 选择行）：

# define vars
fileInput=inputFile.txt
fileRows=subsetRows.txt
fileCols=subsetCols.txt
fileOutput=result.txt

# cut columns and awk rows
cut -f2- $fileInput | cut -f`cat $fileCols` | sed '1d' | awk -v s=`cat $fileRows` 'BEGIN{split(s, a, ","); for (i in a) b[a[i]]} NR in b' > $fileOutput

输出文件：result.txt

1   4   6
3   3   3
7   7   7

问题：
对于小文件，这个解决方案运行良好。但对于更大的文件，50K行和200K列，它需要太长时间，超过15分钟，仍在运行。我认为剪切列可以正常工作，选择行是较慢的部分。

有更好的方法吗？

实际输入文件信息：

# $fileInput:
#        Rows = 20127
#        Cols = 533633
#        Size = 31 GB
# $fileCols: 12000 comma separated col numbers
# $fileRows: 300 comma separated row numbers

有关该文件的更多信息：文件包含GWAS基因型数据。每一行代表一个样本（个体），每一列代表一个SNP。为了进行进一步的区域分析，我们需要对样本（行）和SNP（列）进行子集处理，以使数据作为其他统计软件（如r）的输入更易于管理（减小）。

系统：

$ uname -a
Linux nYYY-XXXX ZZZ Tue Dec 18 17:22:54 CST 2012 x86_64 x86_64 x86_64 GNU/Linux

---

更新：以下解决方案由@JamesBrown提供，该方案导致我系统中的列顺序混乱，因为我使用的是不同版本的awk，我的版本是：GNU Awk 3.1.7

---

- zx8754

Answer 1

即使在如果编程语言是国家，每种语言将代表哪个国家？中，他们说...

Awk：朝鲜。顽固地抵制变化，其用户似乎出于我们只能推测的原因而非常喜欢它。

...每当你看到自己使用sed、cut、grep、awk等命令时，请停下来对自己说：awk可以单独完成！

因此，在这种情况下，重点是提取行和列（调整它们以排除标题和第一列），然后只需缓冲输出以最终打印它。

awk -v cols="1 4 6" -v rows="1 3 7" '
    BEGIN{
       split(cols,c); for (i in c) col[c[i]]  # extract cols to print
       split(rows,r); for (i in r) row[r[i]]  # extract rows to print
    }
    (NR-1 in row){
       for (i=2;i<=NF;i++) 
              (i-1) in col && line=(line ? line OFS $i : $i); # pick columns
              print line; line=""                             # print them
    }' file

使用您的示例文件：

$ awk -v cols="1 4 6" -v rows="1 3 7" 'BEGIN{split(cols,c); for (i in c) col[c[i]]; split(rows,r); for (i in r) row[r[i]]} (NR-1 in row){for (i=2;i<=NF;i++) (i-1) in col && line=(line ? line OFS $i : $i); print line; line=""}' file
1 4 6
3 3 3
7 7 7

使用您的示例文件，并将输入作为变量，在逗号上进行拆分：

awk -v cols="$(<$fileCols)" -v rows="$(<$fileRows)" 'BEGIN{split(cols,c, /,/); for (i in c) col[c[i]]; split(rows,r, /,/); for (i in r) row[r[i]]} (NR-1 in row){for (i=2;i<=NF;i++) (i-1) in col && line=(line ? line OFS $i : $i); print line; line=""}' $fileInput

我相信这肯定会更快。例如，您可以查看Remove duplicates from text file based on second text file以获得一些比较awk和grep等性能的基准测试。

- fedorqui

---

@Braiam 那只是个玩笑...不确定为什么要删除它。 - fedorqui

Answer 2

在Gnu awk 4.0或更高版本中，列排序依赖于for和PROCINFO["sorted_in"]。行号和列号从文件中读取：

$ awk '
BEGIN {
    PROCINFO["sorted_in"]="@ind_num_asc";
}
FILENAME==ARGV[1] {                       # process rows file
    n=split($0,t,","); 
    for(i=1;i<=n;i++) r[t[i]]
} 
FILENAME==ARGV[2] {                       # process cols file
    m=split($0,t,","); 
    for(i=1;i<=m;i++) c[t[i]]
} 
FILENAME==ARGV[3] && ((FNR-1) in r) {     # process data file
    for(i in c) 
        printf "%s%s", $(i+1), (++j%m?OFS:ORS)
}' subsetRows.txt subsetCols.txt inputFile.txt   
1 4 6
3 3 3
7 7 7

将ARGV[3]处理块移动到1和2之前，并在其末尾添加next，可能会带来一些性能提升。

- James Brown

---

不是很在意。如果你要测试的话，我有点想听听性能如何。此外，我在几次会议上写了一些代码，所以很高兴听到它是否有效...质量控制之类的。 :D - James Brown

请注意，使用for (i in c)可能不会提供正确的顺序：为什么awk似乎随机排列数组？ - fedorqui

@JamesBrown 好棒！然后+1，干得好！这个解决方案必须比我的更快，因为它不需要检查所有列是否需要打印，而是直接挑选它们。 - fedorqui

顺便问一下，詹姆斯：PROCINFO["sorted_in"]="@ind_num_asc"; 不应该在 BEGIN 块中设置吗？我很好奇仅在读取第一个文件时设置它是否会影响其他文件。 - fedorqui

显然，该功能仅在较新的版本中可用，4.0 版本在某处被提到过。 - James Brown

Answer 3

不是为了贬低这两个优秀的回答。由于这个问题涉及到大量的数据，我将发布两个答案的组合来加快处理速度。

awk -v cols="$(<subsetCols.txt)" -v rows="$(<subsetRows.txt)" '
BEGIN {
   n = split(cols, c, /,/)
   split(rows, r, /,/)
   for (i in r)
      row[r[i]]
}
(NR-1) in row {
   for (i=1; i<=n; i++)
      printf "%s%s", $(c[i]+1), (i<n?OFS:ORS)
}' inputFile.txt

注意: 这也适用于旧版awk或非gnu awk。

- anubhava

---

哦，很酷，你使用split的返回数字来知道要检查的列数，因此在提取正确的列号时使用它。这非常聪明，干得好。 - fedorqui

Answer 4

为了优化@anubhava的解决方案，我们可以通过利用输入已经排序的事实，摆脱每行搜索10k个值的负担，以确定我们是否在正确的行上。

awk -v cols="$(<subsetCols.txt)" -v rows="$(<subsetRows.txt)" '
BEGIN {
   n = split(cols, c, /,/)
   split(rows, r, /,/)
   j=1;
}
(NR-1) == r[j] { 
   j++
   for (i=1; i<=n; i++)
      printf "%s%s", $(c[i]+1), (i<n?OFS:ORS)
}' inputFile.txt

- tomc

---

我对awk不是很了解，无法确定这篇文章是否值得作为一个独立的答案，对我来说，这些更改似乎太小了，本可以作为@anubhava的评论？ - zx8754

也许你应该对其进行规模化的计时。 - tomc

Answer 5

Python有一个csv模块。你可以将一行读入列表，将所需的列打印到标准输出，反复执行。

这应该切片20,000到30,000列。

import csv
with open('foo.txt') as f:
    gwas = csv.reader(f, delimiter=',', quoting=csv.QUOTE_NONE)
    for row in gwas:
        print(row[20001:30001]

- RonJohn