2 个文件之间的交集

2 个文件之间的交集

我有一个文件包含单核苷酸多态性数据称为snp.bed,看起来像这样:

head snp.bed

    Chr17   214708483   214708484   Chr17:214708484
    Chr17   214708507   214708508   Chr17:214708508
    Chr17   214708573   214708574   Chr17:214708574

我还有一个名为 的文件intersect.bed,如下所示:

head intersect.bed

    Chr17   214708483   214708484   Chr17:214708484 Chr17   214706266   214710783   gene50573
    Chr17   214708507   214708508   Chr17:214708508 Chr17   214706266   214710783   gene50573
    Chr17   214708587   214708588   Chr17:214708580 Chr17   214706266   214710783   gene50573

我想打印出修改后的版本,snp.bed其中包含附加到每行的额外列。如果 中 的一行与snp.bed中 一行的前 4 列匹配intersect.bed,那么我想打印整行,并通过与(基因名称)snp.bed中相应行的最后一列相邻而获得额外的列。intersect.bed或者,如果来自的行snp.bed与来自的任何行都不匹配,intersect.bed则与由字符串“NA”而不是基因名称组成的额外列相邻。

这是我想要的输出:

head snp.matched.bed

    Chr17   214708483   214708484   Chr17:214708484   gene50573
    Chr17   214708507   214708508   Chr17:214708508   gene50573
    Chr17   214708573   214708574   Chr17:214708574   NA

我怎样才能做到这一点?

答案1

该解决方案假设文件的行首没有空格。与您的示例有什么不同,其中有这些空格。

awk '
{
    str = $1$2$3$4; 
}
FNR == NR {
    arr[str] = $NF;
}
FNR != NR {
    gene_name = arr[str] ? arr[str] : "NA";
    print $0, gene_name;
}' intersect.bed snp.bed 

输出

Chr17   214708483   214708484   Chr17:214708484 gene50573
Chr17   214708507   214708508   Chr17:214708508 gene50573
Chr17   214708573   214708574   Chr17:214708574 NA

答案2

这是使用 awk 的解决方案:

$ awk -F '\t' 'BEGIN{while(getline line<"intersect.bed") {N=split(line,a,"\t"); seen[a[1]"\t"a[2]"\t"a[3]"\t"a[4]]=a[N];}} {if(seen[$0]) {name=seen[$0];} else{name="NA"}; print $0 "\t" name}' snp.bed
Chr17       214708483       214708484       Chr17:214708484 gene50573
Chr17       214708507       214708508       Chr17:214708508 gene50573
Chr17       214708573       214708574       Chr17:214708574 NA

我假设单个制表符作为两个输入文件的分隔符。

另请注意,我将“第一个第四列”解释为“前四列”。

答案3

就我个人而言,我认为对于此类任务,最好使用“真正的”编程语言。我喜欢 Python,所以这里有一个 Python 脚本可以完成您想要的操作(它故意冗长,以便您可以理解并轻松修改它):

#!/usr/bin/env python2

# intersect.py

# Read data from the first file
snp_rows = []
with open("snp.bed", 'r') as snp_file:
    for row in snp_file:
        snp_rows.append(row.split())

# Read data from the second file
int_rows = []
with open("intersect.bed", 'r') as int_file:
    for row in int_file:
        int_rows.append(row.split())

# Compare data and compute results
results = []
for row in int_rows:
    if row[:4] in snp_rows:
        results.append(row[:4] + [row[-1]])
    else:
        results.append(row[:4] + ["NA"])

# Print the results
for row in results:
    print(' '.join(row))

将其保存到文件中然后执行:

python2 intersect.py

只是为了好玩,这里是一个使用标准命令(只是grepcut)的 Bash 解决方案:

while read row; do
    match="$(grep -F "${row}" intersect.bed)";
    if [[ -n "${match}" ]]; then
        echo "${row} $(echo ${match} | cut -d' ' -f8)";
    else
        echo "${row} NA";
    fi;
done < snp.bed

请注意,一般不建议使用 Bash 进行严格的文本处理。例如,请参阅以下帖子:

相关内容