上一篇文章《肿瘤变异数据分析和可视化工具maftools：突变数据下载和可视化》主要以TCGA-LUAD为例介绍突变数据和临床数据的下载、处理以及简单的可视化。这篇文章更详细介绍可以利用maftools对肿瘤MAF格式的突变数据做哪些分析。
具体的数据下载和处理方法这里就不再赘述了，请移步上篇文章：《肿瘤变异数据分析和可视化工具maftools：突变数据下载和可视化》。这篇文章里，读入的临床数据终于可以派上用场了。
library(maftools)
luad <- read.maf(maf="TCGA.LUAD.maf", clinicalData="clinical.tsv")本文接上次的内容：《肿瘤变异数据分析和可视化工具maftools：安装和文件格式要求》，本文以TCGA肺腺癌（LUAD）的数据为例介绍突变数据下载和可视化。
在TCGA官网下载TCGA-LUAD项目中mutect2输出的MAF格式的突变以及临床信息：
$ ls -lhrt
total 50M
-rw-r--r-- 1 xiaofei xiaofei 50M May 14 22:44 TCGA.LUAD.maf.gz
-rw-r--r-- 1 xiaofei xiaofei 157K May 15 00:12 clinical.tsvMaftools是一款可以对MAF格式（Mutation Annotation Format）的变异数据进行统计、分析和可视化的R包。除了可以对TCGA来源的MAF文件以外，其他任何变异数据只要是MAF格式都可以使用这款工具进行分析。
Maftools包可主要概括为可视化和分析两大模块，流程和使用方法很简单：通过read.maf读入MAF文件（或者经过格式转换）得到MAF对象，然后将对象传递给对应的分析或者可视化函数就行了。主要模块、函数和主要的分析和可视化功能见下图：
VAF的全称是Variant Allele Frequency（变异等位基因频率）或Variant Allele Fraction（变异等位基因分数）。简单来说就是在基因组某个位点支持alternate/mutant allele的reads覆盖深度占这个位点总reads覆盖深度的比例。以VCF文件中的字段为例，其中DP代表Total Depth，AD代表Allele Depth，因此VAF的计算就是：
$$VAF = \frac{Allele\ Depth}{Total\ Depth} = \frac{AD}{DP}$$
VAF用得比较多的地方是在二倍体germline的genotyping中，杂合位点的VAF在高深度（比如depth>80）情况下应该接近50%；如果VAF接近0.25/0.75说明基因组上可能还有另一份拷贝。另一个应用场景就是癌症基因组的somatic genotyping。肿瘤组织、cfDNA、ctDNA、CTC genotyping的结果中会包含正常的allele（与正常体细胞一致）以及突变的allele，其中突变allele的所占的比例就是VAF。VAF可以用于推断肿瘤的异质性和肿瘤纯度，此外VAF的高低可能会影响癌症的预后。
线性回归问题即已知一系列样本的自变量和因变量的值，求解以下方程中的各θj：
$$h_θ(x) = θ_0 + θ_1x_1 + θ_2x_2 + \cdots + θ_nx_n$$
设样本数量为m，评估拟合的直线与实际样本之间差异的代价函数（Cost Function）为：
$$J(θ_0,θ_1,\cdots,θ_n) = \frac{1}{2m}\sum_{i=1}^{m}(h_θ(x^{(i)})-y^{(i)})^2$$
因此，寻找最佳拟合的线性回归模型则转化为求解该代价函数的最小值，常用方法为最小二乘法（Least Squares Method）和梯度下降法（Gradient Descent Method）。
癌症疫苗靶向肿瘤细胞的抗原可以大致分为两类：肿瘤相关性抗原（tumor-associatied self-antigen）以及肿瘤特异性抗原（tumor-specific antigen）。前者指的是在正常体细胞中也存在，但肿瘤细胞中异常高表达的抗原；后者指的是因肿瘤特异性突变而产生的新抗原/新表位（neoantigen/neoepitope）。新抗原相比肿瘤相关性抗原拥有更强的特异性因此副作用更低，并且不受限于胸腺的中枢耐受。通过高通量测序可以获取大量的肿瘤特异性突变，基于这些突变预测新抗原在癌症的个性化免疫治疗方面有很广阔的应用前景。
I类主要组织相容性复合体（MHC class I）抗原的处理和呈递过程：（1）肿瘤特异性突变产生的突变蛋白会被蛋白酶体降解为8~11aa的肽段；（2）这些肽段被抗原加工相关转运体（TAP）转运进入内质网腔；（3）与新合成的MHC-I结合；（4）最终通过高尔基体转运至细胞膜被CD8+ T细胞识别。
（本文为个人见解，仅供参考。内容主要针对在企业从事生物信息学分析的工作者，不讨论学校和其他科研单位。）
生物信息分析的工作被很多人视为高大上——不需要做实验也能发文章、作为交叉学科门槛比较高、因为高门槛或者沾了IT的光，所以工资待遇相比生物专业其他方向高出不少。但是入行几年或多或少内心都有一些迷茫：这份工作有没有前途、能否做一辈子？自己又如何提高和改变以适应未来的需求变化？我想在此之前应当先思考一下对这份工作是否有正确的认识。
作为一个半路出家的生物信息分析工作者在刚接触这行时真的觉得门槛挺高的：要学编程，要掌握Linux操作系统，还有眼花缭乱的测序数据、文件格式和对应的分析软件（还很难安装）。现在回过头来看，这可能只是一厢情愿的想法。我认为从技术水平上可以分成几个不同的阶段（有别于刘小乐教授对“生物信息学研究”的分级）：
最近听闻前公司和另一家公司A因为一个项目中二代测序数据PCR duplication的问题在扯皮。前公司比较两对Paired-end reads的第8~113bp区间内碱基是否完全相同（PE150）；而A公司把1~150bp都完全相同的才当作PCR duplication。最终参数的分歧使得两家公司的统计结果相差接近一倍，关乎到数据量是否的满足合同要求。
而前公司这部分QC流程的参数就是我设置的。选择截取reads中间的区域进行比较的原因很简单，因为Illumina测序前面几个碱基和末端的区域通常测序质量较差，测序错误较多。如果要求两对reads完全一致，相当于600bp（150*2*2）不能出现任何测序错误。具体为什么设置成8~113bp我并没有详细去分析和论证，只是因为当时测序都外包给N公司做，而N公司内部QC流程去PCR duplication设置的范围设置的就是这个，我也就顺理成章用了同样的参数。
那么怎样去除二代测序数据中的PCR duplication才是科学的姿势呢？虽然对很多人来说应该是一个无聊的论题，但是我还挺好奇，本文姑且在这里分析一二。