如何根据accession code筛选fasta文件中的蛋白序列？

先说说你之前方法失败的原因 #

• pcregrep/grep：默认只逐行匹配，FASTA序列是多行结构，它们只会输出匹配的头部行，而序列行本身不含accession，所以只能抓到头部和首行序列
• awk：大概率是你没处理FASTA的多行结构，或者没正确把accession列表读入内存

方法1：用awk（最可靠，推荐） #

awk可以轻松处理FASTA的多行结构，先把需要保留的accession存进数组，再遍历FASTA文件判断输出。

情况1：FASTA头部是> A0A2K5BUY4 序列描述（accession为第二个字段）

运行这条命令：

awk 'NR==FNR {acc[$1]=1; next} /^>/ {flag=($2 in acc)?1:0} flag' accession.txt sequences.fasta > filtered.fasta

• 命令拆解：
  1. NR==FNR：处理第一个输入文件accession.txt，把每个accession存入acc数组
  2. /^>/：遇到FASTA头部行（以>开头），检查第二个字段是否在acc数组里，设置输出标记flag
  3. flag：只要标记为1，就输出当前行（包括头部和所有后续的序列行，直到下一个头部重新判断标记）

情况2：FASTA头部是>sp|A0A2K5BUY4|序列描述（accession在|分隔的第二个位置）

调整提取accession的逻辑：

awk 'NR==FNR {acc[$1]=1; next} /^>/ {split($0, arr, "|"); acc_id=arr[2]; flag=(acc_id in acc)?1:0} flag' accession.txt sequences.fasta > filtered.fasta

• 用split把头部按|拆分成数组，取第二个元素作为accession再判断是否需要输出

方法2：用grep+sed（适合不想写复杂脚本的情况） #

先找出所有匹配头部的行号，再提取对应行到下一个头部之间的内容：

# 第一步：获取所有匹配头部的行号
grep -n -F -f accession.txt sequences.fasta | cut -d: -f1 > match_lines.txt

# 第二步：提取每个匹配行到下一个头部之前的内容
sed -n -f <(for line in $(cat match_lines.txt); do echo "$line,/^>/-1p"; done) sequences.fasta > filtered.fasta

• 命令拆解：
  1. grep -n -F -f：-n显示行号，-F固定字符串匹配（避免正则冲突），-f用accession.txt作为匹配列表
  2. sed部分：循环每个匹配行号，生成行号,/^>/-1p的命令，意思是从该行输出到下一个>行的前一行

注意事项 #

• 确保accession.txt里每个accession单独占一行，没有多余空格或换行
• 确认FASTA头部的accession和accession.txt里的完全一致（大小写、前缀后缀都要对应）
• 测试时可以先挑1-2个已知的accession，用小样本文件验证命令是否正确

内容的提问来源于stack exchange，提问作者Blumina Romero