sequence_files

# 处理序列文件 (FASTA/FASTQ) ## FASTA 文件 ### 概述 Pysam 提供了 `FastaFile` 类，用于对 FASTA 参考序列进行索引和随机访问。FASTA 文件在使用前必须通过 `samtools faidx` 进行索引。 ### 打开 FASTA 文件 python import pysam # 打开已索引的 FASTA 文件 fasta = pysam.FastaFile("reference.fasta") # 自动查找 reference.fasta.fai 索引 ### 创建 FASTA 索引 python # 使用 pysam 创建索引 pysam.faidx("reference.fasta") # 或使用 samtools 命令 pysam.samtools.faidx("reference.fasta") 这将创建随机访问所需的 `.fai` 索引文件。 ### FastaFile 属性 python fasta = pysam.FastaFile("reference.fasta") # 参考序列列表 references = fasta.references print(f"References: {references}") # 获取长度 lengths = fasta.lengths print(f"Lengths: {lengths}") # 获取特定序列长度 chr1_length = fasta.get_reference_length("chr1") ### 获取序列 #### 按区域获取 使用 **0-基、半开区间** 坐标： python # 获取特定区域 sequence = fasta.fetch("chr1", 1000, 2000) print(f"Sequence: {sequence}") # 返回 1000 个碱基 # 获取整个染色体 chr1_seq = fasta.fetch("chr1") # 使用区域字符串获取 (1-基) sequence = fasta.fetch(region="chr1:1001-2000") **重要：** 数值参数使用 0-基坐标，区域字符串使用 1-基坐标（samtools 约定）。 #### 常见用例 python # 获取变异位置的序列 def get_variant_context(fasta, chrom, pos, window=10): """获取变异位置（1-基）周围的序列上下文。""" start = max(0, pos - window - 1) # 转换为 0-基 end = pos + window return fasta.fetch(chrom, start, end) # 获取基因坐标的序列 def get_gene_sequence(fasta, chrom, start, end, strand): """获取具有链方向意识的基因序列。""" seq = fasta.fetch(chrom, start, end) if strand == "-": # 反向互补 complement = str.maketrans("ATGCatgc", "TACGtacg") seq = seq.translate(complement)[::-1] return seq # 检查参考等位基因 def check_ref_allele(fasta, chrom, pos, expected_ref): """验证位置（1-基 pos）处的参考等位基因。""" actual = fasta.fetch(chrom, pos-1, pos) # 转换为 0-基 return actual.upper() == expected_ref.upper() ### 提取多个区域 python # 高效提取多个区域 regions = [ ("chr1",