大豆MIR156基因家族起源与进化模式研究

朱荣胜，金会会，王志豪，陈庆山（.东北农业大学理学院，哈尔滨　5000；.东北农业大学生命科学学院，哈尔滨　5000；.东北农业大学农学院，哈尔滨　5000）



大豆MIR156基因家族起源与进化模式研究

朱荣胜1，金会会2，王志豪2，陈庆山3
（1.东北农业大学理学院，哈尔滨150030；2.东北农业大学生命科学学院，哈尔滨150030；3.东北农业大学农学院，哈尔滨150030）

摘要：MicroRNA（miRNA）是小分子RNA，参与靶基因转录后调控，长度通常为20~22 nt。MIR156基因在植物中分布广，是目前已知第一个与植物年龄有关的小分子RNA，表达量随植物年龄增长而逐渐减少，对植物生长发育起重要调控作用。研究基于古多倍体大豆基因组共线性分析，结合miRBase收录数据和预测获取的基因，从进化角度对该基因家族进行系统分析。通过对大豆MIR156基因家族成员分布、复制模式、扩张模式和分子系统发育分析，揭示这一古老基因家族基本起源模式。结果表明，该家族由多个祖先基因起源，以全基因组复制及大片段复制等方式发展。研究从方法和结论上对系统分析植物小分子RNA起源和进化具有重要意义。

关键词：基因家族；复制模式；Ks值；大豆

朱荣胜,金会会,王志豪,等.大豆MIR156基因家族起源与进化模式研究[J].东北农业大学学报, 2016, 47(2): 8-15.

Zhu Rongsheng, Jin Huihui, Wang Zhihao, et al. Origin and expansion model research on MIR156 gene family in soybean [J]. Journal of Northeast Agricultural University, 2016, 47(2): 8-15. (in Chinese with English abstract)

豆科植物物种形成于约59 MYA（Million years ago）基因组倍增时期[1]，植物miRNA通过RNA聚合酶转录，识别与其序列互补的靶mRNA并介导翻译抑制或剪切，实现对靶基因的转录后调控[2-3]。由于有些miRNA靶基因具有调控功能的转录因子[4]，miRNA对豆科植物生长发育[5]和胁迫应答[6]具有重要调控功能。

全基因组测序认为大豆是古多倍体[7]。全基因组复制Whole genome duplication（WGD）大片段复制，串联重复被视为进化新成员的重要遗传物质来源方式[8]。WGD由多倍体化引起，是植物界常见现象[9]，极大推动隔离再生及区域物种形成。高等真核生物串联重复发生频繁，表现为一簇基因形成且物种之间拷贝数变异[10]。Lisch等研究表明，大量基因及基因片段可被TEs（Transposable elements）捕获并移动到其他位置，作为产生和扩充新基因机制[11]，WGD对基因进化和变革影响最大。理论上，多拷贝后复制基因的功能冗余会经历几百万年的宽松选择，导致有害突变积累或最终淘汰一个复制拷贝[8]。复制基因的功能分化为：①均保留祖先基因功能；②一个拷贝基因功能保留，另一个基因亚功能化；③两个基因功能发生改变。研究表明调控网络的重构对复制基因进化结果产生重要影响[8]。

MIR156基因家族分别在拟南芥[4]和水稻生长发育中起重要调节作用，主要靶向SPL基因家族[12]。miRNA基因在植物进化进程中是不断丢失与获得的过程[13]。依据家族分布广泛性，可分为三大类：较为古老的一些基因家族；中度保守的基因家族；低度保守的基因家族[14]。植物miRNA基因家族成员在基因扩张后经历明显重排[15]，多呈现分散分布。数据单双子叶植物MIR156家族成员分布，该家族随植物全基因组倍增及其他复制活动，改变成员分布和状况[16]。

目前miRNAs基因起源和进化机制未见系统研究。靶基因可生成miRNA基因[17]；转座子具有反向重复序列，可形成发夹结构，有成为miRNA基因潜在基础[18]；新miRNA基因通过已有miRNA基因复制[19]；也可能源于基因组发夹结构序列的随机突变[20]。来源于靶基因的反向复制是其主要方式[19]。植物miRNA基因每百万年新生2~3个[21]，保守性较高。源于已有的miRNA假说考虑复制之后引起的剂量平衡问题[15]。Nozawa等通过对多个植物基因组miRNA分析，发现miRNA基因主要起源于已存在的miRNA基因或者蛋白编码基因的复制[13]。miRNA随机序列假说表明，这种序列具有发夹结构，中间部分区段丢失使反向复制序列及中间序列符合miRNA前体特征，及后续序列突变，成为有功能的miRNA可能性较大[20]。本研究基于大豆基因组共线性，通过对MIR156基因家族的成员复制、基因家族的扩张模式、倍增时间、系统发育研究，阐明MIR156家族扩张过程及成员间的分化，为进一步研究大豆MIR156基因的起源和功能提供依据。特别是通过对复制基因对和共线性模块挖掘，得到一些候选miRNA基因，为后续的功能验证和试验提供理论基础。

1　材料与方法

1.1数据来源及前体序列补充

大豆基因组数据来源于Phytozome数据库（http: //www.phytozome.com），基因家族数据来源于miRBase数据库（http://www.mirbase.org）。

1.2候选前体序列筛选与获得

从miRBase上得到大豆MIR156成体序列并去重复，得到6条差异成体序列，并以此序列作为种子序列与大豆基因组blast比对，window_size设置为20，word_size设置为12，其他默认。

通过前体序列二级结构预测并结合计算最小自由能，分析大豆miRNA序列及MIR156序列前体长度，通过最小自由能指数（MFEI）等参考指标共得到21个新的高可信序列前体，作为研究数据补充。同时采纳Zhou等部分数据作为补充[22]。

1.3大豆基因组共线性分析

真核生物基因组不同是因为基因同线性和基因排序后共线性不同[23]。尽管植物全基因复制频繁发生，与染色体数目相比基因组规模更趋于稳定，通常都在50倍范围内[23]。本研究以大豆Williams 82基因组v1.0作为研究基因组，结合MCScanX软件[24]及ClustalX对其基因组共线性进行分析。其中MCScanX软件以编码基因序列分析为对象，形成一个共线性区域最低要求5个基因。

1.4基因倍增时间估算

对倍增块中保守的蛋白质编码基因进行氨基酸序列比对，然后编码序列比对；Ks为替换次数/同义替换位点个数，计算方法采用Perl语言编写的Bio::Align::DNAStatistics模块里的calc_KaKs_pair方法[25]。miRNA基因作为非编码基因，采用所在位置临近基因平均Ks值代替计算分化时间[26]。其中单个基因Ks值及基因组共线性模块Ks值作比较[27]，约13 MYA，是大豆形成过程中一次全基因组复制[7]。由于序列经历多次复制，碱基替换率改变，目前采用平均Ks值[25]用于计算分化时间T，计算公式为：T=Ks/2E，E表示分子替换速率[25]。化石证据与进化分析结果表明大豆分子替换速率约为拟南芥的1/5，拟南芥E值为1.5×10-8，大豆E值为6.1× 10-9，E值单位为每个密码子每年突变频率[27]。

图1前体序列筛选流程Fig. 1 Screening process of precursor sequence

1.5家族成员复制

大豆物种形成后伴随多次基因组大片段复制与重组[7]。基因复制通过不等交换、逆转录转座或全基因组重复产生一个与原基因相似的基因或碱基序列[28]。植物基因家族成员数目比动物多[19]，主要因其在植物中特异性扩张造成，物种特异性扩张主要有片段复制、串联重复[26]和逆转录转座[18]，其中片段重复和串联重复为主要方式。重复的基因经历正选择和基因转换等适应性进化的选择压力[29]。

1.6大豆MIR156基因家族进化树构建及分类

对家族成员基因序列进行多序列比对，并用MEGA 6.0软件[30]采用邻接法（Neighbor-Joining method，N-J）（bootstrap[31]设置为1 000）构建系统进化树。最大似然法（Maximum likelihood，M-L）构建系统发生树，以最为合适的基因树作为家族成员进化关系分析树。

2结果与分析

2.1靶序列长度选取

通过大豆在miRNA序列前体GC[35]含量对比分析、最小自由能指数值分布、前体序列长度分布统计，表明大豆miRNA基因GC含量普遍低于40%，序列自由能选取范围小，结构稳定，整体性强。植物前体序列长度及变异容忍范围变化较大，植物前体长度和二级结构有较大差异[32]。通过对大豆miRNA前体序列统计发现，miRNA前体序列多落在离成体序列上下游150 bp以内，因此采取比对靶序列上下游150 bp作为研究序列。图2中为GC含量及最小自由能分布曲线。

2.2作为非编码的miRNA基因的Ks值选取

图3是共线性区块及上下游20 kb包含基因分别求平均值的Ks比较。采用不同范围取值，获得代表miRNA编码基因的近似Ks值。位于基因内取宿主基因值，位于基因间隙的取上下游20 kb距离包含蛋白编码基因的平均Ks值。通过MIR156基因家族17个成员的Ks值两种取值方法比较，位于共线性模块的miRNA基因所在模块整体Ks平均值趋于平稳，而上下游20 kb小范围取值较为合理。结果表明至少家族成员临近有1个基因，最多有8个，平均3.7个。

2.3家族成员相关信息

表1显示大豆MIR156基因家族成员物理位置信息，复制来源信息，位置精细定位信息和估算的分化时间。该家族66%基因位于基因间隙，对位于基因内部miRNA基因，66%位于5'UTR，2个位于exon，1个位于3'UTR。对于新预测部分显示大多位于外显子，有4个成员位于基因间隙。由Ks值可知多数miRNA基因分化于最后一次大豆全基因组复制事件前后；相对较高值拷贝基因趋向于祖先基因。

2.4大豆MIR156位置分布及复制情况

从图4上看，该基因家族成员多呈散状分布，符合植物保守家族特点。相比动物基因家族，其多呈成簇分布，而植物基因家族大多与其他基因家族成员成簇分布[14]。

图2靶序列长度选取参考数据Fig. 2 Target sequence length choose reference data

图3不同方法下Ks值比较Fig. 3 Ks comparation in different methods

表1 MIR156基因家族成员及预测前体信息Table 1 MIR156 gene family members and predict precursor sequence information

续表

由共线性对应关系可知，整个家族可分为A、B、C、D、E五部分。A组成员有n10、n12、n8、n14、x、v、n、y、o、h、w、k、d、i、l、m，结果显示，i、n10与d、n14f反向对应，n14又与n8、n12共线性对应，x、n区域分别与v、h、n8、w、k对应，其中Gm17上的o/y与Gm14上的x/n反向对应，x/n与w/k正向对应，Gm06上的w/k与Gm04上的v/h正向对应。gma-MIR156n/o的Ks值约为1.024，gma-MIR156k/h的Ks值约为1.006，其他成员Ks小于0.3，推测k/h源于n/o复制，之后两个祖先基因先后发生一次拷贝，在大豆形成后各自又发生一次大片段复制，形成目前对应关系。B组成员有r、c、n19、n21，其中n19对应部分基因丢失，相比其他成员，r、c成体序列更接近预测序列，但由分化时间可知，这两个成员是后期分化丢失遗留，并非祖先基因。C组成员为j、e、q、s，由分化时间可知，时间较近，e的分化时间略早。D组成员为u、a、n18、n7、n12。E组是单独的两个共线性对应，分别为p与t、b与f。其中共线性对应区域成员丢失的有n13、g，且多为反向对。由Ks值分布及其对应全基因分化时间可知，A组h、k、o疑似为祖先基因，且来自于一个更早的祖先基因，E组b、f起源较早，但由对应关系及系统发生树聚类可知，对该家族成员增加贡献不大，B、C、D组祖先基因丢失。

2.5 MIR156基因家族系统发生分析

最终采用临近法构建基因树作为分析指导树。结果显示序列分为三个明显分支。由共线性对应关系可知，n4、n11、n6由大片度复制产生。B组r、c相对于其他成员距离较远。g和n21从系统发生树及位置关系判断属于串联复制（见图5）。

图4 MIR156基因家族成员分布及对应关系Fig. 4 MIR156 gene family member distribution and relationship

图5前体序列系统发生树Fig. 5 Precursor sequence phylogenetic tree

3讨论

植物miRNA作用位点无明显特征，匹配程度相对较高[33]。这种作用机制在动物和植物间差异较大，与其miRNA基因主要起源方式有关，不同起源和进化方式最终导致miRNA基因靶作用基因机制不同[13]。在大豆中，保守性高的MIR156基因家族成员散布在整个基因组，成员跨越不同组件[8]。84%基因位于基因间隙，与动物略为不同[13]，基因间隙基因个数与基因组大小无相关性，仅少数部分落在CDS区。由于多数TE-like miRNA基因家族物种特异，TEs不可能是miRNA基因家族主要起源方式[13]。植物miRNA基因家族扩张有共同的祖先[19]。本研究采用编码基因替代非编码基因Ks值，论证取值方式和范围。对复制深度的阈值界定问题存在争议。认定一个起源较早基因规定在相对时间跨度区域。

4结论

本文分析共线性背景下家族成员对应关系及复制时间，结果表明，成员号为f、h、k、b的基因起源较早，认为是祖先基因代表。该基因家族分析阐明大豆基因组发展历程，有利于认识大豆物种形成机制。非编码多拷贝基因起源机制需更优算法及分析模式确定，对编码基因研究具有参考意义。

[参考文献]

[ 1 ] Schlueter J A, Dixon P, Granger C, et al. Mining EST databases to resolve evolutionary events in major crop species[J]. Genome, 2004, 47(5): 868-876.

[ 2 ] Wang L, Cao C, Ma Q, et al. RNA-seq analyses of multiple meristems of soybean: Novel and alternative transcripts, evolutionary and functional implications[J]. BMC plant biology, 2014, 14(6):150-169.

[ 3 ] Bartel D P. MicroRNAs: Genomics, biogenesis, mechanism, and function[J]. Cell, 2004, 116(2): 281-297.

[ 4 ] Wang J W, Czech B, Weigel D. miR156-regulated SPL transcription factors define an endogenous flowering pathway in Arabido-psis thaliana[J]. Cell, 2009, 138(4): 738-749.

[ 5 ] Zhang B, Pan X, Cobb G P, et al. Plant microRNA: A small regulatory molecule with big impact[J]. Developmental Biology, 2006, 289(1): 3-16.

[ 6 ] Sunkar R, Li Y F, Jagadeeswaran G. Functions of microRNAs in plant stress responses[J]. Trends in Plant Science, 2012, 17(4): 196-203.

[ 7 ] Schmutz J, Cannon S B, Schlueter J, et al. Genome sequence of the palaeopolyploid soybean[J]. Nature, 2010, 463(7278): 178-183.

[ 8 ] Zhao M, Meyers B C, Cai C, et al. Evolutionary patterns and coevolutionary consequences of MIRNA genes and microRNA targets triggered by multiple mechanisms of genomic duplications in soybean[J]. The Plant Cell, 2015, 27(3): 546-562.

[ 9 ] Jiao Y, Wickett N J, Ayyampalayam S, et al. Ancestral polyploidy in seed plants and angiosperms[J]. Nature, 2011, 473(7345): 97-100.

[10] Lawrence J. Selfish operons: The evolutionary impact of gene clustering in prokaryotes and eukaryotes[J]. Current Opinion in Genetics & Development, 1999, 9(6): 642-648.

[11] Lisch D. How important are transposons for plant evolution? [J]. Nature Reviews Genetics, 2013, 14(1): 49-61.

[12] Preston J C, Hileman L C. Functional evolution in the plant SQUAMOSA-PROMOTER BINDING PROTEIN-LIKE (SPL) gene family[J]. Frontiers in Plant Science, 2013, 4(4):80-93.

[13] Nozawa M, Miura S, Nei M. Origins and evolution of microRNA genes in plant species[J]. Genome Biology and Evolution, 2012, 4 (3): 230-239.

[14] Zhang B, Pan X, Cannon C H, et al. Conservation and divergence of plant microRNA genes[J]. The Plant Journal, 2006, 46(2): 243-259.

[15] Severin A J, Cannon S B, Graham M M, et al. Changes in twelve homoeologous genomic regions in soybean following three rounds of polyploidy[J]. The Plant Cell, 2011, 23(9): 3129-3136.

[16] Wang S, Zhu Q H, Guo X, et al. Molecular evolution and selection of a gene encoding two tandem microRNAs in rice[J]. FEBS Letters, 2007, 581(24): 4789-4793.

[17] Allen E, Xie Z, Gustafson A M, et al. Evolution of microRNA genes by inverted duplication of target gene sequences in Arabidopsis thaliana[J]. Nature Genetics, 2004, 36(12): 1282-1290.

[18] Piriyapongsa J, Jordan I K. Dual coding of siRNAs and miRNAs by plant transposable elements[J]. Rna, 2008, 14(5): 814-821.

[19] Li A, Mao L. Evolution of plant microRNA gene families[J]. Cell Research, 2007, 17(3): 212-218.

[20] Felippes F F, Schneeberger K, Dezulian T, et al. Evolution of Arabidopsis thaliana microRNAs from random sequences[J]. Rna, 2008, 14(12): 2455-2459.

[21] Fahlgren N, Jogdeo S, Kasschau K D, et al. MicroRNA gene evolution in Arabidopsis lyrata and Arabidopsis thaliana[J]. The Plant Cell, 2010, 22(4): 1074-1089.

[22] Zhou Z, Wang Z, Li W, et al. Comprehensive analyses of microRNA gene evolution in paleopolyploid soybean genome[J]. The Plant Journal, 2013, 76(2): 332-344.

[23] Tang H, Bowers J E, Wang X, et al. Synteny and collinearity in plant genomes[J]. Science, 2008, 320(5875): 486-488.

[24] Wang Y, Tang H, Debarry J D, et al. MCScanX: A toolkit for detection and evolutionary analysis of gene synteny and collinearity[J]. Nucleic Acids Research, 2012, 40(7): 35-49.

[25] Lynch M, Conery J S. The evolutionary fate and consequences of duplicate genes[J]. Science, 2000, 290(5494): 1151-1155.

[26] Maher C, Stein L, Ware D. Evolution of Arabidopsis microRNA families through duplication events[J]. Genome Research, 2006, 16(4): 510-519.

[27] Fan C, Wang X, Hu R, et al. The pattern of phosphate transporter 1 genes evolutionary divergence in Glycine max L. [J]. BMC Plant Biology, 2013, 13(4):48-64.

[28] Magadum S, Banerjee U, Murugan P, et al. Gene duplication as a major force in evolution[J]. Journal of Genetics, 2013, 92(1): 155-161.

[29] Fahlgren N, Howell M D, Kasschau K D, et al. High-throughput sequencing of Arabidopsis microRNAs: evidence for frequent birth and death of MIRNA genes[J]. PloS one, 2007, 2(2): 219-233.

[30] Tamura K, Stecher G, Peterson D, et al. MEGA6: Molecular evolutionary genetics analysis version 6.0[J]. Molecular Biology and Evolution, 2013, 30(12): 2725-2729.

[31] Sanderson M J, Wojciechowski M F. Improved bootstrap confidence limits in large-scale phylogenies, with an example from Neo-Astragalus (Leguminosae) [J]. Systematic Biology, 2000, 49 (4): 671-685.

[32] Axtell M J, Bowman J L. Evolution of plant microRNAs and their targets[J]. Trends in Plant Science, 2008, 13(7): 343-349.

[33]曲波,甄贞,仇有文,等.基于生物信息学方法挖掘奶山羊miRNAs研究[J].东北农业大学学报, 2015, 46(1): 86-93.

Origin and expansion model research on MIR156 gene family in soybean

ZHU Rongsheng1, JIN Huihui2, WANG Zhihao2, CHEN Qingshan3(1. School of Science, Northeast Agricultural University, Harbin 150030, China; 2. School of Life Sciences, Northeast Agricultural University, Harbin 150030, China; 3. School of Agriculture, Northeast Agricultural University, Harbin 150030, China)

Abstract:MicroRNAs (miRNA) were small RNAs participating in target gene transcription regulation which length was usually 20 to 22 nt. MIR156 gene family widely distributed in plants. MIR156 was currently known molecular markers with relation to plant age. The level of expression of this gene family gradually reduces along by age rising. It played an important role in plant growth and development regulation and control. From evolution view, a systemic analysis was implemented to this miRNA family based on soybean genome collinearity analysis. Some of MIR156 family data derived from miRBase database and others were predicted based on computation methods. The ancient family basic origin model was revealed through MIR156 gene family members distribution, reproduction, expansion and molecular phylogenetic analysis. The results showed that the family was originated from several ancestral members and it expanded through the whole gene replication and large fragment duplication. From results and conclusions, this research was significant to analyze systemically the origin and evolution of plant small molecular RNA.

Key words:gene family; replication model; Ks value; soybean

作者简介：朱荣胜（1975-），男，副教授，博士，硕士生导师，研究方向为生物信息学。E-mail: 007zrs@163. com

基金项目：东北农业大学研究生科技创新资助项目（yjscx14053）

收稿日期：2015-07-17

中图分类号：S565.1

文献标志码：A

文章编号：1005-9369（2016）02-0008-08