在unix中,当行名与s不匹配时,将基于一列的两个文件进行内部联接

2024-09-30 01:34:47 发布

您现在位置:Python中文网/ 问答频道 /正文

我在尝试在unix中连接两个数据集时遇到了一个问题,需要您的帮助。我花了很长时间在论坛上寻找解决办法,但却一无所获。在

我在一个数据集中有一个登录号的列表,需要将它们转换成基因符号。为了这样做,我下载了gene2加入.gzfrom NCBI。未压缩的文件大约为7Gb,所以首先我从这个数据集中删除了加入和基因符号

cut -f 2,16 gene2accession > accession2genesymbol

根据wc -l accession2genesymbol有大约7000万行有很多重复的行,所以我用sort accession2genesymbol | uniq删除了这些行,结果得到了~2000万行。在

现在通常我会在R中执行^{} using the ^{} package(返回y中有匹配值的x中的所有行,以及x和y中的所有列);但是,这个数据集太大了,我无法处理。在

下面是一个未排序accession2genesymbol数据集的示例:

^{pr2}$

未排序的访问的一个简短示例如下(对于整个数据集-1576行see the gist):

Accessions
10047140
100913206
10092617
10190704
10190704
103471987
103471997
103472005
103472005
105990514
45006951
45006986
45006986
45007007
45007007
4501883
4501887
94721250
94721261
9558733
9845516
98986457
98986457
98986464
99028871
9910242
9951915
9966805
9966827
9966867
9994185

只有<1100STRONG匹配的<1100STRONG

9951915         LOAG_14435
110047140       LOC110047140

由于没有太多地使用unix进行数据操作/连接,我搜索堆栈溢出以查找类似的问题。在

例如this one。据我所知,unix join函数只能在文件排序的情况下使用,因此我尝试了以下操作:

join -t "\t" <(dos2unix <accession) <(dos2unix <accession2genesymbol.txt)

也许这是不起作用的,因为我需要两个数据集中相同的行号(即,如果数据集的第100行与数据集2的第100行不匹配,那么它就不起作用了),但也许我错了,因为其他原因,这个行不通?在

也许awk是一个更好的解决方案,所以我尝试了this post的建议:

awk '{a[$1]=a[$1] FS $2} END {for (i in a) print i a[i]}' accession accession2genesymbol | sort > file3

这将生成一个包含大约2000万行的文件,由于我的加入只有9000行,所以我希望有9000行(如果这些访问不再存在,可能会更少)。在

我在第一篇文章中尝试了另一个awk解决方案:

awk -F, 'FNR==NR{a[$1];next}($1 in a){print $2}' accession accession2genesymbol > file3
awk: warning: escape sequence `\+' treated as plain `+'

但我最终得到了一个空文件。在

我希望有一个awk(ward)解决方案,python,或者任何可以帮助我解决这个问题的东西。非常感谢你。在


Tags: 文件the数据示例排序基因符号unix
2条回答

join应该适合您的情况。由于您的输入文件没有匹配项,这里是一个虚构的示例,并使用您的映射文件

$ head file
100000009
100000061
100000030

$ join <(sed 1d map) <(sort file)
100000009 sema5bb+
100000030 btr24+
100000061 si:ch211-133n4.9+

假设您的map文件已经排序,您需要删除头sed 1d,并需要对输入file进行排序。注意排序应该是数字排序或词法排序。在

另一种不需要排序的方法是使用grep

^{pr2}$

如果数字和代码的格式不同,就不会有错误匹配。在

我们还没有看到原始的gene2accession文件的示例,但是假设它是一个制表符分隔的字段,第2列是accession,第16列是gene(因为这是您的cut选择的)。我们还假设您的Accessions文件不是绝对巨大的。在

鉴于此,此脚本应该执行您想要的操作:

awk -F'\t' 'NR==FNR{a[$1];next} ($2 in a) && !seen[$2]++{print $2, $16}' Accessions gene2accession

但你可以试试看是否更快:

^{pr2}$

如果是,并且您希望在后续运行中使用sort输出的中间文件:

sort -u -t'\t' -k2,2 gene2accession > unq_gene2accession &&
awk -F'\t' 'NR==FNR{a[$1];next} $2 in a{print $2, $16}' Accessions unq_gene2accession

相关问题 更多 >

    热门问题