生信分析27：基于Ks计算全基因组复制事件

全基因组复制事件（WGD）是某物种的全部染色体由于某种原因整体发生复制而加倍的事件。整体加倍可以是增加一倍，即二倍化事件，也可能增加了两倍，即三倍化事件。

大多数被子植物都经历过多倍化事件。物种在自然选择压力较大的时候可能会发生WGD来丰富自己的DNA原件，如基因、转座子，从而加速自身分化。WGD之后，小部分有利于物种进化的基因被保留，而大多数冗余的基因会消失或变成假基因。

鉴定WGD的方法

Fig 1

Fig 2

Fig 3

1、intragenome duplication：通过严格的比对标准和统计验证将基因组序列与自身进行比对。

（1）Ks 同义替换率（Fig 1）

（2）4DTv (4 fold Degenerate Transversion)，四重简并位点的颠换率（Fig 2）

2、double synteny：直接比较不同物种基因组间的线性关系（Fig 3）。目标物种基因组与参考物种比较，分析两者分化后目标物种特异的全基因组复制事件。需要保证选择的参考物种在与目标物种分化后没有发生自身的WGD，否则共线性关系会很混乱。

Ks计算全基因组复制的原理

Fig 4

理论上，两对由复制产生的基因应该有同样的“年龄”，也就是说由于多倍化产生的复制gene pair间序列的分歧应该是一致的。Ks 常用来作为分子钟计算来计算分化或者复制时间，因此Ks分布产生的峰常用来表示复制（物种内部分析）或者分化（种间分析）。

该分析的可视化通过Fig 4来展示。横坐标代表Ks值，纵坐标是基因对的百分比或者数量。相同物种间的比较检测复制，不同物种间的比较检测分化。

Ks对近期的全基因组复制更敏感，对古老的复制事件灵敏度较差。

Ks是同义替换，不受选择压力，呈现中性进化速率，更符合分子钟的假设。

物种间Ks计算

Fig 5

Fig 6

第一步 提取直系同源基因对（Fig 5）

可以利用wgd软件的dmd功能实现

Cisi.cds.fasta_Lich.cds.fasta.rbh存储了直系同源基因对的列表（Fig 6）。

Fig 7

Fig 8

Fig 9

第二步 计算Ks（Fig 7）

可以通过wgd软件的ksd功能实现，最终生成两个输出文件（Fig 8）

主要结果为tsv文件（Fig 9），第九列为Ks值。

Fig 10

svg文件对Ka、Ks以及两者比例做了最基本的展示（Fig 10）。

物种内部计算Ks值

Fig 11

Fig 12

第一步diamond 比对并进行MCL聚类（Fig 11），输出结果包含两个（Fig 12），Cisi.cds.fasta.mcl包含了聚类结果，每一行是一个基因家族。

Fig 13

第二步 计算Ks（Fig 13）

-mp 1000 用来去除非常大的基因家族

输出为Cisi.cds.fasta.ks.tsv，与前面物种间的分析结果一致。

注意特别大的Ks值要去除，在平时的分析中，一般保留小于3或5的Ks值。

Fig 14

需要注意的是对于物种内的分析，由于串联复制和散在复制的影响，Ks可能存在一个小峰（Fig 14），需要去除，即只保留位于共线性区块内部的基因对（未展示）。

Fig 15

画图（Fig 15）