基因选择压力分析 | Ka/Ks, Dn/Ds 大规模计算,更快更准!
写在前面
来来去去,这个主题相关的,我已经写过很多个推文,甚至还包括与Raindy(福建农林高芳銮老师)合作开发和发表了 EasyCodeML。只能说,确实是一个相对常见但似乎也繁杂的分析。
功能开发
EasyCodeML 主要用于在进化的context上分析选择压力,尤其是做正选择位点分析,亦即,关注到位点。
简单的基因对Ka/Ks计算,如 TBtools 的 Simple Ka Ks Calculator 则关注大规模的基因受选择的情况,亦即,关注到基因。
正选择位点分析起来不容易,而大规模的基因对Ka/Ks计算起来其实也不容易。常见的场景是有数万,甚至十来万个基因对,需要计算KaKs。而做计算的步骤:
以密码子为单位进行序列比对
基于比对结果进行KaKs计算
所以是两部。绝大多数人会使用一些软件,如muscle,mafft等做基因双序列比对,因为他们调用起来放弃,尽管这其实不太合适。因为这些软件本身设计目标是多序列比对。用于CodeML的计算应该使用他们,但是用于两条序列的两两比对,却并不合适。所以最好的方式基本只有:
调用Emboss的Needle程序
使用已有第三方实现
自行做代码实现
其中第一种方法最准但是最慢,第二种方法在Java,python,perl等上面没有良好实现,一般不够准确,甚至不如调用muscle等,第三种方法可以权衡,毕竟是自己做代码实现。
而现在,TBtools用的即第三种方法,准确度上跟Emboss Needle完全一致,但由于是原生代码实现,所以速度很快,同时也支持多线程。
运行速度
针对调用muscle的实现,前述已经提过,大体是一秒钟一对序列,多线程无法提速,甚至会降速(怀疑是进程开销太大),而针对TBtools的原生代码实现(也就是说,我自己coding的),那么单线程一秒钟可以做到200对,如果开四个线程,那么一秒钟可以做到 600对(多线程开销其实并不小)。
实测数据,10000个基因对,muscle调用需要四个多小时。使用TBtools,单线程只需要不到一分钟。如果开四个线程,28秒。
换句话说,常规电脑上去,1w个基因对的计算,一般不需要1分钟。
使用方法
打开 TBtools,跳转到对应功能
使用方式和以前一样,就是记得调整使用的线程数。。。
写在最后
天下武功,唯快不破。