基因組上有很多變異，不同變異的效果不同。
除了wet-lab，現(xiàn)在也有很多軟件預(yù)測。接下來，我打算用SIFT4G、Polyphen-2、SnpEff和provean預(yù)測。
因?yàn)槔习逑敕ǜ淖?，剩下的兩個(gè)軟件就···不弄了
0. vcf文件準(zhǔn)備
msa2vcf可以把多序列比對轉(zhuǎn)為vcf文件。較其他工具，它可以處理gap，但無法轉(zhuǎn)換為reference genome上的位置，需要自己寫程序轉(zhuǎn)換。
1. SnpEff
1.1 檢查database是否正確
軟件里存了很多database，大部分物種都在，用之前可以先檢查一下感興趣物種在不在。
# 列出已構(gòu)建的database
java -jar snpEff.jar databases
不過，即使有，也需要檢查database創(chuàng)建是否有誤，比如genetic code是否正確
java -Xmx4g -jar snpEff.jar -v Lactobacillus_plantarum_gca_001005805
編碼方式居然是standard code，事實(shí)上應(yīng)該是bacteria code ，也就是genetic code 11。
1.2 build自己的database
具體參考：https://pcingola.github.io/SnpEff/se_buildingdb/
如何build，要根據(jù)自己手邊的數(shù)據(jù)。接下來，我采取的是genbank方式。
分三步：1）下載genbank數(shù)據(jù)，2）配置文件，3）run
1.2.1 下載genbank數(shù)據(jù)
拿著chromosome或scaffold的編號(hào)，去ncbi的nucleotide庫中搜，下載數(shù)據(jù)。
下完之后，將所有chromosome的cat起來，命名為genes.gbk。
一定得叫g(shù)enes.gbk
1.2.2 配置文件
在snpEff.config中添加配置信息：
# Lactobacillus plantarum ps128
Lactobacillus_plantarum_ps128.genome : Lactobacillus plantarum
Lactobacillus_plantarum_ps128.chromosomes : NZ_LBHS01000003.1, NZ_LBHS01000004.1, NZ_LBHS01000008.1, NZ_LBHS01000009.1, NZ_LBHS01000010.1, NZ_LBHS01000011.1, NZ_LBHS01000005.1, NZ_LBHS01000006.1, NZ_LBHS01000007.1, NZ_LBHS01000001.1, NZ_LBHS01000002.1
Lactobacillus_plantarum_ps128.NZ_LBHS01000003.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000004.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000008.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000009.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000010.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000011.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000005.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000006.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000007.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000001.1.codonTable : Bacterial_and_Plant_Plastid
Lactobacillus_plantarum_ps128.NZ_LBHS01000002.1.codonTable : Bacterial_and_Plant_Plastid
每個(gè)chromosome都要弄，chromosome的id不能搞錯(cuò)。比如NZ_LBHS01000003.1不能寫成NZ_LBHS01000003。
在snpEff.config目錄下做以下操作:
mkdir -p data/Lactobacillus_plantarum_ps128
mv genes.gbk data/Lactobacillus_plantarum_ps128/
1.2.3.run
java -jar snpEff.jar build -genbank -v Lactobacillus_plantarum_ps128
檢查日志文件，是否有報(bào)錯(cuò)。
如果沒有什么問題，在data/Lactobacillus_plantarum_ps128文件夾下會(huì)出現(xiàn)snpEffectPredictor.bin文件，創(chuàng)建成功。
1.3 預(yù)測
java -Xmx4g -jar snpEff.jar -v Lactobacillus_plantarum_ps128 -ud 0 0_align.formatted.vcf | ./scripts/vcfInfoOnePerLine.pl> 0_align.prediction.vcf
-Xmx4g：給程序分配4G內(nèi)存，視自己情況而定
-ud: 設(shè)置upstream和downstream interval size。如果變異與gene A的距離小于interval size，預(yù)測時(shí)會(huì)包括對gene A的影響。
2. Provean
注意點(diǎn)1：
與SnpEff相比，Provean只能預(yù)測部分蛋白質(zhì)編碼基因的variation。建議Provean放在SnpEff之后，這樣子可以借用SnpEff的結(jié)果作為Provean的輸入，比較簡單。
提取SnpEff結(jié)果中的非intergenic region的變異，作為Provean輸入數(shù)據(jù)。
注意點(diǎn)2：
因?yàn)檠芯课锓N是細(xì)菌，所以只能選用PROVEAN Protein工具，缺點(diǎn)是不能批量，每次只能預(yù)測一個(gè)蛋白質(zhì)序列。
2.1 預(yù)測
為了方便整理，我為每個(gè)蛋白質(zhì)都建了一個(gè)文件夾，里面放了兩個(gè)文件：
- sequence.fa：存儲(chǔ)一條蛋白質(zhì)序列
- var：存放該蛋白質(zhì)上的變異（必須是HGVS表示，SnpEff結(jié)果文件中有）
provean.sh -q sequence.fa -v var
SnpEff中的HGVS表示是三字母氨基酸縮寫，需要轉(zhuǎn)為單字母氨基酸縮寫。
