如何用Bash循环处理不同物种不同样本量的BAM文件并执行bcftools命令？

你可以通过Bash循环实现按物种分组处理样本，核心思路是先提取所有唯一的物种名称，再针对每个物种收集对应的BAM文件并执行bcftools mpileup命令。

完整脚本代码 #

# 获取所有唯一的物种标识（从BAM文件名的下划线前缀提取）
species_ids=$(ls *.bam | sed 's/_.*\.bam$//' | sort -u)

# 遍历每个物种
for sp in $species_ids; do
    # 匹配当前物种的所有BAM文件
    bam_list="${sp}_*.bam"
    # 定义输出BCF文件名（用物种名命名，便于区分）
    output_file="${sp}.bcf"
    
    # 执行bcftools命令
    bcftools mpileup -Ob -o "$output_file" -f ref.fa $bam_list
    
    # 可选：打印处理完成提示
    echo "已完成物种 ${sp} 的分析，输出文件：${output_file}"
done

关键步骤说明 #

• 提取物种标识：ls *.bam | sed 's/_.*\.bam$//' | sort -u

  • ls *.bam：列出当前目录下所有BAM文件
  • sed 's/_.*\.bam$//'：通过正则匹配，去掉每个文件名中_之后的内容以及.bam后缀，只保留物种名
  • sort -u：对提取的物种名去重，得到唯一的物种列表

• 循环处理每个物种：

  • ${sp}_*.bam：通配符匹配当前物种的所有BAM文件
  • 输出文件用物种名命名（如107.bcf、M82.bcf），避免不同物种的输出文件混淆
  • 替换命令中的ref.fa为你实际的参考基因组文件路径

注意事项 #

• 确保当前目录下只有目标BAM文件，避免无关文件被误匹配
• 如果BAM文件名包含特殊字符（如空格），需要调整脚本的文件匹配方式（比如用find命令结合while read循环），但生物信息学场景下这类情况很少见

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