基因家族生信分析

1、基因家族生信分析一、 什么是基因家族概念：是來源于同一個祖先，有一個基因通過基因重復而產(chǎn)生兩個或更多的拷貝而構成的一組基因，他們在結構和功能上具有明顯的相似性，編碼相似的蛋白質(zhì)產(chǎn)物。劃分：按功能劃分：把一些功能類似的基因聚類，形成一個家族。按照序列相似程度劃分:一般將同源的基因放在一起認為是一個家族。1.常見基因家族：WRKY基因家族：是植物前十大蛋白質(zhì)基因家族之一，大量研究表明，WRKY基因家族的許多成員參與調(diào)控植物的生長發(fā)育，形態(tài)建成與抗病蟲。NBS-LRR抗病基因家族：是植物中最大類抗病基因家族之一。MADS-BOX基因家族：是植物體內(nèi)的重要轉(zhuǎn)錄因子，它們廣泛地調(diào)控著植物的生長、發(fā)育和生

2、殖等過程。在植物中參與花器官的發(fā)育，開花時間的調(diào)節(jié)，在果實，根，莖，葉的發(fā)育中都起著重要的作用。熱激蛋白70家族（HSP70）是一類在植物中高度保守的分子伴侶蛋白，在細胞中協(xié)助蛋白質(zhì)正確折疊。二、 基因家族分析流程：l 利用蛋白保守域結構提取號在Pfam數(shù)據(jù)庫提取其隱馬爾科夫模型矩陣文件（*.hmm）l 在數(shù)據(jù)庫（Ensemble 、JGI、NVBI）下載你所需要的物種的基因組數(shù)據(jù)（*.fa,*.gff）l 在虛擬機中Bio-Linux中的hummsearch程序，用隱馬爾科夫模型矩陣文件在蛋白序列文件中搜索含有該保守結構域的蛋白l 將蛋白序列導入MEGA軟件構建進化樹（可以闡明成員之間系統(tǒng)進

3、化關系，從進化關系上揭示其多樣性）l 利用MEME搜索蛋白質(zhì)的保守結構域利用MEME搜索基因家族成員的motif可以揭示基因家族在物種內(nèi)的多樣化及其功能，如果他們都含有相同的motif表明其功能具有相似性，如果部分家族成員含有其他不同的motif,很可能這些成員有其他特異功能，或者可以歸分為一個亞族l 繪制基因染色體位置圖從*.gff文件中抽取我們搜索到的基因位置信息， 在線繪制基因染色體位置圖通過染色體位置分布，可以了解基因主要分布字哪條染色體上，及是否能形成基因簇（被認為是通過重組與錯配促進基因交流）l 基因結構分析從gff文件中抽取基因的結構信息，繪制轉(zhuǎn)錄本結構圖。l 計算串聯(lián)重復基因的

4、Ka,Ks1. 首先將篩選到的基因的cds序列進行多序列對比，篩選identity > 75%，tength大于對比的兩條序列中較長的那條的長度的75%，將篩選到的基因分別用clustalw進行比對，比對結果導入KsKs_Calculster計算Ka，Ks、Ka/ks比，計算核苷酸的非同義替代（ka）與核苷酸的同義替代（ks）的平均速率。2. Ka/ks比值<1表明：通過純化選擇降低了氨基酸變化的速率；比值=1表示中性選擇；比值>1,表明這些基因可能已經(jīng)收到積極選擇，有利于適應性遺傳，這些受正向選擇的基因?qū)⒆鳛橐院蟮难芯恐攸c。 軟件的安裝從圖片中獲得進入NCBI-blast官

5、網(wǎng)復制blast-linux版本的鏈接在Linux終端1. blast的安裝#wget blast鏈接#tar xvfz 文件名 解壓縮文件# mv 解壓縮文件 /root/local/app# mv 解壓縮文件 blast# vi .bashrc#在最后一行添加export $PATH=/root/local/app/blast/bin:$PATH 并保存退出#source .bashrc 運行 #blastp -version 查看是否安裝成功。2hummer的安裝#yum install -y wget /安裝wget#wget hmmer源碼鏈接 #tar -zxvf hmmer-3.

6、2.1# vi .bashrc#(在最末端添加的語句) PATH=$PATH:/biosoft/ /binaries#yum install -y gcc #./configure#make #make check#make install#which hmmsearch 查看是否安裝成功。3.perl的安裝#wget 源代碼鏈接# tar xvfz perl-5.28.1.tar.gz 解壓縮#./configure#make#make install 安裝完成。3.bioperl 的安裝#wget -O - https:/install.perlbrew.pl | bash#perlbre

7、w install-cpanm#/root/perl5/perlbrew/bin/cpanm Bio:Perl三、 具體操作：1.保守域結構分析下載蛋白保守結構域文件、cds、cDNA、gff注釋文件和隱馬爾科夫矩陣模型。以擬南芥為例：下載完成后，需要將文件傳到Linux系統(tǒng)上進行分析： 打開虛擬機輸入ip a將虛擬機IP連接到Xshell上，在Xshell上進行操作，將文件通過xftp（同樣需要連接IP）傳到Linux系統(tǒng)上，然后進行解壓。 (一個命令解壓多個文件:gunzip*.gz)接下來用hummsearch尋找含有該蛋白保守結構域的蛋白及核酸序列安裝hummsearchyum ins

8、tall -y wget /安裝wget#wget hmmer源碼鏈接 #tar -zxvf hmmer-3.2.1# vi .bashrc#(在最末端添加的語句) PATH=$PATH:/biosoft/ /binaries#yum install -y gcc #./configure#make #make check#make install#which hmmsearch 查看是否安裝成功。解壓文件移動到APP目錄下面在app目錄下面新建文件夾 mkdir hmmer將hmmer-3.2.1移到hmmer目錄下面 mmove -v c:/hmmer-3.2.1 c:/hmmer刪除安裝

9、包打開文字編輯器vi /.bashrc 在文字編輯器里最后一行添加以上內(nèi)容安裝好wget如果make check出現(xiàn)錯誤百度用以下方法解決出現(xiàn)complete安裝完成#source /.bashrc #which hmmsearch至此hmmer安裝完成。虛擬機操作：1.導入下載好的文件；2.hmmsearch -cut_tc domtblout NB-ARC.txt NB-ARC.hmm Arabidopsis_thaliana.TAIR10 pep.all.fa 可以用editplus打開.txt文件3. perl domain_xulie.pl 結果文件.txt 蛋白序列文件 domai

10、n.fa 1e-204.clustalw進行多序列比對，得到aln文件和dnd文件。5.hmmbuild 擬南芥特異的hmm模型文件 domain.aln6.hmmsearch cut_tc domtblout result.txt newhmm文件 蛋白質(zhì)序列文件7.在Excel中，根據(jù)特定的evalue進行篩選，并對第一列進行去重復，得到第一列去重復的id，保存為id.txt8. 用perl腳本根據(jù)id提取序列Per get_fa_by_id.pl id.txt 蛋白序列wenjain >結果輸出文件可以根據(jù)的得到的序列文件進行后續(xù)的構建進化樹、motif分析等。2.搜索基因家族成員

11、的MOTIF2.1 需要準備的文件 1. 擬南芥NBS基因蛋白質(zhì)序列 2. 蛋白保守結構域的隱馬爾科夫模型矩陣文件2.2 MOTIF的搜索使用meme軟件命令：meme nbs_pep.fa -protein -oc nbs_motif -nostatus -maxsize 600000 -mor anr -nmotifs 10 -minw 6 -maxw 50 搜索結果存放在nbs_motif文件夾中。文件夾中的eps文件可以用AI打開編輯，可以另存為png或jpg格式，也可打開網(wǎng)頁版，也可用tbtools軟件打開，下載motif在基因上的位置信息。3. 繪制基因在染色體上的位置圖3.1 需

12、要準備的文件1. 擬南芥NBS基因id2. 擬南芥的注釋文件（gff3文件）3. 擬南芥基因組長度4.1在線繪圖工具：MapGene2Chrom4.2 samtools faidx 擬南芥.dna.toplevel.fa 可得到擬南芥.文件 該文件包括 各個染色體，染色體長度。4.3 對基因的id文件在Excel中進行分列，去重復處理。4.4 使用處理過的id文件，對擬南芥的注釋文件進行篩選使用perl腳本得到基因在染色體上的位置。命令： perl get_gene_gff.pl -in1 基因的id文件 -in2 擬南芥gff3文件 -out 新文件名稱 4.5 新文件存放的是基因在染色體上

13、的位置4.6在在線文件MapGene2chrom 中，將基因在染色體上的位置信息文件復制到，input1框中，在input2中粘入samtools得到的fai文件。4.繪制轉(zhuǎn)錄本的結構圖4.1 需要準備的文件1. 擬南芥NBS基因轉(zhuǎn)錄本id（通過家族成員鑒定得到的蛋白id文件）2.擬南芥基因的注釋文件（gtf文件）3.在線繪圖工具：Gene Structure Display Server 2.0http：/gsds.cbi4.2 具體方法1. 準備gtf文件：輸入命令：gffread gff3注釋文件 -T -o 輸出文件（gtf文件）2.editplus 打開gtf文件，去除”transc

14、ript:”3.使用perl腳本提取擬南芥轉(zhuǎn)錄本結構信息：命令： perl get_gtf.pl -in1 擬南芥轉(zhuǎn)錄本id文件 -in2 gtf文件 -out 輸出文件（nbs_gtf.txt）4.通過在線繪圖工具，進行繪圖。5.篩選出串聯(lián)重復基因5.1準備文件1.擬南芥NBS基因CDS序列串聯(lián)重復基因篩選標準【（a）length of alignable sequence covers>75% of longer gene,and(b) similarity of aligned regions >75%】參考文獻：Extent of gene duplication in t

15、he genomes of Drosophila, nematode, and yeast.2.由于篩選時產(chǎn)生的文件較多，因此創(chuàng)建新的目錄：mkdir 新目錄3.用editplus 打開家族成員的id文件，對轉(zhuǎn)錄本id進行處理，使一個基因只拿一個轉(zhuǎn)錄本。4.把id復制到Excel，首先排序處理，然后進行分列，然后以第一列刪除重復值。最后將第一列和第二列進行合并。將處理好的id導入Linux。5.使用perl腳本提取cds序列：命令： perl get_fa_by_id.pl id文件 擬南芥cds序列文件 > cds.fa6.使用blast軟件篩選串聯(lián)重復基因6.1建立目標序列的數(shù)據(jù)庫：

16、makeblastdb -in cds.fa -dbtype nucl -title cds.fa6.2 進行多序列比對：blastn -query cds.fa -db cds.fa -evalue 1e-20 -outfmt 9 -out result.txt6.3 用editplus打開6.4 得到cds序列的長度，使用samtools工具建立索引：命令 ：samtools faidx cds.fa6.5 用perl腳本對result.txt進行篩選，perl KAKS_SHAIXUAN.pl -in1 cds.fa.fai -in2 result.txt -out cleanresul

17、t.txt6.6 用editplus打開，將內(nèi)容復制到Excel，在id后插入一列用公式：if （A1>B1,A1&B1,B1&A1）。然后全選，以第C列刪除重復值。并保存到新的文件中，并導入到Linux中7. 計算串聯(lián)重復基因的KaKs。7.1準備文件1.串聯(lián)重復基因的CDS序列文件7.2 將成對的串聯(lián)重復序列保存在一個文件中。、方法一：復制需要找到的序列的id，在editplus中按ctrl+F搜索，找到后復制粘貼到一個文件中。方法二：首先將成對的id保存在同一個文件，導入到Linux中，在Linux中，利用perl腳本提取序列：perl get_fa_by_id.p

18、l 新的id文件 cds.fa文件 >id1.fa7.3 計算KaKS1.計算之前需要使用CLUSTAW對序列進行比較?？色@得id1.aln。2.使用KaKs_calculator工具將id1.aln文件轉(zhuǎn)換成id1.axt文件命令：axtvenvertor id1.aln id1.axt3.計算KaKs，輸入命令： KaKs_calculstor -i id.axt -o id1_kaks.txt4.如果報錯，則把兩條序列長度保持一致。依此將所有的串聯(lián)重復基因?qū)?，進行計算。四基因家族成員的鑒定（未知隱馬爾科夫模型）1.鑒定測略在NCBI數(shù)據(jù)庫中盡量多下載幾個物種的需要鑒定的蛋白保守結構

19、域序列，以及所要研究物種的所有蛋白序列在虛擬機中本地建庫，并進行blast建庫命令：makeblastdb -in 研究物種的蛋白序列文件 -dbtype prot -title 庫名稱進行序列比對 命令： blastp -query 下載的多個物種序列文件 -db 庫名稱 -evalue 1e-10 -outfmt 6 -out 結構域.blast使用sed 命令去除表頭和結尾得到新的 new結構域. blastAwk print$1 new結構域. Blast |less 可查看打印的結果 Awk print$1 new結構域. Blast > id.txtcat id.txt|so

20、rt|uniq > idd.txt去重復Perl get_fa_by_id.pl idd.txt去重復 蛋白序列wenjain >結果輸出文件在Pfam或者NCBI的cdd中搜索檢查是否有相關蛋白結構域。之后再進行motif分析一些命令及軟件應用說明（參考一些視頻資料）hmmsearch使用說明用途：利用蛋白保守結構域的隱馬爾科夫模型搜索蛋白序列中具有該保守結構域的蛋白用法：說明：result.txt是輸出的結果文件，*.hmm在pfam數(shù)據(jù)庫下載的模型， pep.all.fa物種基因組所有蛋白序列hmmbuild使用說明用途：利用clustalw比對生成的aln文件構建蛋白保守結

21、構域的隱馬爾科夫模型用法：hmmbuild new.hmm domain.aln 說明：new.hmm 是結果文件也就是構建的蛋白保守結構域的隱馬爾科夫模型，domain.aln是clustalw比對生成的aln文件domain_xulie.pl腳本使用說明用途：提取hmmsearch搜索結果中蛋白序列中保守結構域的序列，用于構建新的物種特異的蛋白保守結構域的隱馬爾科夫模型用法：perl domain_xulie.pl（腳本不在使用目錄下要寫全路徑） hmmoutfile pep.all.fa domain.fasta E-value說明：hmmoutfile 是hmmsearch搜索結果文件

22、 pep.all.fa物種基因組所有蛋白序列 domain.fasta是結果存放文件也就是蛋白序列中保守結構域的序列，E-value 是提取序列時設定的E值get_fa_by_id.pl使用說明用途：通過ID號獲取其相應的基因或蛋白序列用法：perl perl get_fa_by_id.pl id.txt cds.fastat >id_cds.fasta 說明：id.txt是包含你的ID的文件，cds.fasta是你叢數(shù)據(jù)庫中下載的包含所有cds序列的文件，id_cds.fasta是輸出文件內(nèi)容是ID對應的序列samtools faidx 用途：提取fasta文件信息用法：samtool

23、s faidx *.fa說明：輸入文件是fasta文件，自動生成輸出目錄*.fa.fai,結果的fai文件第一列是你輸入的fasta文件的ID第二列是其序列長度Gffread 使用說明用途：將基因組注釋文件gff3轉(zhuǎn)化成基因的注釋文件gtf用法：gffread my.gff3 -T -o my.gtf說明：my.gff3是輸入文件基因組注釋文件，my.gtf是輸出文件是基因的注釋文件Get_gene_gff.pl用途：想要繪制基因的染色體位置圖必須要拿到基因在染色體上的具體信息，該腳本就是從總的gff文件抽去你需要的基因的信息如：所在染色體，起始終止位置等信息用法：Perl Get_gene_

24、gff.pl -in1 gene_id.txt -in2 my.gff3 -out gene_location.txt說明：gene_id.txt 是第一個輸入文件基因的ID文件，my.gff3是第二個輸入文件是物種基因組所有蛋白序列 ，gene_location.txt 是結果輸出文件Get_gtf.pl腳本使用說明用途：從基因注釋文件gtf文件中提取轉(zhuǎn)錄本的結構信息用法：perl get_gtf.pl -in1 id.txt -in2 gene.gtf -out structure.txt說明：id.txt 是第一個輸入文件是id文件 ，gene.gtf是第二個輸入文件是基因注釋文件gtf

25、文件 structure.txt是結果輸出文件存放著轉(zhuǎn)錄本的結構信息KaKs_shaixuan.pl 使用說明用途：在多序列比對結果文件中篩選identity大于75%，比對上的序列長度大于對比的兩條序列中最長序列的長度的75%用法：perl KaKs_shaixuan.pl -in1 cds.fai -in2 result.txt -out shaixuan.txt說明：cds.fai 是samtools faidx 對cds的序列文件fasta作用，生成的文件，result.txt 是拿cds進行多序列比對得到的結果，shaixuan.txt 是篩選后的結果存放的文本aln文件轉(zhuǎn)化 axt文件命令：/home/manag