歡迎關(guān)注"生信修煉手冊"
除了filter功能外,snpsift還提供了許多實用的功能,本篇文章主要看下注釋相關(guān)的幾個功能。
1. annotate
annotate命令對輸入的VCF文件進行注釋,需要一個數(shù)據(jù)庫的VCF文件,比如dbsnp等,通過在數(shù)據(jù)庫中查找,將數(shù)據(jù)庫VCF文件中的ID和INFO字段信息輸出出來, 用法如下
java -jar SnpSift.jar annotate dbSnp132.vcf variants.vcf > variants_annotated.vcfdbSnp132.vcf就是數(shù)據(jù)庫對應(yīng)的VCF文件,variants.vcf就是需要注釋的VCF文件,輸出結(jié)果通過重定向符號保存到variants_annotated.vcf文件中。
variants.vcf的內(nèi)容如下
#CHROM ?POS ID REF ?ALT ?QUAL FILTER ?INFO
22 16157571 . T G 0.0 FAIL NS=53
22 16346045 . T C 0.0 FAIL NS=244
22 16350245 . C A 0.0 FAIL NS=192
22 17054103 . G A 0.0 PASS NS=404
22 17071906 . A T 0.0 PASS NS=464
22 17072347 . C T 0.0 PASS NS=464
22 17072394 . C G 0.0 PASS NS=463
22 17072411 . G T 0.0 PASS NS=464默認情況下,會將數(shù)據(jù)庫中的ID和INFO兩個字段的信息都注釋上去,輸出結(jié)果如下
#CHROM POS ID REF ALT QUAL FILTER INFO
22 16157571 . T G 0.0 FAIL NS=53
22 16346045 rs56234788 T C 0.0 FAIL NS=244;RSPOS=16346045;GMAF=0.162248628884826;dbSNPBuildID=129;SSR=0;SAO=0;VP=050100000000000100000100;WGT=0;VC=SNV;SLO;GNO
22 16350245 rs2905295 C A 0.0 FAIL NS=192;RSPOS=16350245;GMAF=0.230804387568556;dbSNPBuildID=101;SSR=1;SAO=0;VP=050000000000000100000140;WGT=0;VC=SNV;GNO
22 17054103 rs4008588 G A 0.0 PASS NS=404;RSPOS=17054103;GMAF=0.123400365630713;dbSNPBuildID=108;SSR=0;SAO=0;VP=050100000000070010000100;WGT=0;VC=SNV;SLO;VLD;G5A;G5;KGPilot123如果你只需要ID或者INFO中的一個,可以通過參數(shù)-id和-info選擇你需要的注釋信息,示例如下
java -jar SnpSift.jar annotate -id dbSnp132.vcf variants.vcf > variants_annotated.vcf添加了-id參數(shù)之后,就只會注釋ID字段的信息了。
2. varType
varType命令會在INFO字段中加入突變類型的信息,共有以下5種突變類型
SNP
MNP
INS
DEL
-
MIXED
如果輸入的VCF文件中只包含了一個樣本的突變信息,這個命令還會給出HOM/HET的信息,用法如下
java -jar SnpSift.jar varType ?test.vcf | grep -v "^#" | head
20 10469 . ?C ?G ?100.0 PASS ?SNP;HOM ?GT:AP 0|0:0.075,0.060
20 10492 . ?C ?T ?100.0 PASS ?SNP;HET ?GT:AP 0|1:0.180,0.345
20 10575 . ?C ?CG 100.0 PASS ?DEL;HET ?GT:AP 0|1:0.000,0.000
20 10611 . ?CG C ?100.0 PASS ?INS;HET ?GT:AP 0|1:0.000,0.010
20 10618 . ?GT TA 100.0 PASS ?MNP;HET ?GT:AP 0|1:0.020,0.0303. gwasCat
gwasCat命令使用GWAS catalog數(shù)據(jù)庫對輸入的VCF文件進行注釋,用法如下
java -jar SnpSift.jar gwasCat gwascatalog.txt test.vcf | tee test.gwas.vcfGWAS Catalog數(shù)據(jù)庫存儲了相關(guān)疾病和性狀的SNP關(guān)聯(lián)分析結(jié)果,輸出結(jié)果如下
1 1005806 rs3934834 ? C ? T ? . ? PASS ? ?AF=0.091;GWASCAT=Body_mass_index ? ?
1 2069172 rs425277 ? ?C ? T ? . ? PASS ? ?AF=0.400;GWASCAT=Height
1 2069681 rs3753242 ? C ? T ? . ? PASS ? ?AF=0.211;GWASCAT=Reasoning ?
1 2392648在輸出結(jié)果中,會給出患病的風險率和相關(guān)疾病或性狀的描述信息
4. dbnsfp
dbnfsp 是一個綜合性的數(shù)據(jù)庫,收錄了SIFT, Polyphen2, LRT, MutationTaster等多種算法預(yù)測的基因組變異對功能的影響。通過snpsfit的dbnsfp命令,可以用這個數(shù)據(jù)庫對變異位點進行注釋,用法如下
java -jar SnpSift.jar dbnsfp -v myFile.vcf > myFile.annotated.vcf第一次運行時,會自動下載dbnsfp對應(yīng)的數(shù)據(jù)庫文件。自動下載的數(shù)據(jù)庫文件是snpSift官方整理的,如果你發(fā)現(xiàn)不是最新版,也可以自己去dbnsfp的官網(wǎng)下載最新版,然后自己整理,命令如下
wget http://dbnsfp.houstonbioinformatics.org/dbNSFPzip/dbNSFP2.9.zip
unzip dbNSFP2.9.zip
(head -n 1 dbNSFP2.9_variant.chr1 ; cat dbNSFP2.9_variant.chr* | grep -v "^#" ) > dbNSFP2.9.txt
bgzip dbNSFP2.9.txt
tabix -s 1 -b 2 -e 2 dbNSFP2.9.txt.gz5. phastCons
phastCons是一款計算保守區(qū)域的軟件,通過snpsift的phastCons命令可以對變異位點進行保守區(qū)域的注釋。
第一步下載數(shù)據(jù)庫,命令如下
mkdir -p ~/snpEff/db/phastCons
cd ~/snpEff/db/phastCons
wget http://hgdownload.soe.ucsc.edu/goldenPath/hg19/phastCons100way/hg19.100way.phastCons/chr1.phastCons100way.wigFix.gz
...
wget http://hgdownload.soe.ucsc.edu/goldenPath/hg19/phastCons100way/hg19.100way.phastCons/chrY.phastCons100way.wigFix.gz第二步,進行注釋,命令如下
java -jar SnpSift.jar phastCons ~/snpEff/db/phastCons file.vcf > file.phastCons.vcf6. geneSets
采用MSigDb數(shù)據(jù)庫提供的基因集注釋信息,對突變位點所在基因進行注釋,用法如下
java -jar SnpSift.jar geneSets -v db/msigDb/msigdb.v3.1.symbols.gmt test.ann.vcf > test.eff.geneSets.vcf輸出結(jié)果會在INFO中新增一個MSigDb注釋信息,示例如下
MSigDb=ACEVEDO_METHYLATED_IN_LIVER_CANCER_DN,CHR1P36,KEGG_OLFACTORY_TRANSDUCTION,REACTOME_GPCR_DOWNSTREAM_SIGNALING,REACTOME_OLFACTORY_SIGNALING_PATHWAY,REACTOME_SIGNALING_BY_GPCR掃描關(guān)注微信號,更多精彩內(nèi)容等著你!