对遗传进化过程进行模拟是研究包括自然选择在内的进化过程的重要工具和方法。常用的模拟工具,有溯祖模拟ms, msms, discoal等,向前模拟工具SLiM等。本文主要介绍discoal的基本使用。
discoal
可以进行溯祖模型的模拟,它是discrete
和coalecent
两个词的简写。在模拟的过程中,可以添加重组、群体大小变化和自然选择等因素。它和ms
等模拟工具比较相近。
下载安装
# https://github.com/kr-colab/discoal
cd <PATH TO DISCOAL>
make discoal
可以在Linux和OSX系统上安装使用。
基本使用
./discoal sampleSize numReplicates nSites -t theta
sampleSize: 抽样或者样本的数量大小,注意不是有效群体数量 numReplicates: 进行模拟的次数 nSites: 模拟序列的长度(bp)
-t
用来指定群体多态性水平 的大小,其中是当下群体数量,表示突变率(突变事件每世代每位点),L为序列长度。
例如:3个样本,进行2次模拟,序列长度为100
./discoal 3 2 100 -t 2
输出结果如下
./discoal 3 2 100 -t 2
1665047201 686400060 #随机数种子
//
segsites: 4
positions: 0.01835 0.09557 0.46556 0.72880
1000
0110
0111
//
segsites: 1
positions: 0.07594
0
1
0
每次模拟的输出用//
分割开,第一次模拟产生了4个多态位点,第二次模拟有一个多态位点。
重组和基因转换
通过-r
可以用来指定模拟过程中的重组率: ,其中表示序列交叉互换的概率每碱基位点。
通过-g
可以用来指定基因转换的参数,(一般不常用)。
群体大小变化
通过-en
可以设定某个群体在某个时间点的群体大小。可以重复多次使用
-en [时间点] [群体ID] [群体大小]
其中时间是以 世代数 为单位,群体大小初始群体数量的比值。如下:
./discoal 3 2 100 -t 2 -r 2.4 -en 0.5 0 0.1 -en 1.2 0 0.8
该群体(ID=0) 在0.5时间单位的时候,群体数量变为了初始(当下)群体数量的0.1倍;其后,在1.2时间单位的时候,群体大小恢复到了初始时的80%。
在有多个群体的时候,可以使用-p
来指定群体数量和每个群体的样本数。
-p [群体数量] [样本数] [样本数] ...
如果群体之间有基因流存在,可以使用-M
来指定他们相互之间的迁徙率,以为单位。-m
可以单独指定某两个群体之间的迁徙率。
./discoal 6 2 100 -t 2 -r 2.4 -p 3 2 2 2 -M 0.05
上例中,有6个样本,分别来自3个群体,每个群体2个样本,这3个群体相互之间的迁徙率为0.05.
如果有群体分化,可以用-ed
来指明,比如:
./discoal 6 2 100 -t 2 -r 2.4 -p 3 2 2 2 -ed 1.0 0 1 -ed 5.0 1 2
第一个群体(0)和第二个群体(1)在1.0个时间单位分化,第二个群体(1)和第三个群体(2)在5.0个时间单位分化。
如果群之间有有混合,可以用-ea
,如下:
-ea [混合时间] [混合群体ID] [两个混合前群体ID] [混合群体中来自第一个群体的比例]
./discoal 6 2 100 -t 2 -r 2.4 -p 3 2 2 2 -ea 0.002 0 1 0 0.15
在0.002时间单位的时候,第二个群体(ID=1)和第一个群体(ID=0)进行了混合,其中混合群体中有15%的个体来自第二个群体(ID=1)。
自然选择
通过-ws
可以来说明一个自然选择清除的发生,
-ws [时间]
比如-ws 0
说明选择位点在采样前即刻在群体中固定,频率达到1;如果时间大于0, 说明采样前t个时间单位()已经固定。
-x
可以用来指定受选择位点在序列上的位置,【0-1】之间,默认为0.5,即在序列的中间位置。
-a
可以用来指定受选择位点的选择压力大小,
比如:
./discoal 6 2 100 -t 2 -r 2.4 -ws 0.05 -a 1000 -x 0.5
表示一个100bp的序列,在第50bp的地方有一个受选择位点,选择压力, 该位点在采样前0.05个时间单位的时候在群体中固定。
如果选择性清除是一个软性选择性清除,可以通过-f
来自指定受选择时,起始等位基因频率。比如-f 0.1
,说明选择发生时,受选择位点的频率为0.1.
对于未固定的选择性清除partial sweeps
,可以使用-c
来指定选择结束时受选择位点的频率。如下一个为固定的软性选择性清除
./discoal 3 2 100 -t 2 -r 2.4 -ws 0.05 -a 1000 -x 0.5 -f 0.1 -c 0.8
一个起始频率为0.1的位点,在采样前0.05个时间单位的时候,其频率达到了0.8。
另外,-T
参数可以将模拟结果以newick进化树的格式输出。
参数区间
上述的例子中,参数都是以单个点估计值给的。如果遇到比较严格的审稿人,可能会怀疑,你是否用其他不同的参数进行了模拟测试。
在discoal
中,参数可以是区间的形式。在相应的参数前加-P
即可以给一个参数区间。具体参数如下:
-Pt low high uniform distribution of θ drawn from the range [low, high]
-Pr low high uniform distribution of ρ (recombination rate for the entire locus)
-Pre mean upperBound exponential distribution of ρ, truncated at upperbound
-Pa low high uniform distribution of α (selection coefficient)
-Pu low high uniform distribution of τ (time since fixation)
-PuA low high uniform distribution of uA- recurrent adaptive mutation rate
-Px low high uniform distribution of sweep position
-Pf low high uniform distribution of f0 (initial selected frequency)
每次模拟,模型将在参数期间内随机取值,最后的模拟结果将代表了整个参数区间。
比如下面例子,给出了突变率的取值期间,在【2-8】内均匀分布,每次模拟随机取值。
./discoal 3 2 100 -Pt 2 8
资料来源:https://github.com/kr-colab/discoal