bioinformatics-fundamentals

生物信息学基础知识

基因组学和生物信息学工作流的基础知识。帮助理解文件格式、测序技术和常见数据处理模式。

何时使用此技能

• 处理测序数据（PacBio HiFi、Hi-C、Illumina）
• 调试 SAM/BAM 比对或筛选问题
• 处理用于基因组组装整理的 AGP 文件
• 验证 AGP 坐标系统和未定位序列的分配
• 了解双端数据与单端数据的区别
• 解读质量指标（MAPQ、PHRED 分数）
• 排查输出为空或读段对损坏的问题
• 访问 GenomeArk QC 数据（GenomeScope、BUSCO、Merqury）
• 整理核型数据或分析染色体数量
• 开展常规生物信息学数据分析

SAM/BAM 格式要点

SAM 标志位（按位）

标志位可以叠加——一个读段可以同时设置多个标志位。

常见标志位：

• 0x0001 (1)：读段为双端测序数据
• 0x0002 (2)：每个片段均正确比对（proper pair）
• 0x0004 (4)：读段未比对
• 0x0008 (8)：配对读段未比对
• 0x0010 (16)：读段比对到反向链
• 0x0020 (32)：配对读段比对到反向链
• 0x0040 (64)：一对中的第一个读段（R1/正向）
• 0x0080 (128)：一对中的第二个读段（R2/反向）
• 0x0100 (256)：次级比对
• 0x0400 (1024)：PCR 或光学重复
• 0x0800 (2048)：补充比对

标志位组合：

• 正确配对的 R1：99 (0x63 = 1 + 2 + 32 + 64)
• 正确配对的 R2：147 (0x93 = 1 + 2 + 16 + 128)
• 未比对的读段：4
• 配对读段未比对：8

完整的标志位表、CIGAR 操作、可选标签和 SAM 必需字段请参阅 reference.md。

正确配对标志位（0x0002）

“正确配对”的含义：

• R1 和 R2 均已比对
• 比对方向正确（通常 R1 为正向，R2 为反向）
• 插入片段长度符合文库预期
• 配对方式符合比对工具的预期

**重要提示：**不同的比对工具对正确配对有不同的判定标准！

MAPQ（比对质量）

公式： MAPQ = -10 * log10(P(mapping is wrong))

常见阈值：

• MAPQ >= 60：置信度高（错误概率 < 0.0001%）
• MAPQ >= 30：质量良好（错误概率 < 0.1%）
• MAPQ >= 20：可接受（错误概率 < 1%）
• MAPQ >= 10：置信度低（错误概率 < 10%）
• MAPQ = 0：多重比对或未比对

**注意：**MAPQ=0 可能表示未比对，也可能表示存在多个同样好的比对位置。

CIGAR 字符串

表示读段与参考序列之间的比对：

• M：匹配或错配（比对匹配）
• I：读段相对于参考序列的插入
• D：读段相对于参考序列的缺失
• S：软剪切（读段中未比对的碱基）
• H：硬剪切（不在读段序列中的碱基）
• N：跳过的区域（用于 RNA-seq 剪接）

示例：100M = 完美匹配 100 bp 示例：50M5I45M = 匹配 50 bp、插入 5 bp、匹配 45 bp

测序技术

PacBio HiFi（高保真）

特点：

• 长读段：通常为 10–25 kb
• 高准确度：>99.9%（Q20+）
• 环状一致性测序（CCS）
• 单端数据（但来自环状分子）
• 非常适合从头组装

推荐比对工具：

• minimap2 预设：map-pb、map-hifi
• BWA-MEM2 也可使用，但针对短读段进行了优化

典型用途：

• 从头基因组组装
• 结构变异检测
• 异构体测序（Iso-Seq）
• 单倍型定相

Hi-C（染色质构象捕获）

特征：

• 双端短读段（通常为 100-150 bp）
• 读段对捕获染色质相互作用
• R1 和 R2 通常比对到不同的 scaffold/染色体
• 需要谨慎处理有效读段对
• 用于 scaffold 构建和三维基因组结构分析

推荐比对工具：

• BWA-MEM2（双端模式）
• BWA-MEM（双端模式）

关键概念： Hi-C 读段对会有意比对到相距较远的位点。区域筛选很容易破坏读段对！

典型用例：

• 基因组 scaffold 构建（连接 contig）
• 三维染色质结构分析
• 单倍型分相
• 组装质量评估

Illumina 短读段

特征：

• 短读段：50-300 bp
• 双端或单端
• 高通量
• 质量评分体系成熟

推荐比对工具：

• BWA-MEM2、BWA-MEM（通用）
• Bowtie2（速度快，局部比对）
• STAR（RNA-seq 剪接比对）

有关各技术的详细技术规格、错误特征和质量指标，请参阅 reference.md。

常用工具及其行为

samtools view

用途： 筛选、转换和查看 SAM/BAM 文件

主要标志：

• -b：输出 BAM 格式
• -h：包含头信息
• -f INT：要求包含指定标志（保留带有这些标志的读段）
• -F INT：筛除指定标志（移除带有这些标志的读段）
• -q INT：最低 MAPQ 阈值
• -L FILE：保留与 BED 文件中区域重叠的读段

重要行为：

• -L（区域筛选）会逐条检查读段，而非检查读段对
• 如果读段对的两个 mate 比对到不同区域，可能会破坏读段对
• 标志筛选（-f、-F）会在区域筛选（-L）之前执行

示例——保留区域内的有效读段对（正确顺序）：

samtools view -b -f 2 -L regions.bed input.bam > proper_pairs_in_regions.bam

bamtools filter

用途： 根据复杂条件进行高级筛选

常用筛选条件：

• isPaired: true - 读段来自双端测序
• isProperPair: true - 读段属于有效读段对
• isMapped: true - 读段已比对
• mapQuality: >=30 - 比对质量阈值

与 samtools 的重要区别：

• isProperPair 比 samtools -f 2 更严格
• 对读段对有效性的检查更全面

samtools fastx

用途： 将 SAM/BAM 转换为 FASTQ/FASTA

关键点： 使用合适的筛选条件，确保 R1/R2 文件相互匹配！

有关包含所有选项和示例的完整工具命令参考，请参阅 reference.md。

常见模式和最佳实践

模式 1：按区域筛选双端数据

错误方式（会破坏读段对）：

# Region filter first -> breaks pairs when mates are in different regions
samtools view -b -L regions.bed input.bam | bamtools filter -isPaired -isProperPair
# Result: Empty output (all pairs broken)

正确方式（保留读段对）：

# Proper pair filter FIRST, then region filter
samtools view -b -f 2 -L regions.bed input.bam > output.bam

模式 2：从筛选后的 BAM 提取 FASTQ

双端：

samtools fastx -1 R1.fq.gz -2 R2.fq.gz \
  --i1-flags 2 \  # Require proper pair
  input.bam

单端：

samtools fastx -0 output.fq.gz input.bam

模式 3：质量筛选

保守（高质量）：

samtools view -b -q 30 -f 2 -F 256 -F 2048 input.bam
# MAPQ >= 30, proper pairs, no secondary/supplementary

宽松（适用于低覆盖度数据）：

samtools view -b -q 10 -F 4 input.bam
# MAPQ >= 10, mapped reads

常见问题摘要

问题 1：区域筛选后输出为空（Hi-C 数据）

区域筛选（samtools view -L）会破坏读段对：一个 mate 在区域内，另一个在区域外，从而丢失有效读段对标志。请在区域筛选之前应用有效读段对筛选：

samtools view -b -f 2 -L regions.bed input.bam > output.bam

问题 2：R1 和 R2 文件的读段数不同

筛选不当破坏了部分读段对。提取时请要求保留有效读段对：

samtools fastx -1 R1.fq -2 R2.fq --i1-flags 2 input.bam

问题 3：Hi-C 数据的比对率偏低

由于 Hi-C 数据存在嵌合读段，比对率偏低是正常的。请使用 Hi-C 专用流程（HiC-Pro、Juicer）。不要过度按 MAPQ 筛选。

问题 4：比对后丢失有效读段对

请检查插入片段长度分布、参考序列不匹配或方向标志错误。

samtools stats input.bam | grep "insert size"
samtools flagstat input.bam

有关全面的故障排除指南和详细解决方案（包括 AGP 处理问题、HiFi 特有问题及诊断命令），请参阅 common-issues.md。

质量指标

N50 及相关指标

N50： 至少 50% 的组装结果位于长度不小于该值的 contig 中时，最短 contig 的长度

相关指标：

• L50： 达到 N50 所需的 contig 数量
• N90： 比 N50 更严格（覆盖率为 90%）
• NG50： 相对于基因组大小的 N50（更适合用于比较）

覆盖率和深度

覆盖率： 至少被一条 reads 覆盖的参考碱基所占的百分比 深度： 覆盖每个碱基的 reads 平均数量

推荐深度：

• 基因组组装（HiFi）：30-50x
• 变异检测：最低 30x
• RNA-seq：20-40 百万条 reads
• Hi-C scaffold 构建：50-100x 基因组覆盖度

完整的覆盖度计算、BUSCO 解读、QV 分数和组装质量指标，请参阅 reference.md。

文件格式速查

FASTA

>sequence_id description
ATCGATCGATCG

• 标题行以 > 开头
• 不包含质量分数

FASTQ

@read_id
ATCGATCGATCG
+
IIIIIIIIIIII

• 每条 read 占四行
• 包含质量分数（通常采用 Phred+33 编码）

BED

chr1    1000    2000    feature_name    score    +

• 坐标从 0 开始
• 半开区间 [start, end)

AGP

chr1    1    5000    1    W    contig_1    1    5000    +
chr1    5001 5100    2    U    100    scaffold    yes    proximity_ligation

• 以制表符分隔的基因组组装格式
• 闭区间坐标，从 1 开始 [start, end]
• 对象和组件的长度必须一致：obj_end - obj_beg + 1 == comp_end - comp_beg + 1

完整的 AGP 规范、坐标系统、验证规则和处理模式，请参阅 reference.md。 AGP 坐标调试和 unloc 处理问题，请参阅 common-issues.md。

坐标系统

从 1 开始（SAM、VCF、GFF、AGP）： 第一个碱基的位置为 1。区间 [2,5] 包含位置 2、3、4、5。 从 0 开始（BED、BAM 二进制格式）： 第一个碱基的位置为 0。区间 [2,5) 包含位置 2、3、4（不包括 5）。

转换： BED_start = SAM_start - 1；BED_end = SAM_end。

最佳实践

通用

1. 始终检查数据类型： 双端还是单端决定了过滤策略
2. 了解测序技术： Hi-C 的行为与 HiFi 不同
3. 按正确顺序过滤： 先筛选正确配对的 reads，再进行区域过滤
4. 验证输出： 检查文件大小、read 数量和 flagstat
5. 使用合适的 MAPQ 阈值： 过于严格会丢失数据，过于宽松会引入噪声

Hi-C 数据

1. 预期 reads 配对距离较远： 不要对它们位于不同 scaffold 感到意外
2. 保留正确配对的 reads： 这对下游 scaffold 构建至关重要
3. 使用支持配对的工具： 标准过滤器可能会破坏配对关系
4. 不要过度按 MAPQ 过滤： Hi-C 的 MAPQ 通常低于 DNA-seq

HiFi 数据

1. 单端处理： 无需考虑配对问题
2. 预期数据质量高： 可以使用严格的过滤条件
3. 使用合适的预设： minimap2 map-hifi 或 map-pb
4. 考虑 read 长度分布： HiFi reads 的长度各不相同

工具测试

1. 创建自包含的数据集： 所选区域中应包含两端 reads
2. 保持正确配对： 这是进行真实测试的必要条件
3. 使用有代表性的数据： 按比例抽样，而不是随机抽样
4. 验证文件大小： 文件过小表示过滤过度

相关技能

• vgp-pipeline - VGP 工作流处理 Hi-C 和 HiFi 数据
• galaxy-tool-wrapping - Galaxy 工具可处理 SAM/BAM 和测序数据格式
• galaxy-workflow-development - 工作流处理测序数据

支持文档

• **reference.md：**详细格式规范（SAM/BAM 完整参考、CIGAR 操作、AGP 格式、FASTQ 编码、工具命令参考、测序技术规范、组装质量指标、覆盖度计算、坐标系统）
• **common-issues.md：**全面的故障排查指南（输出为空、双端测序问题、质量/比对问题、格式转换、Hi-C 和 HiFi 特定问题、AGP 处理错误、诊断命令）
• **genomeark-data-access.md：**GenomeArk AWS S3 数据访问模式（目录结构演变、GenomeScope/BUSCO/Merqury 的 QC 数据位置、获取策略、路径规范化、组装日期提取）
• **genomic-analysis-patterns.md：**领域特定分析模式（核型数据整理、单倍体与二倍体染色体计数、系统发育树物种映射、BED/端粒分析、NCBI 数据集成策略）

版本历史

• **v1.2.0：**拆分为 SKILL.md 和支持文件，以便维护；将 GenomeArk 访问、核型/染色体分析、系统发育映射、AGP 详情以及端粒/NCBI 模式移至支持文件
• **v1.1.1：**新增用于端粒分析的 BED 文件处理模式和 NCBI 数据集成策略
• **v1.1.0：**新增全面的 AGP 格式文档，包括坐标验证、unloc 处理和常见错误模式
• **v1.0.0：**初始版本，包含 SAM/BAM、Hi-C、HiFi 和常见过滤模式