前文回顾
1. GATK官方教程 / 概述及工作前的布置
2. GATK教程 / 体细胞短变异检测 (SNV+InDel)流程概览
3. GATK教程 / 变异检测前的数据预处理
4. GATK / 体细胞短变异检测工具Mutect2的使用
5. Mutect2案例 / 有或没有正常样本配对的肿瘤体细胞变异检测结果对比
体细胞突变结果一般记录在VCF文件中，如果不可视化和绘图，人是没办法阅读的。
突变注释格式 (Mutation Annotation Format, MAF)，通常用于存储检测到的体细胞突变，并可用于后续maftools R包。从而在基因水平上观察不同肿瘤样本的突变图景，结果类似这样：
https://qiita.com/kojix2/items/738ade193e0718509c46https://bioconductor.org/packages/release/bioc/html/maftools.html 这就需要将上游分析获得的VCF文件，注释、整理成标准的Maf格式。需要注意： ① VCF注释 与 ② 注释后整理成maf文件格式，二者可能涉及两类不同的工具(那么能否完美衔接，就是一个潜在的问题)。另外，VCF注释有时需要区分其记录的变异是胚系变异，还是体细胞变异。例如GATK Funcotator注释时，对这两类VCF文件是严格区分的 (调用了完全不同的数据库，当然基因组版本是一样的，只是涉及肿瘤驱动突变的功能注释时，调用的数据库可能有所区别)。 ① VCF注释。由于VCF文件注释工具很多，如：SnpEff、ANNOVAR、GATK Funcotator、VEP、CADD等，不同软件依赖不同的运行环境和注释数据库。 注释后产生的、新的VCF文件，虽然格式上大部分一致，但在FORMAT列中 (如：AD、DP、GQ、AF)，这些子类的个数及排序各有各的不同。这种不同其实不是VCF注释带来的 (不同软件的注释结果一般都会放在INFO列，且注释内容的格式编排不同)，而是由于VCF在产生时使用的不同程序所带来的，如：GATK HaplotypeCaller、GATK Mutect2、vcftools、bcftools)。这导致已注释的VCF文件，发现其FORMAT列可能非常不同，如：GATK HaplotypeCaller (Germline) + bcftools/vcftools GT:AD:DP:GQ:PL 1/1:0,64:64:99:2117,192,0 GT:AD:AF:DP:F1R2:F2R1:FAD:SB 0/1:0,19:0.953:19:0,8:0,6:0,18:0,0,18,1 GT:PL:DP:AD 0/0:0,255,255:598:596,1 ② 注释后整理成maf文件格式。在整理成maftools所要求的Maf格式时，需要使用不同的工具，如snpeffToMaf.pl、GATK Funcotator (没错，此工具既能注释VCF，也能将注释结果转为Maf文件；ANNOVAR应该也可以)。VCF to Maf的技术路径
ANNOVAR的VCF注释+Maf文件整理，网上有很多资料，我个人以前常用，后来舍弃 (非黑，原因很复杂，主要是个人习惯问题)，所以这里不再介绍。
SnpEff对VCF文件注释后，需要配合使用snpeffToMaf.pl转为MAF格式，再使用maftools绘图。但上面介绍过，不同软件获得的VCF文件的FORMAT列的内容、排序很不同，比如AF值一般只有Mutect2程序提供。因此，snpeffToMaf.pl脚本需要不断被修改，以适应不同形式的FORMAT，生成标准MAF文件。
用户自己改就比较麻烦，且容易出错，因此不建议自行修改snpeffToMaf.pl。最好是使用VCF注释工具自动输出的Maf文件。
下面是GATK Funcotator (注释数据库下载方式见文末)注释VCF、并自动输出Maf格式的参考代码 (使用了rush并行计算)：
# 注释数据库：# funcotator_dataSources.v1.7.20200521scat ${result}/metadata.txt | cut -f 1 | sort -u | \ rush -j 6 "gatk Funcotator -R ${ref} \ -V ${result}/{}.mt2.clean.vcf.gz \ --output-file-format MAF \ -O ${result}/{}.mt2.clean.funcotator.maf \ --data-sources-path ~/gatk_ex/funcotator/funcotator_dataSources.v1.7.20200521s/ \ --ref-version hg38"Maf格式长什么样子？
① snpeffToMaf.pl输出的Maf文件 (表头/列的名称)：
Hugo_Symbol Entrez_Gene_Id Center CDS_Change
Chromosome Start_Position End_Position StrandVariant_Classification Variant_Type Reference_Allele Tumor_Seq_Allele1 Tumor_Seq_Allele2 Tumor_Sample_BarcodeProtein_Change i_TumorVAF_WU i_transcript_name下面这3列，不在GATK Funcotator注释结果中；审查其内容后发现是有的，只是列的名称不同：
CDS_Change | i_TumorVAF_WU | i_transcript_name
i_TumorVAF_WU i_transcript_name1.00000 ENST00000335137.30.98824 ENST00000335137.30.35156 ENST00000635509.10.64286 ENST00000635509.10.34286 ENST00000429505.5VAF：变异等位基因频率 (Variant allele frequency, VAF)
② gatk Funcotator输出的Maf文件 (表头/列的名称)：
Genome_Change/Annotation_Transcript/cDNA_Change/Codon_Change/Protein_Changeg.chr1:3728151C>T ENST00...4.5 c.1008C>T c.(1006-1008)gcC>gcT p.A336Ag.chr1:3728190T>C ENST00...4.5 c.1047T>C c.(1045-1047)caT>caC p.H349H注意：下划线数字的倍数关系 (约3倍，即3联体密码子)
gatk Funcotator注释的Maf文件一共217列：
17 Matched_Norm_Sample_Barcode 18 Match_Norm_Seq_Allele1 19 Match_Norm_Seq_Allele2 20 Tumor_Validation_Allele1 21 Tumor_Validation_Allele2 22 Match_Norm_Validation_Allele1 23 Match_Norm_Validation_Allele2 34 Matched_Norm_Sample_UUID 52 UniProt_Natural_Variations 53 UniProt_Experimental_Info 57 COSMIC_overlapping_mutations 59 COSMIC_tissue_types_affected 60 COSMIC_total_alterations_in_gene 61 Tumorscape_Amplification_Peaks 62 Tumorscape_Deletion_Peaks 63 TCGAscape_Amplification_Peaks 64 TCGAscape_Deletion_Peaks 68 CCLE_ONCOMAP_overlapping_mutations 69 CCLE_ONCOMAP_total_mutations_in_gene 71 CGC_Translocation_Partner 72 CGC_Tumor_Types_Somatic 73 CGC_Tumor_Types_Germline 75 DNARepairGenes_Activity_linked_to_OMIM 76 FamilialCancerDatabase_Syndromes 77 MUTSIG_Published_Results 85 Gencode_34_secondaryVariantClassification 92 ClinVar_VCF_CLNDISDBINCL 96 ClinVar_VCF_CLNREVSTAT 98 ClinVar_VCF_CLNSIGCONF 99 ClinVar_VCF_CLNSIGINCL 111 CosmicFusion_fusion_id 112 Familial_Cancer_Genes_Synonym 113 Familial_Cancer_Genes_Reference 114 Gencode_XHGNC_hgnc_id 119 HGNC_Previous_Symbols 125 HGNC_Date_Symbol_Changed 126 HGNC_Date_Name_Changed 127 HGNC_Accession_Numbers 133 HGNC_Gene_Family_Name 136 HGNC_OMIM_ID(supplied_by_OMIM) 137 HGNC_RefSeq(supplied_by_NCBI) 138 HGNC_UniProt_ID(supplied_by_UniProt) 139 HGNC_Ensembl_ID(supplied_by_Ensembl) 140 HGNC_UCSC_ID(supplied_by_UCSC) 142 Simple_Uniprot_alt_uniprot_accessions 190 HGNC_Entrez_Gene_ID(supplied_by_NCBI) 195 AS_UNIQ_ALT_READ_COUNT 217 TLOD
结论：gatk Funcotator 算出来的 tumor_f 列不完全等于 snpeffToMaf.pl 算出来的 i_TumorVAF_WU 列；虽然都可能是由AD (Allelic depths for the ref and alt alleles in the order listed)推导而来 (但在肿瘤体细胞突变中不完全等于：alt/(alt+ref)，如下图)。6/36=0.1666667 | 6/42=0.1428571
4/39=0.1025641 | 4/43 [1] 0.09302326
9/70=0.1285714 | 9/79=0.1139241
158/(158+254)=0.3834951 | 9/14=0.6428571(附)GATK Funcotator注释数据库下载 GATK包含两组预打包的数据源，允许在无需 (大量)额外配置的情况下使用Functator。这些数据源包分别对应于胚系和体细胞突变情形。
广义地说，如果你有胚系VCF，胚系数据源就是你想要开始使用的。相反，如果您有一个体细胞VCF，那么体细胞数据源就是您想要开始使用的数据源。
注释数据库下载方式 (2选1)：
1.ftp://gsapubftp-anonymous@ftp.broadinstitute.org/bundle/funcotator/
2.gs://broad-public-datasets/funcotator/
gsutil ls gs://broad-public-datasets/gsutil ls gs://broad-public-datasets/funcotator/gsutil -m cp -r gs://broad-public-datasets/funcotator/funcotator_dataSources.v1.7.20200521g/ \ funcotator/# annotation# --output-file-format The output file format. Either VCF, MAF, or SEG. # Please note that MAF output for germline use case VCFs is unsupported. # SEG will generate two output files: a simple tsv and a gene list. Required. Possible values: {VCF, MAF, SEG}下载结果 (3.6 GiB)
Operation completed over 81 objects/34.4 GiB.另外，gsutil也可以使用一些bash命令，例如：du -sh
https://gatk.broadinstitute.org/hc/en-us/articles/5358824293659--Tool-Documentation-Index#Funcotatorhttps://gatk.broadinstitute.org/hc/en-us/articles/5358830252187-Funcotatorhttps://gatk.broadinstitute.org/hc/en-us/articles/360035889931-Funcotator-Information-and-Tutorial
https://mp.weixin.qq.com/s/Y-_aDGTNUZVW-QMbhZDHXAhttps://qiita.com/kojix2/items/738ade193e0718509c46https://bioconductor.org/packages/release/bioc/html/maftools.html往期精品(点击图片直达文字对应教程)
机器学习
后台回复“生信宝典福利第一波”或点击阅读原文获取教程合集