基于基因坐标文件提取FASTA序列指定区域时KeyError报错求助

这个错误的核心是你存储FASTA序列的字典里没有'MN908947.3'这个键——大概率是读取FASTA文件时，你把整条标题行（比如>MN908947.3 Severe acute respiratory syndrome...）当成了键，而不是只取空格前的ID部分（也就是'MN908947.3'）。

先修正FASTA读取逻辑 #

把原来的FASTA读取代码改成只提取标题行里第一个空格前的内容作为键：

def read_fasta(fasta_file):
    fasta_dict = {}
    current_id = None
    current_seq = []
    with open(fasta_file, 'r') as f:
        for line in f:
            line = line.strip()
            if not line:
                continue
            if line.startswith('>'):
                # 保存上一条序列
                if current_id is not None:
                    fasta_dict[current_id] = ''.join(current_seq)
                # 只取>后第一个空格前的内容作为ID
                current_id = line[1:].split()[0]
                current_seq = []
            else:
                current_seq.append(line)
        # 处理最后一条序列
        if current_id is not None:
            fasta_dict[current_id] = ''.join(current_seq)
    return fasta_dict

# 读取FASTA文件并验证键
fasta_dict = read_fasta('GENE.fasta.txt')
print(fasta_dict.keys())  # 看看输出里有没有'MN908947.3'

再处理坐标文件提取序列 #

注意GenBank的坐标是1-based，而Python字符串切片是0-based，所以起始位置要减1：

def extract_gene_seqs(coords_file, fasta_dict):
    with open(coords_file, 'r') as coords, open('extracted_seqs.fasta', 'w') as out:
        for line in coords:
            line = line.strip()
            if not line:
                continue
            gene_name, start, end = line.split()
            # 转换为0-based索引
            start_idx = int(start) - 1
            end_idx = int(end)
            # 检查序列ID是否存在
            if 'MN908947.3' not in fasta_dict:
                print("FASTA文件中未找到MN908947.3序列")
                continue
            # 提取序列片段
            target_seq = fasta_dict['MN908947.3'][start_idx:end_idx]
            # 写入结果
            out.write(f">{gene_name}\n{target_seq}\n")

# 执行提取
extract_gene_seqs('你的坐标文件.txt', fasta_dict)

调试小技巧 #

• 运行读取FASTA的代码后，查看打印的键列表，确认'MN908947.3'是否存在，如果不存在，检查FASTA文件的标题行格式有没有特殊字符。
• 确认坐标文件里的起始/终止位置是有效数字，没有额外的空格或符号。

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