自然選擇是五種遺傳力(突變、重組、選擇、基因流、漂變)之一,選擇壓力分析更是進化分析中不可或缺的一項重要內容。
自然選擇的檢測
非同義替代與同義替換的比值,即:ω 值,也就是通常所說的dN/dS(或Ka/Ks)。
(1)當ω =1時,中性進化 (Neutral selection),即不受選擇:
(2)當ω >1時,正選擇(Positive selection);
(3)當 0 <ω <1時,負選擇(Negative selection,也叫凈化選擇或純化選擇 Purifying selection)
不導致氨基酸改變的核苷酸變異我們稱為同義突變,反之則稱為非同義突變。
一般認為,同義突變不受自然選擇,而非同義突變則受到自然選擇作用。
在進化分析中,了解同義突變和非同義突變發生的速率是很有意義的。常用的參數有以下幾種:同義突變頻率 (Ks)、非同義突變頻率(Ka)、非同義突變率與同義突變率的比值(Ka/Ks)。如果Ka/Ks>1,則認為有正選擇效應。如果Ka /Ks=1,則認為存在中性選擇。如果Ka/Ks<1,則認為有純化選擇作用。
dn/ds又叫ka/ks dN>dS 表示基因受到正選擇:
楊子恆老師在他所著的《Computational Molecular Evolution》一書第8章第1節中這樣寫道:
With the synonymous rate used as a benchmark, one can infer whether fixation of nonsynonymous mutations is aided or hindered by natural selection. The nonsynonymous/synonymous rate ratio, ω = dN/dS, measures selective pressure at the protein level. If selection has no effect on fitness, nonsynonymous mutations will be fixed at the same rate as synonymous mutations, so that dN = dS and ω = 1. If nonsynonymous mutations are deleterious, purifying selection will reduce their fixation rate, so that dN < dS and ω < 1. If nonsynonymous mutations are favoured by Darwinian selection, they will be fixed at a higher rate than synonymous mutations, resulting in dN > dS and ω > 1. A significantly higher nonsynonymous rate than the synonymous rate is thus evidence for adaptive protein evolution.
翻譯如下:
以同義突變率作為判定標准,可以推斷非同義突變的保留是受到自然選擇的支持還是阻礙。非同義突變率與同義突變率的比率,ω = dN/dS,可用以測量在蛋白質水平上的選擇壓力。
如果自然選擇對基因的適合度沒有影響,非同義突變將以和同義突變相同的速率被保留,即dN = dS,ω = 1。
如果非同義突變是有害的,凈化選擇將減小它們的保留速度,那么dN < dS,ω < 1。
如果非同義突變更受達爾文選擇偏愛,它們將以大於同義突變的速度被保留,導致dN > dS,ω > 1。
非同義突變速度顯著高於同義突變速度因此成為適應蛋白進化的證據。
CodeML中四個常見的模型假設及其工作流
在CodeML中,考慮不同序列間(考慮位點)或系統發育上的支系間(考慮支系Lineage)的ω 值不同,主要有以下四類常見的模型:
1. 枝模型(Branch model)主要用於對系統發育樹中不同支系 ω值差異性進行界定,主要有三個模型:
(1)One-ratio model:假設系統發育樹中所有支系的 ω 值相等;
(2)Free-ratio model:假設系統發育樹中所有支系的 ω 值不相等;
(3)Two-ratio model:假設前景枝和背景枝的ω 值不同;
2. 位點模型 (Site model)主要假設數據集中不同氨基酸位點受的選擇壓力不同(而不考慮不同支系間受的選擇壓力差異)。
該模型主要用於檢測正選擇( ω >1)作用,共有8個不同假設的模型:
(1)M0(單一比率),即:One-ratio model,假設所有位點具有相同的 ω 值;
(2)M1a(近中性),假設僅有保守位點(0<ω <1)和中性位點( ω =1)而沒有正選擇位點( ω >1)存在,這兩類位點的比率分別為p0和p1,其對應的ω 值分別為ω0、ω1;
(3)M2a(正選擇),該模型在M1基礎上增加了第三類ω值,即假設除了保守位點和中性位點外,還存在處於正選擇壓力下的位點( ω >1),這三類位點的比率分別為p0、p1和p2,其對應的ω 值分別為ω0、ω1和ω2;
(4)M3(離散),假設所有的位點ω 值呈簡單的離散分布趨勢;
(5)M7(beta),假設所有位點的 ω 屬於矩陣(0, 1)並呈beta分布;
(6)M8(beta & ω ) ,該模型在M7基礎上增加另一類ω 值(ω >1);
(7)M8a(beta & ω =1),與M8模型類似,但將ω 值固定為1(ω =0);
3. 枝位點模型 (Branch site model):主要假設不同氨基酸位點的和不同支系間受的選擇壓力均存在差異(既考慮位點間也考慮支系間的 ω 值存在差異),共有四個模型Model A、Model B、Model C和Model D,主要參數如下:
(1)Model A (Model 2, NSites=2, ncatG=ignored)
(2)Model B (Model 2, NSites=3, ncatG=ignored)
(3)Model C (Model 3, NSites=2, ncatG=ignored)
(4)Model D (Model 3, NSites=3, ncatG=2 or 3)
4. 進化枝模型 (Clade Model):與枝位點模型類型,能同時檢測多個進化枝(Clade),共有CmC和 CmD 兩種模型,主要參數如下:
(1)CmC (Model 3, NSites=2, ncatG=2 or 3)
(2)CmD (Model 3, NSites=3, ncatG=ignored)
CodeML 分析大致工作流簡要如下:

基因選擇壓力

將序列基於密碼子比對后,提交序列分析時,先選擇一個最適的進化模型,再選擇一個壓力檢測方法。
前者內置四種模型,包括Branch model、Branch-site model、Site model和Clade model,每個模型至少都有一組成對模型:零假設模型(Null model )和備選假設模型(Alternative model),其中Site model還內置三組成對模型(M0 vs. M3,M1a vs. M2a, M7 vs. M8)
楊子恆PAML選擇壓力分析之codeml(dN dS)
Codeml是PAML軟件包下的一個程序,廣泛應用於估算蛋白編碼序列的同義替換和非同義替換速率,及檢測序列是否受到正選擇。
我們通過實例來分析某個編碼序列氨基酸位點是否存在正選擇作用。
安裝好PAML, 准備序列比對文件,核酸序列得是3的倍數,序列中不能出現純數字和特殊符號,但_可以,不能有額外的空格換行符等,
可用bioedit等工具將編碼序列翻譯后再進行比對,然后保存成phylip格式、Paml可識別phylip格式
樹文件准備-任何可以正確描述序列文件中各序列文件的tree都可以,如果沒有可以借鑒使用的樹,可以用序列文件中的序列用Mrbayes, Phyml, PAUP等軟件建樹,先建立無根樹
PhyML產生的樹,格式可被Paml識別,Paup產生的樹Nexus格式可能不能直接被Paml識別,可用figtree軟件轉存為Newick格式。
運行codeml程序,將序列比對文件和樹文件同時拷貝到paml/bin文件夾下,或者確保codeml程序文件、codeml.ctl配置文件、比對好的XX.phy文件和序列文件、XX.tree樹文件這四個文件處於同一文件夾中。
codeml.ctl配置文件的參數設置

seqfile 序列比對文件 treefile 樹文件 outfile 輸出文件 (任意命名,建議txt后綴便於打開)seqtype=1 clock=0 (使用無根樹需設置為0) model=0
Nsites= 07 8 (0為單一參數模型,7 8 為相似模型,均允許局部替換率可以變化,兩者有一個參數差異,其二者結果將用於LRST檢驗以保證結果的可靠性) 其余參數不做改變。
Model 0的結果里查找omega (dN/dS)的值,這里是0.29437 。(或者在分支模型計算結果里是幾個omega的均值是最大且大於1的,且p值是顯著的)
kappa是ts/tv
考慮單個氨基酸的dN/dS值,要先比較 Model7和model 8的Lnl (lnL )值,以確定哪一個模型更適合序列,然后再選擇該模型的結果。
LRT = 2dl = abs(2 X (Lnl7-Lnl8)) (abs=絕對值)
第一步: 打開輸出文件,找到Model7 和Model8兩個Lnl 值相減,取絕對值,乘以2
第二步: 打開Paml附帶的chi2程序,分別輸入自由度和第一步的結果,自由度直接取1,因為M7,M8一個參數差異 輸入:chi2 1 35,
即出現: df=1 prob = 0.02534=2.535e-02 得到的p小於0.05 (這里在branch model里最終結果序列的模型的p小於0.05 )
此處Paml用Naive Empirical Bayes (NEB) 分析檢測顯著性
CODONML in paml version 4.9f,
做選擇壓力分析的序列文件,需要先clustal和剔除終止密碼子,這步可以在MEGA里完成。
1.MEGA7
MEGA7,是通過分別求取dn,ds的值,然后得到dn/ds。
將fasta格式的序列文件導入MAGA,
然后選擇Distance——computer overall mean distance進入圖一頁面:仔細看下面的選項卡,在substitutions type 中選擇 Syn-Nonsynonymous;Genetic code table 按照自己的序列選擇;modle/method 選擇 Nei-Gojobori method (No. of Differences);Substitutions to include 選擇要計算的dn或者ds。下一步,就能得到dn或者ds,兩者相比得到結果。
如果在distance——computer pairwise#####,然后按照后面步驟操作,結果會得到一個兩兩比較的矩陣(三角),我還不知道這個要怎么用。
如果只計算dn/ds,第一種應該夠用了。
作者覺得視頻教程里up主說做的是branch,然后取model=0是零假設,與我理解的有出入,我以為model=0是site,NSsite可以決定哪個是零假設哪個是備擇假設。
