实验篇——Ka/Ks分析
文章目录
- 前言
- 一、名词解释
- 二、实操
- 1. 安装软件
- 2. 准备文件
- 3. 使用ParaAT2.0 + KaKs_Calculator2.0
- 4. 使用TBtools软件
- 三、额外
- 总结
前言
鉴定不同基因的复制模式
本文得到的共线性基因对文件 来自于上一篇文章中的LIN.collinearity共线性文件
参考文章：
https://www.jianshu.com/p/a3a39b2f341b
一、名词解释
Ka/Ks分析是一种用于评估基因进化速度的方法。较高的Ka/Ks比值可能表明基因在进化过程中经历了功能上的重要变化，而较低的Ka/Ks比值可能表明基因在进化过程中具有保守的功能。
更多相关的解释信息
二、实操
1. 安装软件
安装ParaAT、KaKs_Calculator2.0以及muscle（或者其它序列比对软件） 软件，并配置环境
或者在TBtools软件中操作
2. 准备文件
准备存储共线性基因对的文件 homolog、cds文件、pep文件（注意，cds及pep为fasta格式，“>”后面只接基因ID号）
以LIN.collinearity共线性文件为例
得到homolog文件
cat LIN.collinearity | grep "Lagg" | awk '{print $3"\t"$4}' >LIN.homolog
得到cds文件
原来的文件格式：
sed 's/^>\([^\t]*\)\t.*/>\1/' Lagg.gene.cds >Lagg.cds
要求的文件格式：
得到pep文件
原来的文件格式
sed 's/^\(>[^ ]*\).*/\1/' Lindera_aggregata.gene.pep >Lagg.pep
sed -i 's/\*//g' Lagg.pep
要求的文件格式：
或者也可以用别的编程语言实现，如python
import re
with open("D:\yuceji\shili.pep", 'r') as file:
lines = file.readlines()
output_file = 'output.fasta'
with open(output_file, 'w') as file:
for i in range(len(lines)):
line = lines[i].strip()
if re.match(r'^>', line):
gene_id = line.split()[0]
seq = ''
j = i + 1
while j < len(lines) and not re.match(r'^>', lines[j].strip()):
seq += lines[j].strip()
j += 1
file.write(gene_id + '\n')
file.write(seq + '\n')
3. 使用ParaAT2.0 + KaKs_Calculator2.0
ParaAT.pl -h test.homologs -n test.cds -a test.pep -p proc -m muscle -f axt -g -k -o output 2> ParaAT.log &
-h test.homologs：指定包含同源序列的文件
-n test.cds：指定包含核苷酸序列的文件
-a test.pep：指定包含蛋白质序列的文件
-p proc：指定并行处理的进程数, (proc 可以自行设置）
-m muscle：指定用于多序列比对的软件，这里使用的是muscle
-f axt：指定输出文件的格式，这里使用的是axt格式
-g：启用gap stripping功能，即移除比对序列中的缺失片段。
-k：启用Ka/Ks比值的计算
-o output：指定输出文件的名称
2> ParaAT.log：将标准错误输出重定向到名为"ParaAT.log"的日志文件中。
&：将命令放入后台运行。
4. 使用TBtools软件
基于TBtools软件，可以使用" Sinple Ka/Ks Calculator" 程序
同样依然要准备那三个文件
从左到右每一列的介绍：
-
Seq_1：表示序列对中第一个序列的标识符或名称。
-
Seq_2：表示序列对中第二个序列的标识符或名称。
-
Ka：表示非同义突变（氨基酸替换）的数量，也称为Ka值。
-
Ks：表示同义突变（氨基酸保守替换）的数量，也称为Ks值。
-
Ka/Ks：表示Ka值除以Ks值得到的比值。（Ka/Ks比值用于衡量非同义突变和同义突变的相对丰度，从而推断基因或序列的选择压力）
-
EffectiveLen：表示序列对的有效长度，即用于计算Ka和Ks的比对序列的长度。
-
AverageS-sites：表示平均同义突变位点的数量。
-
AverageN-sites：这一列表示平均非同义突变位点的数量。
-
cN：这一列表示非同义突变位点的校正计数。
-
cS：这一列表示同义突变位点的校正计数。
-
pN：表示非同义突变位点的概率。
-
pS：这一列表示同义突变位点的概率。
其实得到的文件中除了这12列外，还有一列”Note",来记录一些“high sequence divergence value (ps>=0.75)” 这类的信息。
三、额外
总之，最主要的是得到了Ka/Ks的计算结果，达到了我们的目标，其余的也没必要全都了解。
另外，对于结果文件的每一列，我是这样记忆的：
对于Ka与Ks 的区分，可以Ks中的“s"来记忆，因为这可以很快与"same"这个单词联想起来。即相同的意思，故Ks就是同义突变。至于另一个Ka就是非同义突变了。
AverageS-sites 与AverageN-sites 的记忆也是如此。
cN 与 cS ，也是如此，另外"c"可表示"count",即计数
pN与pS，也是如此，另外"p" 可表示"probablity",即概率的意思
总结
本文主要讲述的是Ka/Ks 值的计算，其实无论是使用ParaAT2.0 + KaKs_Calculator2.0来计算，还是使用TBtools软件中的" Sinple Ka/Ks Calculator" 程序。都是为了得到结果。我个人还是比较推荐使用TBtools软件的，因为足够简单（不需要在Linux中再下载一些什么软件后再配置环境了，只要输入三个文件，就能得出结果文件）。当然或许它也会有一些其它问题。这样我们可以具体情况再具体分析。
无丝竹之乱耳，无案牍之劳形。
–2023-8-21 实验篇
