全外显子组测序（WES）分析全流程¶

一句话概述¶

全外显子组测序（WES）分析从FASTQ文件出发，经过质控、比对（BWA-MEM2）、去重、碱基质量校正、变异检测（GATK HaplotypeCaller）、过滤和注释，最终获得高置信度的编码区变异（SNV/InDel）列表。

核心知识点表格¶

知识点	说明
WES vs WGS	WES只测编码区（~1-2%基因组），成本低但覆盖面有限
BWA-MEM2	BWA-MEM的SIMD加速版本，速度提升2-3倍
GATK Best Practices	Broad Institute推荐的变异检测标准流程
HaplotypeCaller	GATK核心变异检测器，基于局部重组装和PairHMM
BQSR	Base Quality Score Recalibration，校正系统性碱基质量偏差
VQSR/Hard filtering	变异质量过滤：机器学习（VQSR）或固定阈值（Hard filter）
gVCF	Genomic VCF，包含所有位点信息（含非变异位点）
Mark Duplicates	标记PCR重复reads，防止假阳性
Capture Kit/BED	WES使用的外显子捕获区域定义文件
VEP/ANNOVAR	变异功能注释工具

各步骤详解¶

第一步：原始数据质控¶

白话解释： 先看看测序数据质量如何——有没有接头污染、碱基质量是否达标、是否有过多的N碱基。

技术细节： - WES数据通常为paired-end 150bp - 平均测序深度目标：100-200x（临床）或50-100x（科研） - 外显子区域覆盖度是核心指标 - 接头序列需要trimming

代码示例：

# 1. FastQC质控
fastqc -t 8 -o fastqc_results/ sample_R1.fastq.gz sample_R2.fastq.gz

# 2. MultiQC汇总
multiqc fastqc_results/ -o multiqc_report/

# 3. 接头去除和质量修剪（fastp）
fastp --in1 sample_R1.fastq.gz \
      --in2 sample_R2.fastq.gz \
      --out1 sample_clean_R1.fastq.gz \
      --out2 sample_clean_R2.fastq.gz \
      --json sample_fastp.json \
      --html sample_fastp.html \
      --thread 8 \
      --qualified_quality_phred 20 \
      --length_required 50 \
      --detect_adapter_for_pe \
      --correction

第二步：序列比对¶

白话解释： 把测序reads比对到人类参考基因组上，确定每条read来自哪个基因组位置。BWA-MEM2是目前最快的短reads比对工具。

技术细节： - 参考基因组：GRCh38/hg38（推荐）或GRCh37/hg19 - BWA-MEM2比BWA-MEM快约2-3x（利用AVX2/AVX512指令集） - 输出BAM文件包含比对坐标、MAPQ、CIGAR等信息 - Read Group信息必须正确设置（后续GATK需要）

代码示例：

# 1. 建立参考基因组索引（只需一次）
bwa-mem2 index /ref/GRCh38.fa

# 2. 比对（设置正确的Read Group）
bwa-mem2 mem -t 16 \
  -R "@RG\tID:sample1\tSM:sample1\tPL:ILLUMINA\tLB:lib1\tPU:unit1" \
  /ref/GRCh38.fa \
  sample_clean_R1.fastq.gz \
  sample_clean_R2.fastq.gz | \
  samtools sort -@ 8 -o sample1_sorted.bam -

# 3. 建立BAM索引
samtools index sample1_sorted.bam

# 4. 比对统计
samtools flagstat sample1_sorted.bam > sample1_flagstat.txt
samtools stats sample1_sorted.bam > sample1_stats.txt

# 5. 查看on-target率（需要capture BED文件）
# Picard CollectHsMetrics
picard CollectHsMetrics \
  I=sample1_sorted.bam \
  O=sample1_hs_metrics.txt \
  BAIT_INTERVALS=capture_targets.interval_list \
  TARGET_INTERVALS=capture_targets.interval_list \
  R=/ref/GRCh38.fa

第三步：去除PCR重复¶

白话解释： PCR扩增过程中，同一条DNA片段可能被复制多次，形成"PCR duplicates"。这些重复reads如果不去除，会让同一个位置的变异信号被人为放大，产生假阳性。

技术细节： - 识别标准：起始位置+终止位置+UMI（如有）完全相同的reads - GATK MarkDuplicates标记而非删除（设置FLAG=1024） - Optical duplicates（同一cluster的重复）也被处理 - WES数据典型duplication rate：10-30%

代码示例：

# 使用GATK MarkDuplicates
gatk MarkDuplicates \
  -I sample1_sorted.bam \
  -O sample1_markdup.bam \
  -M sample1_markdup_metrics.txt \
  --REMOVE_DUPLICATES false \
  --CREATE_INDEX true \
  --VALIDATION_STRINGENCY SILENT

# 或使用sambamba（更快）
sambamba markdup -t 8 \
  sample1_sorted.bam \
  sample1_markdup.bam

# 查看duplication rate
grep "PERCENT_DUPLICATION" sample1_markdup_metrics.txt

第四步：碱基质量校正（BQSR）¶

白话解释： 测序仪报告的碱基质量分数可能存在系统偏差（比如某个测序周期、某个碱基对总是被高估或低估质量）。BQSR通过已知的SNP位点数据库来"校准"这些质量分数，使后续变异检测更准确。

技术细节： - 需要已知变异位点数据库（dbSNP、Mills、1000G等） - 第一步BaseRecalibrator：建立校正模型 - 第二步ApplyBQSR：应用校正 - BQSR对WES影响通常较小但仍推荐执行 - 校正基于：read group、位置、上下文序列、测序周期

代码示例：

# 下载必要的已知位点数据库（从GATK resource bundle获取）
# https://console.cloud.google.com/storage/browser/genomics-public-data/resources/broad/hg38/v0
# dbSNP: Homo_sapiens_assembly38.dbsnp138.vcf.gz
# Mills: Mills_and_1000G_gold_standard.indels.hg38.vcf.gz
# 1000G: 1000G_phase1.snps.high_confidence.hg38.vcf.gz

# 1. 建立校正表
gatk BaseRecalibrator \
  -I sample1_markdup.bam \
  -R /ref/GRCh38.fa \
  --known-sites /ref/Homo_sapiens_assembly38.dbsnp138.vcf.gz \
  --known-sites /ref/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz \
  --known-sites /ref/1000G_phase1.snps.high_confidence.hg38.vcf.gz \
  -O sample1_recal_data.table \
  -L /ref/exome_capture.bed  # 限制在外显子区域加速

# 2. 应用校正
gatk ApplyBQSR \
  -I sample1_markdup.bam \
  -R /ref/GRCh38.fa \
  --bqsr-recal-file sample1_recal_data.table \
  -O sample1_recal.bam \
  -L /ref/exome_capture.bed

# 3. （可选）评估校正效果
gatk BaseRecalibrator \
  -I sample1_recal.bam \
  -R /ref/GRCh38.fa \
  --known-sites /ref/Homo_sapiens_assembly38.dbsnp138.vcf.gz \
  --known-sites /ref/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz \
  -O sample1_recal_data_post.table

gatk AnalyzeCovariates \
  -before sample1_recal_data.table \
  -after sample1_recal_data_post.table \
  -plots recalibration_plots.pdf

第五步：变异检测（HaplotypeCaller）¶

白话解释： HaplotypeCaller是GATK的核心工具。它在每个活跃区域（有变异信号的地方）进行局部重新组装，生成所有可能的单倍型（haplotype），然后用PairHMM算法计算每条read支持各单倍型的概率，最终确定最可能的基因型。

技术细节： - 输出gVCF（含reference blocks）便于多样本联合分析 - ERC（Emit Reference Confidence）模式记录非变异位点信息 - 对于群体研究：HC→gVCF→GenomicsDBImport→GenotypeGVCFs - 对于单样本/肿瘤：直接出VCF或使用Mutect2

代码示例：

# === 单样本模式（直接出VCF） ===
gatk HaplotypeCaller \
  -I sample1_recal.bam \
  -R /ref/GRCh38.fa \
  -O sample1_raw.vcf.gz \
  -L /ref/exome_capture_padded.bed \
  --native-pair-hmm-threads 4 \
  -D /ref/Homo_sapiens_assembly38.dbsnp138.vcf.gz

# === gVCF模式（多样本联合分析推荐） ===
gatk HaplotypeCaller \
  -I sample1_recal.bam \
  -R /ref/GRCh38.fa \
  -O sample1.g.vcf.gz \
  -L /ref/exome_capture_padded.bed \
  -ERC GVCF \
  --native-pair-hmm-threads 4 \
  -D /ref/Homo_sapiens_assembly38.dbsnp138.vcf.gz

# === 多样本联合基因型分析 ===
# 1. 合并gVCF（GenomicsDBImport）
gatk GenomicsDBImport \
  -V sample1.g.vcf.gz \
  -V sample2.g.vcf.gz \
  -V sample3.g.vcf.gz \
  --genomicsdb-workspace-path genomicsdb \
  -L /ref/exome_capture_padded.bed \
  --tmp-dir /tmp/

# 2. 联合基因型分析
gatk GenotypeGVCFs \
  -R /ref/GRCh38.fa \
  -V gendb://genomicsdb \
  -O cohort_raw.vcf.gz \
  -D /ref/Homo_sapiens_assembly38.dbsnp138.vcf.gz

第六步：变异过滤¶

白话解释： 原始VCF里有很多假阳性变异需要过滤掉。有两种策略：(1) VQSR用机器学习（高斯混合模型）学习真阳性变异的特征分布，适合大样本；(2) Hard filtering用固定阈值，适合小样本或非人物种。

技术细节： - VQSR需要>30个外显子组样本（SNP）或更多样本 - VQSR对外显子组数据效果不如WGS（因为位点数少） - WES推荐对SNP使用VQSR、对InDel使用Hard filtering - Hard filtering GATK推荐阈值见下

代码示例：

# === 方法1：VQSR（大队列推荐） ===
# SNP VQSR
gatk VariantRecalibrator \
  -V cohort_raw.vcf.gz \
  -R /ref/GRCh38.fa \
  --resource:hapmap,known=false,training=true,truth=true,prior=15.0 hapmap_3.3.hg38.vcf.gz \
  --resource:omni,known=false,training=true,truth=false,prior=12.0 1000G_omni2.5.hg38.vcf.gz \
  --resource:1000G,known=false,training=true,truth=false,prior=10.0 1000G_phase1.snps.high_confidence.hg38.vcf.gz \
  --resource:dbsnp,known=true,training=false,truth=false,prior=2.0 Homo_sapiens_assembly38.dbsnp138.vcf.gz \
  -an QD -an MQ -an MQRankSum -an ReadPosRankSum -an FS -an SOR \
  -mode SNP \
  -O snp_recal.recal \
  --tranches-file snp_recal.tranches

gatk ApplyVQSR \
  -V cohort_raw.vcf.gz \
  -R /ref/GRCh38.fa \
  --recal-file snp_recal.recal \
  --tranches-file snp_recal.tranches \
  --truth-sensitivity-filter-level 99.5 \
  -mode SNP \
  -O cohort_snp_recal.vcf.gz

# === 方法2：Hard Filtering（小样本/单样本） ===
# 分离SNP和InDel
gatk SelectVariants -V sample1_raw.vcf.gz -select-type SNP -O snps_raw.vcf.gz
gatk SelectVariants -V sample1_raw.vcf.gz -select-type INDEL -O indels_raw.vcf.gz

# SNP过滤
gatk VariantFiltration \
  -V snps_raw.vcf.gz \
  -R /ref/GRCh38.fa \
  --filter-expression "QD < 2.0" --filter-name "QD2" \
  --filter-expression "FS > 60.0" --filter-name "FS60" \
  --filter-expression "MQ < 40.0" --filter-name "MQ40" \
  --filter-expression "MQRankSum < -12.5" --filter-name "MQRankSum-12.5" \
  --filter-expression "ReadPosRankSum < -8.0" --filter-name "ReadPosRankSum-8" \
  --filter-expression "SOR > 3.0" --filter-name "SOR3" \
  -O snps_filtered.vcf.gz

# InDel过滤
gatk VariantFiltration \
  -V indels_raw.vcf.gz \
  -R /ref/GRCh38.fa \
  --filter-expression "QD < 2.0" --filter-name "QD2" \
  --filter-expression "FS > 200.0" --filter-name "FS200" \
  --filter-expression "ReadPosRankSum < -20.0" --filter-name "ReadPosRankSum-20" \
  --filter-expression "SOR > 10.0" --filter-name "SOR10" \
  -O indels_filtered.vcf.gz

# 合并
gatk MergeVcfs \
  -I snps_filtered.vcf.gz \
  -I indels_filtered.vcf.gz \
  -O sample1_filtered.vcf.gz

# 只保留PASS的变异
bcftools view -f PASS sample1_filtered.vcf.gz -o sample1_PASS.vcf.gz -Oz

第七步：变异注释¶

白话解释： 拿到过滤后的VCF后，需要知道每个变异意味着什么——它在哪个基因？会不会导致蛋白质改变？在人群中多常见？有没有致病性报告？

技术细节： - VEP (Variant Effect Predictor)：Ensembl官方工具，最全面 - ANNOVAR：速度快、数据库丰富 - 注释维度：基因影响（missense/nonsense等）、人群频率（gnomAD）、致病性预测（SIFT/PolyPhen/CADD）、临床数据库（ClinVar）

代码示例：

# === VEP注释 ===
vep --input_file sample1_PASS.vcf.gz \
    --output_file sample1_annotated.vcf \
    --vcf \
    --offline \
    --cache \
    --dir_cache /ref/vep_cache/ \
    --species homo_sapiens \
    --assembly GRCh38 \
    --everything \
    --fork 8 \
    --plugin CADD,/ref/CADD/whole_genome_SNVs.tsv.gz \
    --plugin dbNSFP,/ref/dbNSFP/dbNSFP4.7a.gz,SIFT_pred,Polyphen2_HDIV_pred,MutationTaster_pred,AlphaMissense_pred \
    --custom /ref/clinvar.vcf.gz,ClinVar,vcf,exact,0,CLNSIG,CLNREVSTAT

# === ANNOVAR注释 ===
# 下载数据库
annotate_variation.pl -buildver hg38 -downdb -webfrom annovar refGene humandb/
annotate_variation.pl -buildver hg38 -downdb -webfrom annovar gnomad30_genome humandb/
annotate_variation.pl -buildver hg38 -downdb -webfrom annovar clinvar_20250721 humandb/
annotate_variation.pl -buildver hg38 -downdb -webfrom annovar dbnsfp47a humandb/

# 运行注释
table_annovar.pl sample1_PASS.vcf.gz humandb/ \
  -buildver hg38 \
  -out sample1_annovar \
  -remove \
  -protocol refGene,gnomad30_genome,clinvar_20250721,dbnsfp47a \
  -operation g,f,f,f \
  -nastring . \
  -vcfinput \
  -polish

# === 致病变异筛选（示例） ===
# 筛选条件：
# 1. 功能影响：nonsynonymous/stopgain/splicing
# 2. 人群频率：gnomAD AF < 0.01 (1%)
# 3. 致病性预测：CADD > 20 或 ClinVar Pathogenic
bcftools query -f '%CHROM\t%POS\t%REF\t%ALT\t%INFO/CSQ\n' sample1_annotated.vcf | \
  grep -E "missense|stop_gained|splice_donor|splice_acceptor" > candidate_variants.txt

实战命令（可复制）¶

# ===== 完整WES分析流水线 =====

SAMPLE="patient01"
REF="/ref/GRCh38.fa"
DBSNP="/ref/Homo_sapiens_assembly38.dbsnp138.vcf.gz"
MILLS="/ref/Mills_and_1000G_gold_standard.indels.hg38.vcf.gz"
KG="/ref/1000G_phase1.snps.high_confidence.hg38.vcf.gz"
BED="/ref/exome_capture_padded.bed"
THREADS=16

# 1. 质控 + trim
fastp --in1 ${SAMPLE}_R1.fastq.gz --in2 ${SAMPLE}_R2.fastq.gz \
  --out1 ${SAMPLE}_clean_R1.fastq.gz --out2 ${SAMPLE}_clean_R2.fastq.gz \
  --json ${SAMPLE}.fastp.json --html ${SAMPLE}.fastp.html \
  --thread 8 --detect_adapter_for_pe

# 2. 比对
bwa-mem2 mem -t ${THREADS} \
  -R "@RG\tID:${SAMPLE}\tSM:${SAMPLE}\tPL:ILLUMINA\tLB:lib1\tPU:unit1" \
  ${REF} ${SAMPLE}_clean_R1.fastq.gz ${SAMPLE}_clean_R2.fastq.gz | \
  samtools sort -@ 8 -o ${SAMPLE}_sorted.bam -

samtools index ${SAMPLE}_sorted.bam

# 3. 去重
gatk MarkDuplicates \
  -I ${SAMPLE}_sorted.bam -O ${SAMPLE}_markdup.bam \
  -M ${SAMPLE}_dup_metrics.txt --CREATE_INDEX true

# 4. BQSR
gatk BaseRecalibrator \
  -I ${SAMPLE}_markdup.bam -R ${REF} \
  --known-sites ${DBSNP} --known-sites ${MILLS} --known-sites ${KG} \
  -O ${SAMPLE}_recal.table -L ${BED}

gatk ApplyBQSR \
  -I ${SAMPLE}_markdup.bam -R ${REF} \
  --bqsr-recal-file ${SAMPLE}_recal.table \
  -O ${SAMPLE}_recal.bam -L ${BED}

# 5. 变异检测
gatk HaplotypeCaller \
  -I ${SAMPLE}_recal.bam -R ${REF} \
  -O ${SAMPLE}.g.vcf.gz -ERC GVCF \
  -L ${BED} -D ${DBSNP}

# 6. 基因型分析（单样本简化版）
gatk GenotypeGVCFs \
  -R ${REF} -V ${SAMPLE}.g.vcf.gz -O ${SAMPLE}_raw.vcf.gz

# 7. 过滤（Hard filtering）
gatk SelectVariants -V ${SAMPLE}_raw.vcf.gz -select-type SNP -O ${SAMPLE}_snps.vcf.gz
gatk VariantFiltration -V ${SAMPLE}_snps.vcf.gz -R ${REF} \
  --filter-expression "QD<2.0||FS>60.0||MQ<40.0||SOR>3.0" \
  --filter-name "SNP_FILTER" -O ${SAMPLE}_snps_filtered.vcf.gz

# 8. 注释
table_annovar.pl ${SAMPLE}_snps_filtered.vcf.gz humandb/ \
  -buildver hg38 -out ${SAMPLE}_annotated \
  -protocol refGene,gnomad30_genome,clinvar_20250721 \
  -operation g,f,f -nastring . -vcfinput -polish

面试常问点¶

Q1: WES和WGS各自的优缺点？¶

A: WES优点：成本低（WGS的1/3-1/5），外显子区深度高（同样cost下覆盖度更高），数据量小便于存储分析。WES缺点：无法检测非编码区变异（调控区、内含子深处的splicing变异）、无法可靠检测结构变异（SV）、捕获偏好性导致某些区域覆盖不均。

Q2: BWA-MEM2相比BWA-MEM有什么改进？¶

A: BWA-MEM2利用CPU的SIMD指令（AVX2/AVX512/SSE）对Smith-Waterman对齐和种子扩展进行向量化加速，速度提升2-3倍。算法逻辑和比对结果与BWA-MEM基本一致（几乎identical的输出），只是底层实现优化。

Q3: 为什么GATK需要Mark Duplicates而不是直接Remove？¶

A: 标记（mark）而非删除的原因：(1) 保留所有数据便于后续重分析；(2) 某些分析可能需要duplicates信息；(3) HaplotypeCaller会自动忽略标记为duplicate的reads；(4) 如果有UMI信息可以后续更精确地去重。

Q4: BQSR的原理是什么？¶

A: BQSR观察到：在已知不应有变异的位点（即非dbSNP位点），如果某read的某碱基与参考不同，那这个碱基可能是测序错误。它按照read group、位置context、碱基对类型、测序周期等维度分组统计错误率，与质量分数比较并建立校正模型。校正后的质量分数更准确反映真实错误概率。

Q5: Hard filtering各参数的含义？¶

A: QD（Quality by Depth）：变异质量除以覆盖深度，<2.0表示质量不够；FS（FisherStrand）：正负链偏好的Fisher's检验p值，>60表示严重链偏好；MQ（Mapping Quality）：比对质量均值，<40表示映射不可靠；SOR（Symmetric Odds Ratio）：链偏好的另一个指标；ReadPosRankSum：变异在reads中位置的偏差。

Q6: 什么时候用VQSR，什么时候用Hard filtering？¶

A: VQSR：(1) 人类WGS/WES大队列（>30样本）；(2) 有足够的训练数据（known variants）。Hard filtering：(1) 小样本（<30）；(2) 非模式生物；(3) WES数据中InDel（位点数可能不够训练VQSR）；(4) 靶向测序（panel）。

Q7: gVCF模式的意义？¶

A: gVCF不仅记录变异位点，还以reference block形式记录非变异位点的信息（包括该位点的置信度）。好处：(1) 多样本联合分析时能区分"没有变异"和"没有数据"；(2) 可以灵活地随时加入新样本重新分析而无需重跑HaplotypeCaller；(3) 便于大队列分析（GenomicsDBImport）。

Q8: 如何评估WES分析结果的质量？¶

A: 关键指标：(1) 比对率>95%；(2) On-target率>70%（target区域reads占比）；(3) 目标区域平均深度>50x；(4) 目标区域>80%的碱基覆盖度≥20x；(5) Duplication rate<30%；(6) Ti/Tv比（全外显子组应为~3.0-3.3）；(7) Het/Hom比（非近亲群体应为~1.5-2.0）。

Q9: ClinVar和gnomAD在变异解读中的作用？¶

A: gnomAD提供人群频率信息——如果一个变异在健康人群中频率>1%，通常不是致病变异（常见多态性）。ClinVar提供临床注释——收集了已报道的变异-疾病关联和致病性分类（Pathogenic/Likely pathogenic/VUS/Benign）。两者结合用于ACMG变异分类。

易错点¶

1. Read Group信息设置错误或缺失¶

问题： GATK HaplotypeCaller要求BAM文件必须有正确的@RG header，缺失会报错。 正确做法： 在BWA比对时通过-R参数指定完整RG信息：ID（唯一标识）、SM（样本名）、PL（测序平台）、LB（文库）、PU（测序单元）。

2. 未使用capture BED文件限制分析区域¶

问题： WES数据在非捕获区域覆盖很低，对全基因组做变异检测会产生大量假阳性。 正确做法： 在HaplotypeCaller中使用-L指定capture BED（最好左右各padding 100-150bp），只在目标区域检测变异。

3. BQSR使用了不匹配的参考基因组版本¶

问题： 参考基因组是hg38但known sites文件是hg19，坐标完全不匹配。 正确做法： 确保所有文件（参考基因组、known sites、BED）使用同一基因组版本。必要时用LiftOver转换坐标。

4. Hard filtering阈值不区分SNP和InDel¶

问题： 对SNP和InDel使用相同过滤阈值，但它们的质量指标分布不同。 正确做法： 先用SelectVariants分离SNP和InDel，分别应用不同的阈值（GATK推荐值不同）。

5. 注释时忽略了变异的多等位基因分裂¶

问题： VCF中一行可能包含多个ALT allele（multi-allelic），注释工具可能只处理第一个。 正确做法： 用bcftools norm -m -将multi-allelic拆分为bi-allelic后再注释。同时用bcftools norm -f ref.fa进行左对齐。

6. 混淆germline和somatic变异检测¶

问题： 用HaplotypeCaller检测肿瘤somatic突变，或用Mutect2检测遗传变异。 正确做法： Germline变异用HaplotypeCaller；肿瘤somatic突变用Mutect2（需要tumor-normal配对）。

补充知识¶

GATK资源下载¶

# Google Cloud GATK resource bundle
# https://console.cloud.google.com/storage/browser/genomics-public-data/resources/broad/hg38/v0

# 主要文件：
# - Homo_sapiens_assembly38.fasta（参考基因组）
# - Homo_sapiens_assembly38.dbsnp138.vcf.gz（注意：bundle中的dbSNP版本为Build 138，
#   最新的dbSNP已到Build 157(2025.3)，但GATK VQSR流程中dbSNP仅用作"known"标记，
#   Build 138已足够；如需最新频率信息请单独下载dbSNP 157）
# - Mills_and_1000G_gold_standard.indels.hg38.vcf.gz
# - 1000G_phase1.snps.high_confidence.hg38.vcf.gz
# - hapmap_3.3.hg38.vcf.gz
# - 1000G_omni2.5.hg38.vcf.gz

WES质控指标参考¶

指标	合格标准	优秀标准
比对率	>95%	>98%
On-target率	>60%	>75%
目标区域平均深度	>50x	>100x
≥20x覆盖率	>80%	>95%
Duplication rate	<30%	<15%
Insert size均值	150-300bp	200±30bp

加速建议¶

BWA-MEM2：比BWA-MEM快2-3x
GATK Spark模式：分布式计算（不太稳定）
DeepVariant：Google的深度学习变异检测器，准确度高，支持GPU加速（也可CPU运行）
Sentieon：GATK的商业加速版（需license）
按染色体并行：将BED按染色体拆分，各染色体独立运行HC后合并

本教程涵盖WES分析从原始数据到注释变异的完整GATK Best Practices流程，适合临床和科研WES数据分析。

全外显子组测序（WES）分析全流程¶

一句话概述¶

核心知识点表格¶

各步骤详解¶

第一步：原始数据质控¶

第二步：序列比对¶

第三步：去除PCR重复¶

第四步：碱基质量校正（BQSR）¶

第五步：变异检测（HaplotypeCaller）¶

第六步：变异过滤¶

第七步：变异注释¶

实战命令（可复制）¶

面试常问点¶

Q1: WES和WGS各自的优缺点？¶

Q2: BWA-MEM2相比BWA-MEM有什么改进？¶

Q3: 为什么GATK需要Mark Duplicates而不是直接Remove？¶

Q4: BQSR的原理是什么？¶

Q5: Hard filtering各参数的含义？¶

Q6: 什么时候用VQSR，什么时候用Hard filtering？¶

Q7: gVCF模式的意义？¶

Q8: 如何评估WES分析结果的质量？¶

Q9: ClinVar和gnomAD在变异解读中的作用？¶

易错点¶

1. Read Group信息设置错误或缺失¶

2. 未使用capture BED文件限制分析区域¶

3. BQSR使用了不匹配的参考基因组版本¶

4. Hard filtering阈值不区分SNP和InDel¶

5. 注释时忽略了变异的多等位基因分裂¶

6. 混淆germline和somatic变异检测¶

补充知识¶

GATK资源下载¶

WES质控指标参考¶

加速建议¶

📚 相关文章推荐