Python随机核酸序列GC含量计算bug修复及碱基统计字典实现问题

1. GC含量计算错误原因 #

• 核心问题：GCcontent变量定义在for循环内部，每次迭代都会被当前序列的GC值覆盖，最终打印的仅为最后一条（第2条）序列的GC含量，并非两条序列的总和。
• 若要统计所有序列的GC总含量，需要提前初始化全局累计变量，每次循环中累加当前序列的GC数量。

2. 单条序列碱基统计字典实现思路 #

直接通过字符串count方法统计4种碱基的数量，生成固定键的字典即可，不存在的碱基会自动统计为0，完全符合要求的输出格式。

修正后完整代码 #

import random

def randseq(abc, length):
    return "".join([random.choice(abc) for i in range(random.randint(1, length))])

N = 2
longest_seq = ""
shortest_seq = randseq("ATCG", 10)
total_GC = 0  # 初始化GC总含量累计变量

for i in range(N):
    print(f'Sequence {i + 1}:')
    seq = randseq("ATCG", 10)
    if len(seq) > len(longest_seq):
        longest_seq = seq
    if len(seq) < len(shortest_seq):
        shortest_seq = seq
    # 计算当前序列GC含量并累加到总含量
    current_G = seq.count("G")
    current_C = seq.count("C")
    total_GC += current_G + current_C
    # 生成碱基统计字典
    base_stats = {
        "A": seq.count("A"),
        "T": seq.count("T"),
        "C": current_C,
        "G": current_G
    }
    print(seq)
    print(f"Sequence {i+1}: {base_stats}")

print("所有序列总GC含量是:", total_GC)

运行效果示例 #

Sequence 1:
TCGGTG
Sequence 1: {'A': 0, 'T': 2, 'C': 1, 'G': 3}
Sequence 2:
GCATCGTCAA
Sequence 2: {'A': 3, 'T': 2, 'C': 3, 'G': 2}
所有序列总GC含量是: 9

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