距离准备工作笔记过去了近一周,我的准备工作仍然没能完成。
这个过程中我基本上把前辈的笔记浏览了一遍,把B站的视频也看了一遍;但是由于网络原因我无法重现相关的内容。
这些都不是目前最重要的看过的这些内容帮助我理解了我手头的数据的相关信息有助于后面的工作,重要的是下游分析,因为我手头上已经有了一个上游分析的结果。只不过缺点是没有经过原始数据的清洗过程,不知道我得到的数据进一步分析会不会有别的影响。
上游的分析最终得到变异的VCF文件进而可以进行下游的分析。
曾老师在视频的最后讲到内容涉及到下有分析的是更高级的课程所涉及的,他的博客曾经也有相关的内容。
链接:https://www.jianshu.com/p/49d035b121b8
VCF下游分析
主要是:注释和过滤
注释
VEP:http://grch37.ensembl.org/info/docs/tools/vep/index.html

snpEFF:http://snpeff.sourceforge.net/

ANNOVAR:https://doc-openbio.readthedocs.io/projects/annovar/en/latest/

以下是曾老师之前的笔记:
1.Annovar使用记录 (http://www.bio-info-trainee.com/641.html)
2.用annovar对snp进行注释 (http://www.bio-info-trainee.com/441.html)
3.对感兴趣的基因call variation(http://www.bio-info-trainee.com/2013.html)
4.WES(六)用annovar注释(http://www.bio-info-trainee.com/1158.html)
通过搜索我还找到了些其他人的介绍:
https://blog.goldenhelix.com/the-sate-of-variant-annotation-a-comparison-of-annovar-snpeff-and-vep/
突变注释工具SnpEff,Annovar,VEP,oncotator比较分析(https://www.jianshu.com/p/6284f57664b9)
博客是金螺旋(Golden Helix)的人所写,另一个是引用者简书作者生信杂谈,他讲到金螺旋开发了一个variants注释工具Varseq。自诩Varseq比其他的分析方法更好,之后很多人问Annovar的作者Kai Wang对这篇博客的意见,然后Kai Wang在博客下的评论里说:金螺旋用这篇博客攻击Annovar/VEP/snpEff是不合适的。并且给研究者们发了大量邮件来给这篇博客打广告,这可能给许多科学家造成误导。然后科学家们对博客中对Annovar解释不当的地方一一反驳,最后说Annovar注释的没什么不对,最多就是有些地方不够精确。
反正具体谁好用当然是数据说话了,我没有通过实例验证,不过从博文中的比较图可以看到三个软件对splicing variant的注释结果差异较大,很大一部分原因是由于对同一个variant的描述使用了不同的术语所造成的。从作者的韦恩图中看到SNPeff略胜一筹。

不过不同的数据选用的工具自然有不同的选择。我手头上有一个复杂疾病的靶向捕获测序的结果,想筛选出调控区的SNP可能的功能影响。因此需要进一步的功能注释。
首先这个数据测序完成与2013年,当时的主流数据库可能是hg19,目前我需要更新版本到hg38。那么首先第一步应该进行坐标转换。其实这个在学习ChIPseq的分析的时候有学到过,当时曾老师演示了UCSC 的 liftOver将果蝇的基因组dm3转换为dm6,但是只看过一遍,当时没有完整的操作也没有记录的很清楚,我给忘了怎么操作。然后朋友比较熟悉CrossMap.py(http://crossmap.sourceforge.net/)给我用这个转换了。我在询问的过程中才想起来之前学习过。

使用代码:
CrossMap.py vcf hg19ToHg38.over.chain.gz pooling_variants_all_variants.hg19.vcf ~/db/gencode/GRCh38/GRCh38.primary_assembly.genome.fa pooling_variants_all_variants.hg19-hg38.vcf

这个工具核心就是chain文件,有了这个就可以做对应想要的版本转换了。
第一步任务完成得到了新的VCF文件。
但是进行注释需要相应的软件和数据库我也没有,朋友那里有Annovar注释相关的软件但是他们的数据库还是hg19的所以没办法得到相应的注释
听外显子视频的最后一章节曾老师讲到使用annovar可以一次性注释多个数据库,曾老师还讲到例子中利用 SIFT ,PolyPhen-2 以及 PROVEAN 软件, 预测 SNV 对蛋白质功能的影响程度,仅当 3 种软件均预测同一遗传变异对蛋白质 的功能影响较大时,才认定该遗传变异具有高危害性。利用 PROVEAN 软件预测 Indel 对蛋白质功能的影响。但是其实dbNSFP数据库,就注释了这些变异对蛋白功能的影响。刚好我问了朋友他以前下载过dbNSFP,但是用了如下代码完成了注释后
table_annovar.pl pooling_variants_all_variants.hg19.avinput ~/db/annovar_humandb/ -buildver hg19 -out hg19_anno -remove -protocol ensGene,gnomad211_exome,esp6500siv2_all,EAS.sites.2015_08,ALL.sites.2015_08,dbnsfp35a -operation g,f,f,f,f,f -nastring . -csvout -polish --thread 8
得到的数据疾病上外显子上的突变位点都有了相应的打分情况。但是我想要看的调控区就比较少。
后面他提议说用CADDhttps://cadd.gs.washington.edu/score
来注释,这个可以提供全基因组的致病性预测,然后因为这个软件支持的检测版本比较高,刚好之前转换了坐标的信息就有了用武之地。

参考资料:https://github.com/kircherlab/CADD-scripts/blob/master/README.md
打分越高提示致病性越高,感觉还是挺有用的,这个结果再进一步分析应该能达到老板想要的目的。在这个过程中复习了以前的知识,如果不复习永远没有进步。