如何从批量Prokka注释结果的ffn文件中自动提取TcdA基因

所有方案默认文件结构符合描述：所有菌株的独立文件夹存放在同一根目录下，每个菌株文件夹内的.ffn文件和文件夹同名（如Strain1文件夹对应Strain1.ffn），无需手动进入任何单菌株文件夹操作。

方案1：Shell脚本（Linux/macOS/WSL环境直接运行，无额外依赖） #

• 终端进入存放所有菌株文件夹的根目录
• 直接运行以下脚本，提取结果会统一汇总到all_TcdA_sequences.fasta文件中，每条序列自动附带对应菌株名前缀，方便后续溯源分析

# 清空/初始化输出文件
> all_TcdA_sequences.fasta

# 遍历所有子目录下的ffn注释文件
for ffn in */*.ffn; do
    strain=$(echo "$ffn" | cut -d'/' -f1)
    awk -v s_id="$strain" '
        BEGIN{keep=0}
        /^>/{
            keep=0
            if($0 ~ /glycosylating toxin TcdA/){
                keep=1
                print ">"s_id"|"substr($0,2)
            }
        }
        keep{print}
    ' "$ffn" >> all_TcdA_sequences.fasta
done

脚本逻辑说明：

• 自动识别所有子文件夹下的.ffn格式文件，跳过其他后缀的Prokka结果文件
• 仅提取序列头（>开头的注释行）包含目标标识的完整fasta条目，不会错提其他基因序列
• 自动给提取到的序列加菌株名前缀，避免不同菌株的同源序列头重复

方案2：Python脚本（全系统通用，适配无Shell环境的Windows系统） #

将以下脚本保存为extract_TcdA.py，放到所有菌株文件夹的根目录下，确保本地已安装Python3环境，直接运行python extract_TcdA.py即可，输出结果和Shell脚本完全一致：

import glob
import os

output_file = "all_TcdA_sequences.fasta"
with open(output_file, "w", encoding="utf-8") as outf:
    # 遍历所有子目录下的ffn文件
    for ffn_path in glob.glob("*/*.ffn"):
        strain_id = os.path.basename(os.path.dirname(ffn_path))
        keep_seq = False
        with open(ffn_path, "r", encoding="utf-8") as inf:
            for line in inf:
                if line.startswith(">"):
                    keep_seq = False
                    if "glycosylating toxin TcdA" in line:
                        keep_seq = True
                        outf.write(f">{strain_id}|{line[1:]}")
                elif keep_seq:
                    outf.write(line)

结果校验提示 #

提取完成后可快速校验提取数量：

• Shell环境下运行grep -c ">" all_TcdA_sequences.fasta，返回的数字即为提取到的TcdA序列总条数
• Windows环境下可直接打开输出的fasta文件，统计>开头的行数即可
  如果存在注释大小写不统一的情况，可修改脚本匹配逻辑为不区分大小写：Shell中awk匹配规则改为~ /glycosylating toxin TcdA/I，Python中判断条件改为"glycosylating toxin tcda" in line.lower()。

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