脚本分享—快速统计基因组组装结果

我们小编欢乐豆有个压箱底的 perl 脚本，由于编程语言"洁癖"，想要彻底抛弃 perl 语言转向 python，于是他使用 AI 辅助下进行了转换，由于脚本相对简单，转换竟然就成功了。中间发现四种碱基含量百分比和原脚本统计有出入，检查确认是序列大小写没有注意的原因，修改后就完美运行了，这里分享给大家！

安装python模块

# 使用pip安装
pip install biopython

查看脚本参数

python N50Stat.py  -h
usage: N50Stat.py [-h] -i INPUT_FILE [-o OUTPUT_FILE]

Calculate statistics for DNA sequences in FASTA format.

optional arguments:
  -h, --help            show this help message and exit
  -i INPUT_FILE, --input_file INPUT_FILE
                        Input file containing DNA sequences in FASTA format.
  -o OUTPUT_FILE, --output_file OUTPUT_FILE
                        Output file to store the statistics.

脚本获取方法

关注下方微信公众号【微因】，后台回复关键字【脚本】 （不含中括号哟），建议粘贴复制，避免出错，获取脚本与测试文件。

脚本输出结果

脚本输出结果如下：

代码解释说明

先来用 AI 对脚本进行下解释说明：

• 导入模块: argparse：用于解析命令行参数的模块。Bio 中的 SeqIO：Biopython 库的一部分，用于读取和写入生物学序列文件格式。
• 函数: calc_n50(seq_lengths, percentile)：计算给定序列长度列表和指定百分位数的 N50 长度。calc_median(arr)：计算给定列表的中位数。base_count(seq, counters)：计算序列中核苷酸碱基（A、T、G、C、N）的出现次数。write_statistics(stats, output_file)：将计算得到的统计信息写入输出文件。print_statistics(stats)：将计算得到的统计信息打印到控制台。calculate_statistics(file_path, output_file)：处理 FASTA 文件，计算各种统计信息，并可选择将其写入输出文件。
• 主要部分: 使用 argparse 模块处理命令行参数。调用 calculate_statistics 函数，并提供输入文件路径和输出文件路径作为参数。
• 执行: 如果脚本作为主模块运行（name == "main"），它会调用 main 函数，而后者又调用带有命令行参数的 calculate_statistics 函数。
• 命令行用法: 脚本至少需要指定输入文件（-i 或 --input_file）。可选地，可以提供一个输出文件（-o 或 --output_file）以保存统计信息。例如，要运行脚本：python script_name.py -i input.fasta -o output_statistics.txt此脚本计算各种统计信息，如总序列数、总碱基数、最小和最大序列长度、平均和中位数序列长度，以及 N25、N50、N75、N90、N95 长度。此外，它计算每个核苷酸碱基的百分比，以及(A + T)和(G + C)的组合百分比。结果可以打印到控制台或保存到输出文件。

怎么样，有没有用，要不要收藏或者用起来呀？ 如果在学习过程中有什么疑问欢迎大家留言讨论。