我有 2 个文件，文件 1 的第 1 列必须替换为文件 2 的第 2 列，文件 1 的第 2,3,4-5 或 5-4 列（交叉匹配）与第 1,4,5 列匹配文件 2 的 -6 或 6-5。
文件1
SNP Chr Pos EA NEA EAF Beta SE Pvalue Neff
1:79137 1 79137 A T 0.25 -0.026 0.0073 4.0e-04 231420
1:79033 1 79033 A G 0.0047 -0.038 0.056 4.9e-01 225429
1:118630 1 118630 C T 0.99 -0.033 0.055 5.5e-01 226311
1:533179 1 533179 A G 1 -0.098 0.19 6.1e-01 185906
文件2
1 1:79033_A_G 0 79033 A G
1 1:79137_A_T 0 79137 T A
1 1:118630_C_T 0 118630 T C
1 1:533179_A_G 0 533179 G A
我需要输出如下所示：
SNP Chr Pos EA NEA EAF Beta SE Pvalue Neff
1:79137_A_T 1 79137 A T 0.25 -0.026 0.0073 4.0e-04 231420
1:79033_A_G 1 79033 A G 0.0047 -0.038 0.056 4.9e-01 225429
1:118630_C_T 1 118630 C T 0.99 -0.033 0.055 5.5e-01 226311
1:533179_A_G 1 533179 A G 1 -0.098 0.19 6.1e-01 185906
这些文件没有确切的行数，并且文件不是制表符分隔的。我尝试了下面的代码，但它不起作用，你能更正我的代码吗？
awk 'NR==FNR{chr[$1]=$1;snp[$2]=$2;pos[$4]=$4;a1[$5]=$5;a2[$6]=$6;next} ($1 in chr)&&($4 in pos)&& ((($5 in a1) && ($6 in a2)) || (($6 in a1) && ($5 in a2))) {$2==snp[$2]}' file 2 file1
编辑1：
下面的 Perl 代码会出现一些错误，并产生大约 20 000 行的重复行，一个例子是，
文件1
SNP Chr Pos EA NEA EAF Beta SE Pvalue Neff
7:10100610 7 10100610 A G 0.0002 0.13 0.58 8.2e-01 120658
7:10100610 7 10100610 C G 0.0013 0.1 0.13 4.4e-01 139170
10:1006107 10 1006107 C G 1 -0.11 0.42 7.9e-01 152016
文件 2
7 7:10100610_G_A 0 10100610 A G
7 7:10100610_G_C 0 10100610 C G
10 10:1006107_C_G 0 1006107 G C
这些行的预期输出：
7:10100610_G_A 7 10100610 A G 0.0002 0.13 0.58 8.2e-01 120658
7:10100610_G_C 7 10100610 C G 0.0013 0.1 0.13 4.4e-01 139170
10:1006107_C_G 10 1006107 C G 1 -0.11 0.42 7.9e-01 152016
但是 perl 代码给出的输出
7:10100610_G_A 7 10100610 A G 0.0002 0.13 0.58 8.2e-01 120658
10:1006107_C_G 7 10100610 C G 0.0013 0.1 0.13 4.4e-01 139170
10:1006107_C_G 10 1006107 C G 1 -0.11 0.42 7.9e-01 152016
答案1
该join命令将完成连接多个文件中的匹配行的工作。但它对其输入文件有一些要求，因此您需要在此过程中创建一些临时文件，以及一些额外的字段。
awk '{printf $2" "$3" "$4" "$5"%"$1"%"; $1="";print $0 "%" NR }' < file1 | sort > 1.tmp
awk '{print $1" "$4" "$5" "$6"%"$2} $5 != $6 {print $1" "$4" "$6" "$5"%"$2}' < file2 | sort > 2.tmp
join -a 1 -t % -o 1.4 2.2 1.2 1.3 1.tmp 2.tmp | sort -t % -n | awk -F % '!$2{$2=$3}{print $2" "$4}'
一步步
预处理第一个文件：
awk '{printf $2" "$3" "$4" "$5"%"$1"%"; $1="";print $0 "%" NR }''
输出示例：
1 118630 C T%1:118630% 1 118630 C T 0.99 -0.033 0.055 5.5e-01 226311%4
这 4 个字段由 分隔%，分别是：
- 必须匹配的“键”（输入字段 2-5）
- 原始第一列（如果没有匹配则需要）
- 原行的剩余部分
- 原始行号（这样我们就可以恢复之后的文件顺序
sort）
该输出通过管道传输sort到临时文件中，因为join需要对其输入进行排序。
对于第二个文件：
awk '{print $1" "$4" "$5" "$6"%"$2} $5 != $6 {print $1" "$4" "$6" "$5"%"$2}'
输出示例：
1 118630 C T%1:118630_C_T
1 118630 T C%1:118630_C_T
当您指定字段 5 和 6 应该匹配时，将打印第二行并将它们交换（前提是它们不相同）。这里的 - 分隔字段%是
- 需要匹配的“key”
- 第2栏
同样，输出通过管道传输sort到另一个临时文件中。
然后是主要的“加入”步骤：
join -a 1 -t % -o 1.4 2.2 1.2 1.3 1.tmp 2.tmp
当第二组中没有匹配项时， -a 1指示保留第一组中的行。将分隔符设置为（而不是空格）。该参数产生以下四个输出字段：join-t %%-o
- 文件 1，第 4 列：行号
- 文件2，第2列：替换自
file2（当没有匹配时，这将为空） - 文件 1，第 2 列：原始的第 1 列
file1 - 文件 1，第 3 列：该行的其余部分来自
file1
输出行示例：
4%1:118630_C_T%1:118630% 1 118630 C T 0.99 -0.033 0.055 5.5e-01 226311
然后sort可以恢复原来的文件顺序（按数字排序，字段分隔符%）
sort -t % -n
最后awk检查“替换”字段是否为空（因为未找到匹配项），如果是，则使用原始的column1。它还会丢弃行号和所有这些%s。
awk -F % '!$2{$2=$3}{print $2" "$4}'
最终输出线：
1:118630_C_T 1 118630 C T 0.99 -0.033 0.055 5.5e-01 226311
答案2
我会在 Perl 中这样做，因为它有一个sort函数可以让我们轻松地将A T和T A视为同一件事。例如（显示所有示例文件的输出组合）：
$ perl -lane 'if(!$k){$name{join("","chr".$F[0],$F[3],sort($F[4],$F[5]))}=$F[1]; }else{$var=join("", "chr".$F[1],$F[2],sort($F[3],$F[4])); $F[0]=$name{$var} if $name{$var};print join "\t", @F; } $k++ if eof' file2 file1
SNP Chr Pos EA NEA EAF Beta SE Pvalue Neff
7:10100610_G_A 7 10100610 A G 0.0002 0.13 0.58 8.2e-01 120658
7:10100610_G_C 7 10100610 C G 0.0013 0.1 0.13 4.4e-01 139170
10:1006107_C_G 10 1006107 C G 1 -0.11 0.42 7.9e-01 152016
1:79137_A_T 1 79137 A T 0.25 -0.026 0.0073 4.0e-04 231420
1:79033_A_G 1 79033 A G 0.0047 -0.038 0.056 4.9e-01 225429
1:118630_C_T 1 118630 C T 0.99 -0.033 0.055 5.5e-01 226311
1:533179_A_G 1 533179 A G 1 -0.098 0.19 6.1e-01 185906
或者，稍微更清晰一些：
$ perl -lane 'if(!$k){
$name{join("","chr".$F[0],$F[3],sort($F[4],$F[5]))}=$F[1];
}
else{
$var=join("", "chr".$F[1],$F[2],sort($F[3],$F[4]));
$F[0]=$name{$var} if $name{$var};
print join "\t", @F;
}
$k++ if eof' file2 file1
SNP Chr Pos EA NEA EAF Beta SE Pvalue Neff
7:10100610_G_A 7 10100610 A G 0.0002 0.13 0.58 8.2e-01 120658
7:10100610_G_C 7 10100610 C G 0.0013 0.1 0.13 4.4e-01 139170
10:1006107_C_G 10 1006107 C G 1 -0.11 0.42 7.9e-01 152016
1:79137_A_T 1 79137 A T 0.25 -0.026 0.0073 4.0e-04 231420
1:79033_A_G 1 79033 A G 0.0047 -0.038 0.056 4.9e-01 225429
解释
perl -lane：这-a使得 perl 的行为就像 awk 一样，自动将其输入分割到空格上的数组中@F。由于 Perl 数组从 开始0，$F[0]将是第一个字段，$F[1]将是第二个字段，依此类推。字段 N 是$F[N-1]。使-nperl 将其参数作为文本文件读取，并将 给出的脚本应用于-e其中的每一行。只是-l从每个输入行中删除尾随换行符，并向每个print调用添加一个换行符。$k++ if eof$k：如果到达文件末尾 (eof)，变量就会加1。然后我们可以使用if(!$k)（如果未定义 $k）作为NR==FNRawk 中的等效项。if(!$k){$name{join("","chr".$F[0],$F[3],sort($F[4],$F[5]))}=$F [1];}
: if this is the first file,file2, join fields 1, 4, and the sorted fields 5 and 6, into a string and use that string as the key in the hash (associative array)name. Then, save the variant's name from file2 as the value associated with that key. The sorting lets us treatA TandT Aas equivalent. I use"chr".$F[0]to deal with cases like1 123and11 23`，其中染色体和位置的串联给出相同的数字，即使染色体实际上不同。else{：如果我们现在正在读取第二个文件，file1.$var=join("", $F[1],$F[2],sort($F[3],$F[4]));: 构建密钥。这次使用字段 2、3，并对 4 和 5 进行排序。$F[0]=$name{$var} if $name{$var};name：如果该键有值，则将第一个字段设置为存储在哈希中的值。需要if确保我们不会更改标头或可能存在file1但不存在的任何其他变体file2。print join "\t", @F;：打印字段，包括上面刚刚所做的更改。