就基因家族工作做一簡單介紹
基本思路
數(shù)據(jù)準(zhǔn)備
確定好研究的基因家族后（比如：NBS廓鞠，MADS-box etc.）帚稠，下面就可以下載相關(guān)數(shù)據(jù)。
- 所研究物種的基因組序列床佳； genome.fa
- 所研究物種蛋白序列滋早；pep.fa
- 所研究物種gff文件
- 目標(biāo)基因家族的隱馬科夫模型
- or RefSeq 對(duì)應(yīng)基因家族的蛋白序列
對(duì)應(yīng)基因組信息可根據(jù)發(fā)表文章中提供的路徑進(jìn)行下載即可；
對(duì)于隱馬科夫模型可以從Pfam進(jìn)行下載（比如：NBS）
點(diǎn)擊BROWSE后砌们，進(jìn)行搜索NBS即可杆麸；本次下載PF00931
下載RefSeq對(duì)應(yīng)抗病序列如下所示搁进；
或者下載RefSeq植物所有信息，根據(jù)注釋信息選取相關(guān)序列
見基因功能注釋
鑒定家族成員
-
基于保守結(jié)構(gòu)模型進(jìn)行鑒定
所用工具為hmmer
hmmbuild/hmmsearch/hmmscan/hmmalign 這幾個(gè)功能是主要用于蛋白質(zhì)結(jié)構(gòu)與分析和注釋的hmmer中小工具
軟件安裝
conda install -c bioconda hmmer=3.3
在鑒定基因家族時(shí)昔头，常用到的工具是hmmsearch饼问，里面常用的算法有三種。一般我們使用--cut_tc算法對(duì)隱馬可夫模型進(jìn)行搜索揭斧，tc算法是使用pfam提供的hmm文件中trusted cutoof的值進(jìn)行篩選莱革，相對(duì)比較可靠。
以2015年發(fā)表文章Identification and distribution of the NBS-LRR gene family in the Cassava genome為例讹开，解析其步驟盅视。
分析流程如下圖所示
首先根據(jù)NB-ARC 模型，利用hmmersearch 篩選質(zhì)量較高的基因（E-value?<?1 × 10?20）旦万。 而后使用clustalw2對(duì)高質(zhì)量的序列進(jìn)行多序列比對(duì)闹击，多序列比對(duì)后，對(duì)這些置信的序列進(jìn)行隱馬可夫模型的構(gòu)建（使用hmmbuild）纸型，最后用新建的模型拇砰，重新對(duì)所需基因進(jìn)行篩選。
第一步 篩選高質(zhì)量NBS基因
## hmmersearh進(jìn)行搜索
hmmsearch --cut_tc --domtblout NBS-ABC.out NBS-ARC.hmm \
Arabidopsis_thaliana.TAIR10.pep.all.fa
# --cut_tc: use profile's TC trusted cutoffs to set all thresholding
# --domtblout <f> : save parseable table of per-domain hits to file <f>
## 根據(jù)e-value 篩選高質(zhì)量基因
grep -v "#" NBS-ABC.out|awk '($7 + 0) < 1E-20'|cut -f1 -d " "|sort -u > NBS-ARC_qua_id.txt
第二步：使用clustalw2進(jìn)行多序列比較,構(gòu)建新模型
比對(duì)具體步驟請(qǐng)看序列比對(duì)和構(gòu)建進(jìn)化樹（clustalw和phylip）
hmmbuild NBS-ARC.second.out NBS-ARC_qua.aln
hmmsearch --cut_tc --domtblout NBS-ARC.second.out \
NBS-ARC_qua.hmm ../Arabidopsis_thaliana.TAIR10.pep.all.fa
第三步：再次進(jìn)行過濾
grep -v "#" NBS-ABC.second.out|awk '($7 + 0) < 1E-03' | cut -f1 -d " "|sort -u >final.NBS.list
-
基于同源注釋
簡而言之狰腌，利用已知的NBS基因作為ref除破，進(jìn)行blastp比對(duì)
如果要取得可信度較高的結(jié)果，可講上述兩種方法得到的基因取交集即可
com m -12 gene1 gene1 >common.list
