更新20181216
今天又把KaKs_calculator2.0的幫助手冊翻出來看了一遍,找到了一句話:Before calculation, gaps and stop codons between compared sequences will be removed.軟件應(yīng)該是自己會自動刪除gaps和終止密碼子,比對后應(yīng)該不用自己手動刪除終止密碼子
參考 組學(xué)大講堂 公眾號分享的文章 如何使用 KaKs_Calculator 計算Ka/Ks 重復(fù)論文 Chloroplast genome analysis of resurrection tertiary relict Haberlea rhodopensis highlights genes important for desiccation stress response? 中ccsA基因的ka/ks的計算過程:
論文中用到了五個物種ccsA基因,根據(jù)拉丁名搜索NCBI organelle genome數(shù)據(jù)庫找到對應(yīng)的CDS序列,使用 mafft 比對
mafft --auto --clustalout ccsA_KaKs_pra.fasta > ccsA_aligned.aln
(小疑問:mafft比對結(jié)果是默認把氨基酸大寫字母轉(zhuǎn)換成小寫字母嗎?)
KaKs_Calculator2.0 download | SourceForge.net 由此網(wǎng)址下載KaKs_calculator
解壓出來以后分別用到src目錄下的 AXTConvertor 程序和bin/Linux目錄下的 KaKs_Calculator 程序,使用之前使用 chmod u+x?AXTConvertor
chmod u+x?KaKs_Calculator 賦予文件可執(zhí)行權(quán)限
之后用?AXTConvertor 程序?qū)⒈葘Y(jié)果文件轉(zhuǎn)化為axt文件
./src/AXTConvertor ccsA_aligned.aln ccsA.axt
然后用?KaKs_Calculator 計算 Ka/Ks值
./bin/Linux/KaKs_Calculator -i ccsA.axt -o ccsA_2.kaks -m YN -c 11
-o參數(shù)指定輸出文件, -m參數(shù)指定方法,-c參數(shù)指定密碼子表;具體內(nèi)容可以使用 -h 參數(shù)查看幫助文檔
ccsA_fasta_and_kaks_result.rar_免費高速下載|百度網(wǎng)盤-分享無限制
分析所用到的fasta文件和最終輸出結(jié)果的百度云鏈接
問題:輸入5條序列會輸出C52個Ka/Ks值,最終的結(jié)果是取平均值嗎?
最后用ggplot2簡單畫一幅箱線圖

pal2nal.pl腳本將蛋白比對轉(zhuǎn)換成CDS比對命令
perl pal2nal.pl example_pep_aln.fasta example_CDS.fasta -output.fasa > example_CDS_aln.fasta
下載地址?http://www.bork.embl.de/pal2nal/#Download
更新:20181215
比對格式之間互相轉(zhuǎn)化(使用biopython的AlignIO模塊)
from Bio import AlignIO
from Bio import Alphabet
AlignIO.convert("inputfile","fasta","outfile","clustal",alphabet=Alphabet.generic_dna)
參考文獻
微信文章: 《 一個計算同義突變(Ks)與非同義突變(Ka)的方法》(其中有1:先比對蛋白序列,然后利用perl腳本轉(zhuǎn)換成CDS比對;2:將核酸序列聯(lián)配轉(zhuǎn)換成AXT格式的python腳本)