VCF可以通過不同的上游分析獲得剧腻,在格式上又有一些不同裸诽,盲目轉(zhuǎn)化會(huì)帶來很多問題持舆。查看VCF文件內(nèi)容可知是通過GATK上游處理得到的向叉。
less raw.vcf
##fileformat=VCFv4.1
##FILTER=<ID=LowQual,Description="Low quality">
##FORMAT=<ID=AD,Number=.,Type=Integer,Description="Allelic depths for the ref and alt alleles in the order listed">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Approximate read depth (reads with MQ=255 or with bad mates are filtered)">
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Genotype Quality">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=PL,Number=G,Type=Integer,Description="Normalized, Phred-scaled likelihoods for genotypes as defined in the VCF specification">
##GATKCommandLine.HaplotypeCaller=
為了對變異結(jié)果進(jìn)一步進(jìn)行可視化分析,查找了一下網(wǎng)上的相關(guān)資料吠勘,可以通過annovar軟件進(jìn)行轉(zhuǎn)化妆艘,然后通過maftools R包進(jìn)行可視化分析。確定了整體的分析策略之后看幼,先安裝annovar軟件。
1.安裝annovar
官網(wǎng)下載:http://annovar.openbioinformatics.org/en/latest/user-guide/download/
注冊 網(wǎng)站發(fā)送郵件 直接郵件下載 拷貝到服務(wù)器 (需使用機(jī)構(gòu)郵箱)幌陕。
拷貝到服務(wù)器诵姜,在相應(yīng)的文件夾下進(jìn)行解壓即可。
tar -zxvf annovar.latest.tar.gz
2.VCF轉(zhuǎn)換為annovar格式
### 批量轉(zhuǎn)換
for i in *vcf; do perl ~/biosoft/annovar/convert2annovar.pl -format vcf4 $i -outfile $i.annovar -includeinfo; done
3.VCF變異根據(jù)基因進(jìn)行注釋
for i in *vcf; do \
perl ~/biosoft/annovar/table_annovar.pl $i.annovar ~/path/to/humandb/ -buildver hg19
-out variants --otherinfo -remove -protocol ensGene -operation g -nastring NA
這個(gè)時(shí)候會(huì)遇到一個(gè)報(bào)錯(cuò):
hg19_ensGeneMrna.fa does not exist. Please use 'annotate_variation.pl --downdb ensGene ~/biosoft/annovar/humandb' download the database
按照提示下載相應(yīng)的注釋參考數(shù)據(jù)庫搏熄,提示相應(yīng)的ref文件已經(jīng)無法下載棚唆。重新構(gòu)建ref
perl ./annotate_variation.pl --buildver hg19 --downdb seq ~/biosoft/annovar/humandb/hg19_seq
perl ./retrieve_seq_from_fasta.pl /home/data/ssy066/biosoft/annovar/humandb/hg19_ensGene.txt -seqdir ~/biosoft/annovar/humandb/hg19_seq -format ensGene -outfile ~/biosoft/annovar/humandb/hg19_ensGeneMrna.fa
重新進(jìn)行注釋。
4.轉(zhuǎn)化為maf格式文件
目標(biāo)文件夾下會(huì)出現(xiàn)后綴為multianno.txt的文件心例,導(dǎo)入R語言中轉(zhuǎn)化為maf文件宵凌,進(jìn)行可視化即可。止后、
library(maftools)
var.annovar.maf = annovarToMaf(annovar = "~/BN4027P6.raw.vcf.anno.hg19_multianno.txt",
Center = 'CSI-NUS', refBuild = 'hg19', MAFobj = T, table = 'ensGene')
plotmafSummary(maf = var.annovar.maf, rmOutlier = TRUE, addStat = 'median', dashboard = TRUE, titvRaw = FALSE)
參考內(nèi)容
1.maftools : Summarize, Analyze and Visualize MAF Files (bioconductor.org)