pysam 是一个 Python 模块，用于读取、操作和写入基因组比对格式（SAM/BAM/CRAM）以及变异文件（VCF/BCF），它是对 htslib 的封装。该工具主要用于高通量测序流程，能够高效访问比对到参考基因组的数十亿个 DNA 片段。
用途
pysam 适用于以下场景：处理下一代测序（NGS）数据；分析基因组变异（SNP、indel）；从特定基因组区域提取读段；构建自定义生物信息学流程；对测序数据进行质量控制。
核心原则
索引访问：BAM 文件必须建立索引（.bai），才能对基因组区域进行高效的随机访问。
坐标系统：基因组坐标中的位置采用 0 基（Python 风格），但在某些上下文中范围采用 1 基。
读段属性：每个读段包含序列、质量分数、比对位置和标志。
快速参考
标准导入：
import pysam
基本模式：
1. 打开 BAM 文件
samfile = pysam.AlignmentFile("aligned_reads.bam", "rb")
2. 遍历特定基因组区域内的读段
for read in samfile.fetch("chr1", 10000, 10100): print(f"Read: {read.query_name}, Quality: {read.mapping_quality}") print(f"Sequence: {read.query_sequence}") print(f"Position: {read.reference_start}")
3. 变异分析（VCF）
vcf = pysam.VariantFile("mutations.vcf") for rec in vcf.fetch("chr1", 10000, 10100): print(f"Pos: {rec.pos}, Ref: {rec.ref}, Alt: {rec.alts}") print(f"Genotype: {rec.samples['sample1']['GT']}")
4. 写入比对读段
outfile = pysam.AlignmentFile("output.bam", "wb", template=samfile) for read in samfile: if read.mapping_quality > 30: outfile.write(read) outfile.close()
关键规则
应当做的：始终使用带索引的文件——用 pysam.index("file.bam") 创建索引以实现快速访问；检查读段标志——使用 read.is_paired、read.is_unmapped 过滤读段；处理未比对读段——未比对读段的 reference_start = -1；显式关闭文件——使用上下文管理器或 .close() 避免资源泄漏。
不应做的：不要遍历整个 BAM——使用带区域的 fetch() 以提高效率；不要忽略质量分数——低质量碱基可能导致假变异；不要混用坐标系统——在 0 基与 1 基索引之间保持一致。
高级模式
按基因统计读段数：
使用基因注释文件
genes = {} # gene_name -> (chr, start, end) for read in samfile.fetch(): # 检查读段是否与任何基因重叠 for gene, (chr, start, end) in genes.items(): if read.reference_name == chr and start <= read.reference_start < end: genes[gene]['count'] += 1
变异过滤：
按质量和深度过滤变异
for rec in vcf.fetch(): depth = rec.samples['sample1']['DP'] qual = rec.qual if depth > 10 and qual > 20: # 处理高质量变异 pass
pysam 提供了基因组数据处理所需的底层访问能力，使研究人员能够直接处理生命本身的原始数据。