【计算分子进化】基于伽马分布建模位点可变置换率

在讨论核苷酸替换模型时，并不考虑一条序列上的不同位点可能存在不同的置换率的情况，即假定不同位点都以相同的方式和速率进化。

这种假设显然是不符合真实情况的，不同位置的序列受到不同的选择压力、突变偏差的影响。举几个可能的情况：

1.编码区和非编码区，受到不同大小/来源的选择压力影响

2.编码基因的外显子和内含子区，收到不同大小/来源的选择压力影响

3.基因内，基因起始和终止区会受到更强的纯化选择压力

4.密码子内，第一位和第二位也会受到更强的选择压力

5.某些区域的GC含量会显著高于其他区域

因此，需要考虑位点间可变的置换率。高置换率的区域，可能积累了更多的变化，如果不考虑这种情况，就会导致我们低估序列距离。 

---

1.伽马分布

要了解gamma分布，首先要理解泊松分布和指数分布

1.泊松分布，模拟一个离散事件在连续时间或空间内发生次数的概率分布，典型例子有：某网站在某段时间内的点击量，学校在某段时间内进入的学生数等等

2.指数分布，模拟连续时间/空间上离散事件发生的间隔时间的概率分布，常用于可靠性分析和等待时间模型，典型例子有：某网站被访问的时间间隔。

3.伽马分布：模拟未来事件发生前的等待时间，即直到发生第k个事件的等待时间，而事件发生的概率还是遵循泊松分布的。

设位点速率r服从伽马分布，伽马分布的密度函数为：

 

 （这里摘自百科）

 

在用伽马分布建模位点速率时，设beta = alpha，来减少参数的使用，即我们这里可以只考虑形状参数alpha。

同时，在这种假设下均值为1，方差为1/alpha。

 当alpha取不同的值时，实际上描述的是不同类型的序列 

1. alpha＞1,则分布为钟形，意味着大部分位点为 1附近的中等速率，而极少位点以很慢或很快的速率变化。当 alpha一>∞时，分布衰减为所有位点只有一个单一速率的模型。

2. alpha＜1,分布为一个高度倾斜的L形，意味着除一些快速置换“热点” 外，大部分位点的置换率都非常低，或者几乎“不变”。

在实际使用时，需要联合多条真实序列数据来估计合适的alpha

2. 基于K80的伽马距离

ps:原K80模型的P矩阵

 这里给出推导的思路，具体过程和式子详见课本

 课本上还同理推导了JC69可变位点速率模型，其思路都是一致的。

同理可以推广到F84,TN93,HKY85模型等。

可能后续还会更新本节，比如加点推导细节。。。。