直系同源基因ks_【比较基因组】如何利用paml计算kaks

前言：

本篇文章的作者有两人，一者是我，另一者是知乎ID：OrzKill，江湖人称表哥，是重庆师大的一名研究生（研究蜜蜂的，对，没错，蜂蜜什么的就可以找他要一些珍品~嘻嘻嘻）。这篇专栏文章的缘由是因为我在群里说自己是用paml算的kaks，表哥知道后就一直让我教他，而后在两个人互相交流下（其实我只是教会他用paml简单的算kaks而已，获得paml的输入文件100%是他的功劳），整理出一个新的流程出来，简单粗暴，而且高效，舍弃了我之前自己获得paml输入文件时的繁琐流程，故而将这个流程完整分享给大家。

首先，我们要明白kaks是个虾米，这里简单说一下，就是选择压力，即正选择与负选择，可参考下文

dNdS与KaKs的关系，你搞清了吗？​www.360doc.com

本次测试环境为linux系统，blast+版本2.6.0+，python版本3.7，在其他系统下运行本流程的小伙伴可能需要自行修改一些东西，不过大体是没变的。

###本次测试数据
链接：https://pan.baidu.com/s/1OuHduB8w_xuMWv9CAwqivQ 
提取码：ftd3

本次使用到的测试数据为酵母的两个近缘种，数据来源于NCBI，下载完基因组后用transdecoder预测得到cds跟pep文件，数据量不大，在虚拟机或者子系统下均可以跑通。

###首先用conda安装transdecoder
conda install -c bioconda transdecoder
###接着用使用transdecoder预测，测试数据为1.fna跟2.fna
TransDecoder.LongOrfs  -t 1.fna

在1.fna.transdecoder_dir里面获得预测的cds跟pep，正常还需要去冗余序列，但这里仅作为演示方便，所以直接接着往下做

我重新命名了一下，并且去掉序列名后面的orf注释，在测试文件分别对应1.cds跟1_new.cds

去掉前

去掉后

###使用python3去掉，用awk或者sed都可以，只不过我个人熟悉python
import sys
x = sys.argv[1]
file = open(x, "r")
lines = file.readlines()
for line in lines:
    if ">" in line:
       strings = line.split(" ")[0]
       print(strings)
    else:
       print(line.rstrip())
###保存为change_name.py
python change_name.py 1.cds > 1_new.cds