使用 awk 通过文件中的 ID 从 multifasta 文件中提取序列

Posted 2023-02-21

【中文标题】使用 awk 通过文件中的 ID 从 multifasta 文件中提取序列

【英文标题】：extract sequences from multifasta file by ID in file using awk

【发布时间】：2018-09-18 18:14:36

【问题描述】：

我想从 multifasta 文件中提取与由单独的 ID 列表给出的 ID 匹配的序列。

FASTA 文件 seq.fasta：

>7P58X:01332:11636
TTCAGCAAGCCGAGTCCTGCGTCGTTACTTCGCTT
CAAGTCCCTGTTCGGGCGCC
>7P58X:01334:11605
TTCAGCAAGCCGAGTCCTGCGTCGAGAGTTCAAGTC
CCTGTTCGGGCGCCACTGCTAG
>7P58X:01334:11613
ACGAGTGCGTCAGACCCTTTTAGTCAGTGTGGAAAC
>7P58X:01334:11635
TTCAGCAAGCCGAGTCCTGCGTCGAGAGATCGCTTT
CAAGTCCCTGTTCGGGCGCCACTGCGGGTCTGTGTC
GAGCG
>7P58X:01336:11621
ACGCTCGACACAGACCTTTAGTCAGTGTGGAAATCT
CTAGCAGTAGAGGAGATCTCCTCGACGCAGGACT

ID 文件 id.txt：

7P58X:01332:11636
7P58X:01334:11613

我想获取仅包含与 id.txt 文件中的 ID 匹配的序列的 fasta 文件：

>7P58X:01332:11636
TTCAGCAAGCCGAGTCCTGCGTCGTTACTTCGCTTT
CAAGTCCCTGTTCGGGCGCC
>7P58X:01334:11613
ACGAGTGCGTCAGACCCTTTTAGTCAGTGTGGAAAC

我真的很喜欢我在答案 here 和 here 中找到的 awk 方法，但是那里给出的代码对于我给出的示例仍然不能完美运行。原因如下：

(1)

awk -v seq="7P58X:01332:11636" -v RS='>' '$1 == seq print RS $0' seq.fasta

此代码适用于多行序列，但 ID 必须单独插入代码中。

(2)

awk 'NR==FNRn[">"$0];next fprint f ORS $0;f="" $0 in nf=$0' id.txt seq.fasta

此代码可以从 id.txt 文件中获取 ID，但只返回多行序列的第一行。

我想最好的办法是修改代码 (2) 中的 RS 变量，但到目前为止我所有的尝试都失败了。请问有人可以帮我吗？

【问题讨论】：

awk 'NR==FNRids[$0];next /^>/f=($1 in ids) f' file.ids file.data ***.com/a/39898676/6260170

我会使用 bioawk，但我的方法“单独插入变量”，这可能不是最优的：for seq_id in $(cat id.txt); do bioawk -c fastx -v seq_id="$seq_id" '$name == seq_id print ">"$name"\n"$seq' seq.fasta; done

【参考方案1】：

关注awk 可能对您有所帮助。

awk 'FNR==NRa[$0];next /^>/val=$0;sub(/^>/,"",val);flag=val in a?1:0 flag' ids.txt  fasta_file

【参考方案2】：

$ awk -F'>' 'NR==FNRids[$0]; next NF>1f=($2 in ids) f' id.txt seq.fasta
>7P58X:01332:11636
TTCAGCAAGCCGAGTCCTGCGTCGTTACTTCGCTT
CAAGTCCCTGTTCGGGCGCC
>7P58X:01334:11613
ACGAGTGCGTCAGACCCTTTTAGTCAGTGTGGAAAC

【讨论】：

应该如何修改它以接受来自上一步的管道 seq.fasta 文件？

whatever | awk '...' id.txt -