MC1R基因編碼區(qū)的系統(tǒng)發(fā)育分析和正選擇位點(diǎn)檢測(cè)

李麗莎+李祥龍+周榮艷+李蘭會(huì)+張永改+馬夢(mèng)云+張秀玲

摘要： 黑素皮質(zhì)素受體1（melanocortin 1 receptor，MC1R）基因在黑色素細(xì)胞的表達(dá)上發(fā)揮重要作用，并與皮膚色素沉積密切相關(guān)。MC1R基因的選擇壓力分析對(duì)識(shí)別重要的功能位點(diǎn)具有深遠(yuǎn)影響。為研究不同動(dòng)物MC1R基因在進(jìn)化過(guò)程中是否受到正選擇作用，基于位點(diǎn)模型對(duì)MC1R基因進(jìn)行研究。結(jié)果表明，動(dòng)物MC1R基因主要經(jīng)歷中性漂變和純化選擇作用，并且發(fā)現(xiàn)MC1R基因受到正選擇作用的是編碼29、191、253、304氨基酸的位點(diǎn)。本研究從分子進(jìn)化的角度為MC1R基因調(diào)控黑色素的形成提供了一個(gè)新的思路。

關(guān)鍵詞： MC1R；分子進(jìn)化；親緣關(guān)系；系統(tǒng)發(fā)育；正選擇

中圖分類號(hào)： Q756 文獻(xiàn)標(biāo)志碼： A 文章編號(hào)：1002-1302（2016）03-0044-04

黑素皮質(zhì)素受體1（melanocortin 1 receptor，MC1R）基因，也稱為促黑素細(xì)胞激素受體（MSHR）基因，該基因只有1個(gè)編碼區(qū)，由擴(kuò)展位點(diǎn)（extension locus，E）編碼，是控制動(dòng)物黑色素合成的一個(gè)重要基因。哺乳動(dòng)物的MC1R一般由310個(gè)氨基酸組成，有7個(gè)跨膜結(jié)構(gòu)域，是G蛋白耦合受體-黑素皮質(zhì)素受體（melanocortin receptor，MCRs）的家族成員之一[1]，該基因在黑色素細(xì)胞的表達(dá)上發(fā)揮重要作用，并與皮膚色素沉積密切相關(guān)。MC1R與α-促黑素（α-melanocyte stimulating hormone，α-MSH）和腎上腺皮質(zhì)激素（ACTH）結(jié)合，增強(qiáng)了腺苷酸環(huán)化酶的活性，使環(huán)腺苷酸水平提高，進(jìn)而提高酪氨酸激酶（TYR）表達(dá)，從而形成黑色素[2]。它們之間的結(jié)合使動(dòng)物毛色或羽色的調(diào)控得以正常進(jìn)行。目前，已知一些MC1R基因變異引起腺苷酸環(huán)化酶活化能力不同，并且與小鼠、牛、馬、狐貍、家雞等動(dòng)物中的皮毛差異有關(guān)[3-4]。許多哺乳動(dòng)物有著同源特征的MC1R基因已被描述，而對(duì)于鳥類來(lái)講，MC1R基因色素性狀的調(diào)控機(jī)理與哺乳動(dòng)物相同[5]，因此MC1R基因通過(guò)控制真黑色素的數(shù)量來(lái)控制哺乳動(dòng)物毛色及鳥類羽色，某些動(dòng)物疾病的發(fā)生也與MC1R基因有關(guān)[6]。

正選擇（positive selection）也稱為達(dá)爾文選擇（Darwinian selection），是指將具有提高個(gè)體適合度的有利突變的等位基因在某個(gè)群體中固定下來(lái)的選擇作用，從而適應(yīng)不同的生存環(huán)境[7]。對(duì)正選擇的研究能更好地理解生物物種的進(jìn)化過(guò)程。從基因水平上來(lái)檢測(cè)適應(yīng)性進(jìn)化，這既有利于揭示生物的進(jìn)化史，同時(shí)加深對(duì)該基因的結(jié)構(gòu)和功能變異的理解[8]。

隨著對(duì)MC1R基因的進(jìn)一步深入研究，全面檢測(cè)MCIR基因位點(diǎn)將對(duì)動(dòng)物品種多樣性的保護(hù)和利用提供依據(jù)。對(duì)MC1R基因的研究已經(jīng)數(shù)年，傳統(tǒng)研究大多集中在MC1R基因的結(jié)構(gòu)、生物學(xué)功能和作用機(jī)理，往往忽略對(duì)其結(jié)構(gòu)、功能和作用機(jī)理起決定作用的基因序列的研究。在這些序列中蘊(yùn)藏著大量的生物信息，對(duì)今后MC1R基因的深入研究提供重要依據(jù)和參考。本研究對(duì)21種動(dòng)物MC1R基因編碼區(qū)進(jìn)行分子進(jìn)化分析，以期明確該基因在動(dòng)物演化過(guò)程中是否發(fā)生了正選擇作用，從而更加深入了解MC1R基因。

1 材料與方法

1.1 序列數(shù)據(jù)

本研究利用已有的MC1R基因序列數(shù)據(jù)進(jìn)行研究，所有的數(shù)據(jù)來(lái)源于NCBI網(wǎng)站（http：//www.ncbi.nlm.nih.gov/）GenBank數(shù)據(jù)庫(kù)，下載20種哺乳動(dòng)物45條MC1R基因的編碼區(qū)（963 bp），再加上1條原雞（AY220303.1 Gallus gallus）序列作為外群，總共21個(gè)物種46條MC1R基因序列（表1）。

1.2 方法

1.2.1 序列分析 利用Clustal W 2.1軟件對(duì)MC1R基因的蛋白序列進(jìn)行聯(lián)配，再通過(guò)Pal2nal[9]在線工具生成密碼子多序列比對(duì)，生成序列文件。使用Bioedit 7.0軟件[10]對(duì)序列進(jìn)行編輯并分別保存成FASTA格式和PHY格式的文件。用MEGA 6.0軟件[11]將比對(duì)序列轉(zhuǎn)化成NEXUS格式用于貝葉斯算法來(lái)構(gòu)建系統(tǒng)發(fā)育樹。利用DAMBE軟件[12]進(jìn)行堿基替代飽和度檢驗(yàn)（test of substitution saturation），以確定序列是否適合建樹，若建樹所用的序列沒(méi)有達(dá)到飽和，則可以建樹并進(jìn)行系統(tǒng)發(fā)育分析，若達(dá)到飽和則無(wú)建樹的必要。通過(guò)該軟件將比對(duì)后序列轉(zhuǎn)化成PAML格式用于選擇壓力分析。

1.2.2 系統(tǒng)發(fā)育樹的構(gòu)建 通過(guò)PAUP4.0[13]和Modeltest3.7[14]計(jì)算出符合該序列的核苷酸最佳替換模型。基于赤池信息量準(zhǔn)則（Akailke Information Criterion，AIC）選出的模型為TIM3+I+G，其參數(shù)結(jié)果用于MrBayes軟件[15]的計(jì)算，利用貝葉斯算法來(lái)構(gòu)建MC1R基因的系統(tǒng)發(fā)育樹。貝葉斯算法是基于最大似然法（Maximum Likelihood，ML）和貝葉斯推理法（Bayesian Inference，BI）的系統(tǒng)發(fā)生分析程序，其后驗(yàn)概率（posterior probabilities）通過(guò)3條熱鏈（heated chain）和1條冷鏈（cols chain）運(yùn)行1 000 000代估計(jì)，MCMC（Markov Chain Monte Carlo）分析中命令mcmc Samplefreq=100，其含義為每100代保存1棵樹，直至average standard deviation of split frequencies值小于0.01（該值代表2個(gè)獨(dú)立分析當(dāng)前的相似程度）可認(rèn)為達(dá)到平衡。丟棄老化樣本（sump burnin=2 500）后由剩余樣本來(lái)構(gòu)建合意樹（Majority consensus tree，Bayesian tree）。獲得的系統(tǒng)發(fā)育樹可通過(guò)Treeview軟件[16]來(lái)觀看。

1.2.3 正選擇位點(diǎn)的檢測(cè) Codeml程序在估算蛋白質(zhì)編碼序列同義替換和非同義替換速率以及檢測(cè)序列是否經(jīng)受正選方面被廣泛使用。本研究利用Codeml程序來(lái)進(jìn)行位點(diǎn)模型（site models）分析，以檢測(cè)基因在進(jìn)化過(guò)程中所經(jīng)受的選擇壓力的變化[13]。位點(diǎn)模型都是基于選擇壓力參數(shù)ω（編碼蛋白質(zhì)的核苷酸的非同義突變dN與同義突變dS的比值）來(lái)衡量分子進(jìn)化受到選擇壓力的方向和量度。若ω>1，即非同義突變被固定的速率遠(yuǎn)遠(yuǎn)大于同義突變，表明該位點(diǎn)受到正選擇壓力作用，位于正選擇壓力作用下的變異有利于新變異體的生存，同時(shí)該密碼子被稱為正選擇位點(diǎn)[17-19]，而ω=1和 ω<1 分別表示基因在進(jìn)化過(guò)程中經(jīng)歷了中性選擇和負(fù)選擇。其中，位點(diǎn)模型[20-21]假設(shè)不同位點(diǎn)的進(jìn)化速率不同，但在系統(tǒng)發(fā)育樹的不同分枝上沒(méi)有差別。參考以下6種模型：（1）模型M0（單一速率模型）假設(shè)不同位點(diǎn)進(jìn)化速率相同；（2）模型M1a（中性選擇模型）假設(shè)存在一部分保守位點(diǎn)，此時(shí)ω<1且其在序列中所占比例為p0，另一部分為中性選擇位點(diǎn)，此時(shí)ω=1且其在序列中所占比例為p1=1-p0；（3）模型M2a（正選擇模型）是在M1a模型的基礎(chǔ)上加入正選擇位點(diǎn)的比例即p2=1-p1-p0；（4）模型M3（離散模型）假設(shè)p0、p1和p2分別表示負(fù)選擇、中性選擇和正選擇，其對(duì)應(yīng)的選擇壓力參數(shù)分別是ω0、ω1和ω2；（5）模型M7（Beta分布模型）假設(shè)0<ω<1；（6）模型M8（Beta-ω分布模型）假設(shè)一部分位點(diǎn)的ω>1。利用似然比檢驗(yàn)（Likelihood ratio test，LRT）來(lái)比較零假設(shè)和備擇假設(shè)模型：M0與M3、M1a與M2a 和M7 與M8的最大似然值對(duì)數(shù)差的2倍（2ΔlnL），再由改良后的χ2分布進(jìn)行似然比檢驗(yàn)，由P值來(lái)判斷備擇假設(shè)是否成立。

2 結(jié)果與分析

2.1 不同物種MC1R基因的堿基替代飽和度檢驗(yàn)

因?yàn)檫M(jìn)化中轉(zhuǎn)換發(fā)生的頻率要高于顛換，更容易經(jīng)受替換和飽和效應(yīng)的影響，并且飽和效應(yīng)會(huì)直接影響到系統(tǒng)發(fā)育分析結(jié)果的可信度[22]，所以在系統(tǒng)發(fā)育分析中檢測(cè)堿基替代的飽和度是至關(guān)重要的。因?yàn)榫幋a序列中的三聯(lián)體密碼子第2個(gè)編碼位置通常是保守的，而第3個(gè)密碼子位置具有兼并性且對(duì)核苷酸同義替換有很大的貢獻(xiàn)。普遍認(rèn)為這是中性理論的微觀證據(jù)。如果第3個(gè)密碼子位置處替代達(dá)到飽和度，那么序列將會(huì)失去進(jìn)化信號(hào)，則以此編碼序列來(lái)研究基于進(jìn)化機(jī)制的結(jié)論可信度不高[23]。本研究使用DAMBE軟件在F84模型下檢測(cè)堿基替代的飽和性，分別利用密碼子3個(gè)位置處核苷酸轉(zhuǎn)化數(shù)（s）和顛換數(shù)（v）對(duì)遺傳距離作圖，從圖1可以看出，結(jié)果均接近線性回歸。表明堿基替換未達(dá)到飽和，因此可以進(jìn)行系統(tǒng)發(fā)育分析。

2.2 不同物種MC1R基因的系統(tǒng)發(fā)育分析

本研究在GenBank數(shù)據(jù)庫(kù)進(jìn)行搜索和篩選，共獲得21種動(dòng)物46條MC1R基因的編碼區(qū)序列：小家鼠、牦牛、人、野豬、家馬和家馬等等。利用貝葉斯方法構(gòu)建系統(tǒng)發(fā)育樹，利用原雞的MC1R基因編碼區(qū)序列作為外類群，建樹結(jié)果見(jiàn)圖2。

由MC1R基因構(gòu)建的系統(tǒng)發(fā)育樹主要分成2大分支。一支為牦牛、瘤牛、家牛、印度水牛、野豬、家馬、巖羚羊、羊駝、雙峰駝、中華白海豚、亞馬遜白海豚、矮鰭海豚、家貓、細(xì)腰貓14種物種，其中瘤牛、牦牛、家牛、印度水牛、巖羚羊、中華白海豚、亞馬遜白海豚和矮鰭海豚8個(gè)物種之間親緣關(guān)系較近，而家馬、羊駝、雙峰駝、野豬、家貓和細(xì)腰貓6個(gè)物種間親緣關(guān)系較近；另一支為小家鼠、加納小家鼠、亞洲家鼠、褐家鼠、拉布

拉多白足鼠、人6種物種，其中小家鼠、加納小家鼠、亞洲家鼠、褐家鼠和拉布拉多白足鼠5個(gè)物種之間親緣關(guān)系較近。原雞與哺乳動(dòng)物之間區(qū)別明顯。由此可見(jiàn)，動(dòng)物MC1R基因構(gòu)建的系統(tǒng)發(fā)育樹能很好的反應(yīng)不同物種之間的親緣關(guān)系。

2.3 正選擇位點(diǎn)分析

通過(guò)MrBayes軟件構(gòu)建動(dòng)物MC1R基因的系統(tǒng)發(fā)育樹后，利用PAML 4軟件分析MC1R基因在進(jìn)化過(guò)程中所經(jīng)受的選擇壓力（表2）。第1對(duì)模型M0和M3似然比檢驗(yàn)結(jié)果顯示2ΔlnL=169.42，df=4，P<0.01。因此模型M3顯著優(yōu)于模型M0，說(shuō)明位點(diǎn)間承受的選擇壓力具有異質(zhì)性。模型M3的3類密碼子位點(diǎn)所對(duì)應(yīng)的ω值分別為0.035 30，0.279 60和31.853 85，由于ω2對(duì)應(yīng)的p2為0，因此不存在正選擇位點(diǎn)。第2對(duì)模型M1a和M2a似然比檢驗(yàn)結(jié)果顯示 2ΔlnL=0，df=2，P=1.00。因此備擇假設(shè)模型M2a不成立，但并不推薦該模型為參考標(biāo)準(zhǔn)。第3對(duì)模型M7和M8似然比檢驗(yàn)結(jié)果顯示2ΔlnL=6.18，df=2，P=0.045。因此備擇假設(shè)模型M8成立，其ω=1.501 4>1，表明MC1R基因受到選擇壓力，存在正選擇位點(diǎn)（29W，191L，253L，304L）。

本研究用于檢測(cè)正選擇位點(diǎn)的3對(duì)模型發(fā)現(xiàn)MC1R基因編碼的氨基酸位點(diǎn)在進(jìn)化過(guò)程中所經(jīng)受中性選擇和純化選擇作用比例較高，但有4個(gè)氨基酸進(jìn)化受到正選擇作用的驅(qū)動(dòng)。

3 討論與結(jié)論

在進(jìn)化過(guò)程中，基因的復(fù)制過(guò)程所產(chǎn)生的遺傳變異會(huì)造成功能分歧，隨后純化選擇將功能固定下來(lái)。少數(shù)位點(diǎn)甚至單個(gè)位點(diǎn)的核苷酸替代都能改變其編碼的蛋白質(zhì)功能。而整個(gè)基因大多數(shù)位點(diǎn)都處于純化選擇或中性選擇的進(jìn)化狀態(tài)，以至于少數(shù)位點(diǎn)的正選擇信息會(huì)被大多數(shù)位點(diǎn)所稀釋，從而導(dǎo)致正選擇信息被整體湮沒(méi)[24]。現(xiàn)代分子適應(yīng)性進(jìn)化分析多以dN/dS的比值來(lái)衡量選擇壓力。若ω>1，且兩兩模型之間差異顯著，則表明編碼序列在對(duì)應(yīng)的位點(diǎn)經(jīng)受正選擇作用。與傳統(tǒng)方法相比，這種基于最大似然法所估計(jì)的正選擇檢測(cè)方法，可把長(zhǎng)期由純化選擇和中性選擇作用掩蓋下的正選擇位點(diǎn)識(shí)別出來(lái)，反應(yīng)出更高的靈敏度[25]。

本研究采用位點(diǎn)模型來(lái)檢測(cè)選擇壓力，結(jié)果表明，在21種動(dòng)物中只檢測(cè)到了4個(gè)正選擇位點(diǎn)，可見(jiàn)MC1R基因編碼的蛋白從整體水平上主要受中性漂變和純化選擇的作用。推測(cè)可能是由以下2個(gè)原因：（1）MC1R基因所編碼的蛋白功能重要，氨基酸序列較為保守。（2）MC1R基因序列進(jìn)化歷史久遠(yuǎn)，且其較短的序列只能提供較少的可檢測(cè)位點(diǎn)，這樣也難以檢測(cè)出適應(yīng)性進(jìn)化信號(hào)[26]。

適應(yīng)性進(jìn)化從整體水平上了解了動(dòng)物MC1R基因在進(jìn)化過(guò)程中所經(jīng)受的選擇作用，同時(shí)也分析了不同物種之間親緣關(guān)系的遠(yuǎn)近。對(duì)于親緣關(guān)系較近的物種，選擇壓力小；而親緣關(guān)系較遠(yuǎn)的物種選擇壓力較大且基因序列較保守。MC1R基因作為控制動(dòng)物毛色形成的主效基因，經(jīng)證實(shí)其變異可作為許多哺乳動(dòng)物黑、紅毛色性狀的變化的主要原因[27]，該基因功能改變會(huì)表現(xiàn)出褐黑色素的紅色或黃色性狀，反之會(huì)表現(xiàn)出真黑色素的黑色或棕色性狀[28]。

本研究針對(duì)MC1R基因編碼區(qū)進(jìn)行分子進(jìn)化分析，發(fā)現(xiàn)MC1R基因編碼的蛋白從整體水平上主要受中性漂變和純化選擇的作用，對(duì)于正選擇作用的分析，共檢測(cè)到4個(gè)正選擇位點(diǎn)，這些位點(diǎn)對(duì)探討MC1R基因編碼蛋白功能的改變提供了重要依據(jù)，從分子進(jìn)化的角度為MC1R基因調(diào)控黑色素的形成提供了一個(gè)新的思路。

參考文獻(xiàn)：

[1]甘海云，李建斌，王紅梅，等. 黑素皮質(zhì)激素受體1基因與牛的毛色[J]. 家畜生態(tài)學(xué)報(bào)，2007，28（3）：99-102.

[2]楊永升，李 寧，鄧學(xué)梅，等. 黑素皮質(zhì)素受體1——哺乳動(dòng)物黑色素形成中的關(guān)鍵基因[J]. 遺傳，2004，26（4）：544-550.

[3]Sazanov A，Masabanda J，Ewald D，et al. Evolutionarily conserved telomeric location of BBC1 and MC1R on a microchromosome questions the identity of MC1R and a pigmentation locus on chromosome 1 in chicken[J]. Res Dec，1998，8：651-654.

[4]Hunt G，Rot H J，van Oost B A. Identification of a premature stop codon in the melanocyte-stimulating hormone receptor gene（MC1R） in Labrador and Golden ret rievers with yellow coat color[J]. Animal Genetics，2000，31（3）：194-199.

[5]Boissy R E，Nordlund J J，Rheins L A. Molecular basis of congenital hypopigmentary disorders in humans：a review[J]. Pigment Cell Research，1997，10：21-24.

[6]Leal-Klevezas D S，Martínez-Vázquez I O，López-Merino A，et al. Single-step PCR for detection of Brucella spp. from blood and milk of infected animals[J]. Journal of Clinical Microbiology，1995，33（12）：3087-3090.

[7]Shimizu K K，Purugganan M D. Evolutionary and ecological genomics of arabidopsis[J]. Plant Physiology，2005，138（2）：578-584.

[8]Nei M，Kumar S. Molecular evolution and phylogenetics[M]. New York：Oxford University Press Inc，2000.

[9]Suyama M，Torrents D，Bork P. PAL2NAL：robust conversion of protein sequence alignments into the corresponding codon alignments[J]. Nucleic Acids Research，2006，34：W609-W612.

[10]Hall T A. Bioedit：a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT[J]. Nucleic Acids Symp Ser，1999，41：95-98.

[11]Tamura K，Dudley J，Nei M，et al. MEGA4：molecular evolutionary genetics analysis（MEGA）software version 4.0[J]. Molecular Biology and Evolution， 2007，24：1596-1599.

[12]Xia X，Xie Z. DAMBE：data analysis in molecular biology and evolution[J]. Journal of Heredity，2001，92：371-373.

[13]Yang Z H. PAML 4：phylogenetic analysis by maximum likelihood[J]. Molecular Biology and Evolution，2007，24（8）：1586-1591.

[14]Posada D，Crandall K A. MODELTEST：testing the model of DNA substitution[J]. Bioinformatics，1998，14（9）：817-818.

[15]Ronquist F，Huelsenbeck J P. MrBayes 3：bayesian phylogenetic inference under mixed models[J]. Bioinformatics，2003，19（12）：1572-1574.

[16]Page M R D. TreeView：an application to display phylogenetic trees on personal computers[J]. Computer Applications in the Biosciences，1996，12（4）：357-358.

[17]Choisy M，Woelk C H，Gue G J，et al. Comparative study of adaptive molecular evolution in different human immunodeficiency virus groups and subtype[J]. Journal of Virology，2004，78（4）：1962-1970.

[18]Zhang C Y，Ding N，Chen K P，et al. Complex positive selection pressures drive the evolution of HIV-1 with different co-receptor tropisms[J]. Science China：Life Sciences，2010，53（10）：1204-1214.

[19]許慧琳，張文彤，趙耐青，等. H2N2亞型人甲型流感病毒HA1序列進(jìn)化正選擇位點(diǎn)研究[J]. 中華流行病學(xué)雜志，2007，28（4）：385-389.

[20]Nielsen R，Yang Z. Likelihood models for detecting positively selected amino acid sites and application to the HIV-1 envelope gene[J]. Genetics，1998，148（3）：929-936.

[21]Yang Z. Maximum likelihood estimation on large phylogenies and analysis of adaptive evolution in human influenza virus A[J]. Journal of Molecular Evolution，2000，51（5）：423-432.

[22]呂雪梅，王應(yīng)祥，張亞平. 蜂猴線粒體細(xì)胞色素b基因變異特點(diǎn)及系統(tǒng)發(fā)育分析[J]. 動(dòng)物學(xué)研究，2001，22（2）：93-98.

[23]梁國(guó)明，劉桂瓊，姜?jiǎng)灼? 反芻動(dòng)物PRNP基因進(jìn)化機(jī)制研究[C]//中國(guó)畜牧獸醫(yī)學(xué)會(huì)養(yǎng)羊?qū)W分會(huì)全國(guó)養(yǎng)羊生產(chǎn)與學(xué)術(shù)研討會(huì)議論文集，2010.

[24]姬廣超，王明輝，高會(huì)江，等. 哺乳動(dòng)物MT基因的進(jìn)化選擇與功能分歧[J]. 東北農(nóng)業(yè)大學(xué)學(xué)報(bào)，2010（9）：82-88.

[25]Yang Z H. Inference of selection from multiple species alignments[J]. Current Opinion in Genetics & Development，2002，12（6）：688-694.

[26]Yang Z. Computational molecular evolution[M]. New York：Oxford University Press，2006.

[27]Schmutz S M，Moker J S，Berryere T G，et al. An SNP is used to map MC1R to dog chromosome 5[J]. Animal Genetics，2001，32（1）：43-44.

[28]Newton J M，Wilkie A L，He L，et al. Melanocortin 1 receptor variation in the domestic dog[J]. Mammalian Genome，2000，11（1）：24-30.