背 景
基因家族的分類是以其編碼的蛋白質結構、功能及進化關係為基礎的。基因家族成員既可以是不同物種中具有進化相關性的基因（直系同源基因，orthologs），也可以是同一物種中具有相似生物學功能的同源基因（旁系同源基因，paralogs）。家族成員鑑定主要是利用同源基因（蛋白質）之間的序列相似性。那麼研究基因的進化是選擇DNA序列還是蛋白序列呢？一般情況如果序列的相似性≥70%，選擇DNA序列；當<70%時DNA或蛋白序列皆可。BLAST工具可以用來快速搜索相關序列進而進行同源序列的鑑別，然而BLAST是一種局部比對方法，從整體上講它找到的眾多相關序列之間可能關係很遠；而且BLAST並未在其分析中增加查詢序列中保守區域的權重，從局部來說其發現差異較大的同源基因（如特定物種的種內旁系同源基因）的能力也受到了限制。另一方面，Pfam、SMART、InterPro等資料庫則利用序列譜、隱馬爾可夫模型等方式記錄同源基因保守區域的序列特徵，它們都是從同源序列的多序列比對中總結出來的統計學模型，在遠源家族成員的鑑定中表現更為出色。因此利用序列特徵的統計學模型進行家族成員的鑑定是對BLAST方法的一個很好的補充。然而，序列中的保守區域一般不會覆蓋到序列整體，因此對於一些具有複雜結構的基因當使用只針對局部保守序列的統計學模型進行預測時往往產生的結果數量過多，序列之間的差異加大。如果要確認一組同源基因的功能相關性，還需要藉助系統發育分析的方法來幫助判斷。構建系統發育樹有三種主要的建樹方法，分別是距離法、最大節約法(maximumparsimony,MP)和最大似然法(maximum likelihood,ML)。最大似然法考察數據組中序列的多重比對結果，優化出擁有一定拓撲結構和樹枝長度的進化樹，這個進化樹能夠以最大的概率導致考察的多重比對結果；距離樹考察數據組中所有序列的兩兩比對結果，通過序列兩兩之間的差異決定進化樹的拓撲結構和樹枝長度，基於距離的方法有UPGMA、ME（Minimum Evolution，最小進化法）和NJ（Neighbor-Joining，鄰接法）等；最大節約法考察數據組中序列的多重比對結果，優化出的進化樹能夠利用最少的離散步驟去解釋多重比對中的鹼基差異。對於近緣序列，一般用MP，MP用到的假設最少；遠緣序列一般用NJ或者ML。本文將詳細闡述不同物種間的蛋白家族分析以及進化樹的構建
目 的
用擬南芥MAPKK家族的蛋白質序列作為已知參照，找出玉米基因組中可能的 MAPKK家族成員，並進行亞家族的分類，確定序列之間的遠近關係
擬南芥MAPKK蛋白家族介紹：
Arabidopsis MAP Kinase Kinase family
https://molbio.mgh.harvard.edu/sheenweb/mapk_project.html思 路
獲取參照基因家族（擬南芥MAPKK家族）的蛋白質序列，分析其序列特徵（Pfam Search）
下載玉米基因組蛋白質序列，提取其中包含參照家族序列特徵的序列（HMM search），對結果進行分析
在已獲得的序列集合中利用其他序列特徵進行篩選（psiblast phi_pattern），確定玉米中的 MAPKK
將玉米、擬南芥的MAPKK序列一起進行系統發育分析，確定序列之間的相互關係及亞家族分類
步 驟
1. 獲取擬南芥MAPKK家族的蛋白質序列，進行序列特徵分析
1) 進入TAIR網址（www.arabidopsis.org），從Browse-Gene Families中找到MAPKK基因家族，獲取 Genomic Locus ID
2) 進入TAIR網站，選擇Tools-Bulk data retrieval-sequences，利用Genomic Locus ID將相應序列以fasta格式進行下載，注意調整下載序列的類型為蛋白質序列(Araport11 protein sequences)，選擇search against 為 Get one sequence per locus
將下載的fasta文件存為Ara_MAPKK.fasta，各個基因行標可改為以下形式：>ATMKK1_AT4G26070.2
Ara_MAPKK.fasta
3) 進入Pfam網站(https://pfam.xfam.org/)，選擇SEARCH - Batch search對Ara_MAPKK.fasta進行Pfam庫中的序列特徵分析。在Sequences file欄中上載入序列文件，Cut-off選擇Gathering threshold，填入Email address，點擊Submit進行分析
4) 將Email中的結果保存為Ara_MAPKK_Pfam.txt，在pfam網站通過Pfam ID查看此基因家族的序列特徵並在Curation & model處下載這一序列特徵的HMM ，存為 Pkinase.hmm
Ara_MAPKK_Pfam.txt
Pkinase.hmm
2. 玉米基因組蛋白質序列的獲取，特徵分析以及符合特徵條件的序列提取
1) 從maizegdb網站下載玉米的蛋白序列
2) 利用 Hmmer 工具中的 hmmsearch 命令對玉米基因組中的蛋白質序列進行Pkinase序列特徵(Pkinase.hmm)的鑑別，得到文件 Pkinase_Maize_HMM_domtblout
Hmmer連結：
進行Pkinase序列特徵鑑別
hmmsearch --domtblout Pkinase_Maize_HMM_domtblout \--cut_ga ../Pkinase.hmm ./Zea_mays.B73_RefGen_v4.pep.all.faPkinase_Maize_HMM_domtblout
3) 從hmmsearch的結果中獲取非冗餘的基因ID（編碼蛋白可能會有多個區域符合Pkinase_hmm序列特徵），並根據這個ID_list從玉米全基因組蛋白質序列數據中提取相應的蛋白質序列，並以fasta格式保存提取的序列，命名文件為Pkinase_Maize_HMM_domtblout_out.fa
3. 用blast工具中的psiblast對上一步驟獲得的結果進行進一步過濾及提取
1) Pkinase.hmm這個序列特徵能夠把很多蛋白激酶超家族的成員找到，而MAPKK只是超家族中的一個很小的有其特定功能任務的一個小組。因此，我們需要繼續尋找MAPKK的其他序列特徵，並用這個序列特徵對獲得的結果進行進一步的限定。在TAIR - gene family對MAPKK家族的描述中，包含了另一個很小的序列特徵片段。我們把這個特徵片段保存到一個名為Ara_MAPKK_pattern.txt的文件中
Ara_MAPKK_pattern.txt
2) Blast系列比對工具中有一個psiblast，psiblast 有一個參數（-phi_pattern）就是在蛋白質序列比對的結果中利用特定的pattern特徵（Ara_MAPKK_pattern.txt）對獲得的結果進行匹配限定。分析的目標序列數據則要先利用makeblastdb命令進行格式化以便能夠讓blast 系列工具識別
blast工具連結：
https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST利用makeblastdb進行格式化
makeblastdb -in Pkinase_Maize_HMM_domtblout_out.fa –dbtype prot利用psiblast進行序列的進一步匹配限定，得到結果文件Pkinase_Maize_HMM_domtblout_VS_psiblast
psiblast -db Pkinase_Maize_HMM_domtblout_out.fa \-query Ara_MAPKK.fasta \-out Pkinase_Maize_HMM_domtblout_VS_psiblast -evalue 0.001 \-phi_pattern Ara_MAPKK_pattern.txt -outfmt 7Pkinase_Maize_HMM_domtblout_VS_psiblast
3) 獲取psiblast結果文件中的目標基因ID（第二列），並根據這個ID_list從玉米基因組蛋白質序列文件中提取序列並存為fasta格式的序列文件maize_genelocus.fa
4. 家族成員之間的遠近關係分析
1) 利用MEGA-X軟體構建系統發育樹
MEGA是一個功能非常強大的分子進化遺傳分析軟體，可用於序列比對、進化樹的推斷、估計分子進化速度、驗證進化假說等。做系統發育樹先要做多序列比對，然後把多序列比對的結果提交給建樹軟體進行建樹，所以在用MEGA建樹時可以輸入一個已經比對好的多序列比對，也可以輸入一條原始序列，讓MEGA先來做多序列比對，再建樹。MEGA的使用比較簡單，但是參數選擇至關重要，方法如下：
File-Open A File-選擇擬南芥的原始蛋白Ara_MAPKK.fasta為例（如果是要進行物種間的比較需將兩者的蛋白序列合併為一個.fa文件），隨後系統會詢問你選擇Align還是analyze，若是原始序列選擇Align，若是進行過多序列比對的文件選擇analyze
因為我的序列是原始序列，首先選擇Align，接下來需進行多序列比對，選擇Alignment，MEGA提供了ClustalW和Muscle兩種多序列比對方法，ClustalW的基本原理是首先做序列的兩兩比對，根據該兩兩比對計算兩兩距離矩陣，是一種經典的比對方法，使用範圍也比較廣泛。Muscle的功能則僅限於多序列比對，它的最大優勢是速度比ClustalW的速度快幾個數量級，我在這裡選擇Muscle，隨後會出現序列比對參數設置窗口（不知道怎麼改選擇默認參數即可），這就是序列比對的打分矩陣，之前已在序列比對一文說明，原文連結如下：
將多序列比對得到的結果保存為.meg文件用於後續建樹，選擇 DATA -Export Aligment - MEGA Format
從箭頭處打開.meg文件，得到一個 TA 文件
在建樹前要選擇一個最優的模型，提高建樹的精確度。如果想要快速建樹可以省去這一步，直接選擇默認的模型。選擇MODELS中的Find Best DNA/Protein Models(ML) ，軟體就會根據你的數據幫你計算尋找最適合的模型（選擇默認參數計算即可）。得到的結果中具有最低BIC分數（BayesianInformation Criterion）的模型被認為是最好地描述替代模式如JTT+G，但是MEGA提供的模型選擇有限，因此可在結果中找到BIC最低值對應的MEGA有的模型即可（JTT）
2) 利用ggtree可視化樹結構
其實MEGA軟體也能夠美化構建好的進化樹，但是顏值還是不夠，R包ggtree是個更好的選擇。我將玉米maize_genelocus.fa與擬南芥Ara_MAPKK.fasta兩個蛋白序列合併進行如上方法是進化樹構建，輸出newick格式的樹，需要保存其bootstrap值以及branch length值
ggtree 畫圖
https://yulab-smu.github.io/treedata-book/chapter1.htmllibrary(ggplot2)library(ggtree)setwd("C:/Users/fudiy/Desktop")tree <- read.tree("maize_Ara.nwk")pdf(file="file.pdf",width=10,height=6)p <- ggtree(tree) + geom_nodepoint(color="#b5e521", alpha=1/4, size=10)p + geom_tippoint(color="#FDAC4F", shape=8, size=3) + geom_tiplab(size=3, color="purple")dev.off()ggtree的方法很多，可以慢慢學著畫